A numerical investigation of sinkage and trim effects on the resistance of trimaran hull form in calm water

对于压差阻力占总阻力比例较大的高速船舶，细长的船型可以显著降低阻力^[1]。但是，若船体过细，其横稳性会变差。而在主船体两侧各设置一个小的侧体，则可在提高航速的同时大大改进横稳性。近20多年来，备受关注的三体船概念^[2-3]就充分考虑了这一因素^[4-8]，三体船的应用带来了船型设计上的变革。澳大利亚Austal公司的Benchijigua Express三体渡船（船长L=126.7 m）和美国海军“独立”号三体濒海战斗舰（船长L=127.6 m）的投入使用，激起了人们对大型高速三体船研究的兴趣^[2-3]。Benchijigua Express三体渡船和“独立”号濒海战斗舰的速度均高达40 kn以上，长度超过100 m，其显著特征是弗劳德数Fr > 0.6，处于半滑行状态。与排水型船型（一般Fr < 0.3）或滑行艇型（一般Fr > 1.0）相比，半滑行状态由于动升力的作用，静浮力和浮力分布发生了变化，致使船体的航行姿态（升沉和纵倾）发生了显著变化，进而湿表面积和兴波状况也产生了相应的变化，带来的效果就是阻力性能发生显著变化^[9-13]。半滑行船型的阻力性能不同于排水型船型和滑行艇型，对于高速半滑行三体船型，其航行姿态对阻力的影响是阻力预报方面必不可少的一个重要研究内容。

造船界很早就意识到了航行姿态对阻力的影响。Yeung^[14]最早从船舶兴波阻力薄船理论的角度开展了升沉和纵倾对阻力影响的研究。航行姿态对三体船阻力的影响，始于Yang和Löhner^[9]、王中等^[10]及李柯等^[15]的研究，不同于Yeung的早期薄船理论，他们分别使用更加准确的Rankine源法或细长体理论考察了这一问题。数值结果表明，在0.3<Fr<0.6范围内，航行姿态的影响显著。同时，结果也表明，势流理论在更高速时与试验结果存在较大偏差，其原因是势流理论在数学模型和数值模拟等方面均存在误差，因而会产生一些不真实的现象。

随着计算机技术和船舶水动力学的迅猛发展，近年来，较多学者开始倾向于采用CFD方法考察这方面的问题^[11-13]。本文将以一艘方艉折角线型的前三体船为例，采用CFD方法重点探讨航行姿态对阻力成分（摩擦阻力和压差阻力）及不同成分阻力大小的影响规律。

浮力的本质是压力沿船体表面的积分。根据伯努利方程，当航速较大时，压力就会降低，从而沿船体表面积分的浮力也会发生变化。对于排水型船，当航速较低时，航行中压力改变会很小，故浮力的变化通常可忽略不计。因此，排水型船（体积弗劳德数满足${{F}}{{{r}}_\nabla }{\rm{ < 1}}{\rm{.0}}$，一般对应的长度弗劳德数满足Fr < 0.3）的航行姿态影响通常忽略不计；而对于滑行船型（体积弗劳德数满足${{F}}{{{r}}_\nabla }{\rm{ > 3}}{\rm{.0}}$，对应的长度弗劳德数Fr >1.0），因升力是主要的托举力，故航行姿态由升力决定。对于半滑行船型（体积弗劳德数满足$ 1.0 < {Fr}_\nabla < 3.0 $，一般对应的长度弗劳德数满足0.3~0.4< Fr <1.0），由于航速变化带来的压力变化、进而浮力以及航行姿态的变化相当显著，因此必须予以考虑，这也是航行姿态对半滑行船型（或者过渡船型）的阻力预报影响较大的原因。三体船型属高速船型，处于半滑行状态。为分析航行姿态对阻力的影响，本文以某三体船为例，考虑三体船的2种升沉和纵倾状态，即自由航行状态（释放三体船的升沉和纵倾运动，简称“自由模型”）和固定航行姿态的状态（约束三体船的升沉和纵倾运动，简称“约束模型”），采用CFD方法对各自的阻力进行了计算。

绕三体船的流动应满足连续性方程和动量方程（N-S方程）。假定流体为不可压缩流体，流体域内时均形式的连续性方程和动量方程用张量形式表示为^[16-17]：

式中：${u_i},{u_j}$为速度分量的时均值（$i,j = 1,2,3$）；$t$为时间；${x_i},{x_j}$为三维直角坐标系的3个坐标轴（$i,j = 1,2,3$）；$p$为压力时均值；$\rho $为流体密度；$\mu $为流体动力黏性系数；${g_i}$为重力加速度分量（$i = 1,2,3$）；$\rho \overline {u_i'u_j'} $为雷诺应力项。

本文采用标准$k - \varepsilon $湍流模型来封闭RANS方程。湍动能$k$、湍动耗散率$\varepsilon $的输运方程分别如式（3）和式（4）所示^[16-17]：

式中：${\nu _{\rm{t}}}$为湍动黏度，定义为${\nu _{\rm{t}}} = \rho {C_\mu }\dfrac{{{k^2}}}{\varepsilon }$；$\nu $为流体运动黏性系数。通常情况下，标准$k - \varepsilon $湍流模型的系数取值如下^[16-17]：

本次计算采用的是一个前三体构型的三体船模型^[18]。该船模的主尺度如表1所示，横剖面示意图如图1所示。

Figure 1.  Schematic diagram of transverse section of trimaran model

三体船的构型布局如图2所示。图中：$a$为侧体中心位置与主体中心位置的纵向距离，侧体中心位于主体舯前为正；$b$为侧体中心位置与主体中心位置的横向距离，$b$值始终为正；L为主体长度；L_s为侧体长度。本次计算的三体船模型对应的构型参数a=1.0 m，b=0.7 m。

Figure 2.  Schematic diagram of trimaran layout

借助CFD商业软件STAR-CCM+完成数值模拟，采用有限体积法进行空间离散，流体方程则采用非定常求解器求解。所选$k - \varepsilon $湍流模型选择适用于复杂几何外部的流动问题求解。自由液面的捕捉采用流体体积（volume of fluid，VOF）法。均匀来流采用静水VOF法处理。三体船的六自由度运动模拟利用DFBI (dynamic fluid body interaction) 求解完成，自由运动模型释放了纵摇和垂荡这2个方向的自由度^[16-17]。

本文的数值计算采用重叠网格技术，将计算域分成背景区域和重叠网格区域两部分。使用2套网格进行网格划分，计算总网格数240万。其中，在主体和2个侧体周围以及自由表面处均进行网格加密。

本次计算的流体域尺度以主体长度L为基准选取。流体域长4.5L，宽3L，高2L，其中水域深度为L，空气域高度为L。计算流体域的边界条件由5个速度入口和1个压力出口组成。主体的船艏距离速度入口的距离为L。

数值计算模型的具体边界条件设置如图3所示。

Figure 3.  The boundary conditions

本次验证的模拟工况为：Fr = 0.546（航速V = 3.416 m/s）。分别求时间步长Δt = 0.01，0.02和0.04 s这3种条件下，三体船在静水中的总阻力时历曲线。取40~70 s时间段内的结果曲线进行比较，如图4（图中，R_t为总阻力）所示。由图中曲线可见，3种时间步长的总阻力计算平均值与试验值的偏差约为6%，Δt=0.02 s时的总阻力时历曲线与Δt = 0.01 s时的非常接近。

Figure 4.  Time history curves of total resistance for different time steps

取40 ~70 s时间段内的总阻力计算平均值作为偏差衡量标准，以Δt = 0.01 s的计算结果为基准，比较其他时间步长计算结果与其之间的偏差，如表2所示。从表中可看出，Δt = 0.02 s的总阻力计算平均值与Δt = 0.01 s的计算平均值仅相差−0.1%，满足时间步无关性要求。本文数值模拟计算采用的时间步长为Δt = 0.02 s。

通过改变自动网格模块下基础尺度的值来调整网格数量。研究了3种不同网格数下的CFD计算模型，3种计算模型的网格数分别为150万、240万和392万。模拟工况为Fr =0.546（V =3.416 m/s），Δt = 0.02 s。

图5所示为3种网格模型的总阻力时历曲线比较。由图可见，网格数为240万和392万的计算结果曲线非常接近。

Figure 5.  Time history curves of total resistance for different grid sizes

取40~70 s时间段内的总阻力计算平均值作为偏差衡量标准，以网格数392万的计算结果为基准，比较其他网格数计算结果与其之间的偏差，结果如表3所示。由表3中可以看出，网格数为240万和392万的计算结果平均值仅相差−0.2%，这表明该网格数对计算结果几乎无影响，满足网格收敛性要求。本文的数值模拟采用网格数为240万的计算模型。

根据以上验证结果，选择时间步长Δt = 0.02 s、网格数量为240万的数值模型进行计算。数值计算范围为：Fr = 0.1~1.0（${{F}}{{{r}}_\nabla } = 0.28\sim 2.80$），共选取21个速度点。

计算前三体船约束模型的静水总阻力，并将其与自由状态下的静水阻力试验值^[18]进行比较，如图6所示。

Figure 6.  Comparison between total resistance and experimental values of constrained model

由图可见，Fr < 0.3时，约束模型的计算总阻力与自由状态下的试验值非常接近；而在Fr > 0.3的中、高速段，两者相差较大，其中在Fr = 0.5时试验值比约束模型计算值高出了34.4%。

在自由状态下的静水阻力试验中，测量该前三体船模艏、艉测量点的升沉值^[16]。模型试验中艏、艉测量点的位置如图7所示。其中，升沉值沿垂直方向上升为正，下降为负；纵倾值指船体纵向绕质心旋转的角度，艏倾为正，艉倾为负。

Figure 7.  The position diagram of sinkage measurement points at bow and stern

对前三体船自由模型的升沉和纵倾进行数值计算，监测三体船的升沉和纵倾随Fr的变化，并与自由状态下的模型试验值进行比较，结果分别如图8和图9所示。

Figure 8.  The sinkage curves of measurement points at bow and stern of trimaran model

Figure 9.  The trim curves of trimaran model

由图8可见，自由模型艏、艉测量点的升沉计算值与试验数据吻合良好。船体艏部测量点在Fr > 0.4时出现了明显的抬升，而艉部测量点则是从Fr = 0.3开始出现明显的下沉。随着航速的不断增加，艏、艉测量点的升沉值也随之增大，至Fr > 0.6时艏、艉测量点的升沉值均趋于一个比较稳定的值，随航速的变化开始变缓。

由图9可见，自由模型的纵倾计算值与试验数据吻合较好。在整个航速区间内，三体船的纵倾角均为负值，即船体在不同航速下航行时，始终处于艉倾状态。当Fr < 0.3时，三体船的纵倾并不明显；当0.3 < Fr < 0.6时，纵倾角随着航速的增加快速增大；在Fr = 0.6时，纵倾角达到最大值，约1.44°；当0.6 < Fr < 0.8时，纵倾角几乎保持不变，其大小稳定在1.4°左右；而当Fr > 0.8时，纵倾角则随航速的增加缓慢减小，在Fr = 1.0时，其值约为1.25°。

由以上数值计算和模型试验结果可以看出，在Fr > 0.3的中、低速范围内，三体船的航行姿态变化并不明显，而在Fr > 0.3的中、高速段，三体船出现了非常显著的航行姿态变化。

将数值计算得到的自由模型静水总阻力与自由状态下的模型试验总阻力进行对比，如图10所示。由图可见，当Fr < 0.3时，自由模型的计算总阻力与试验值非常接近，这与约束模型的计算结果一致，说明当Fr < 0.3时，航行姿态变化较小，对总阻力的影响甚微，是否考虑航行姿态的影响几乎不影响总阻力计算结果。而当Fr > 0.3时，约束模型阻力计算值与试验值的偏差幅度增大，这说明考虑了航行姿态变化的自由模型计算值精度显著提高。当Fr = 0.5时，约束模型的总阻力计算值与试验值相差34.4%，而自由模型阻力计算值与试验值的偏差值则降低了7.2%，即计算中考虑航行姿态变化后计算结果精度显著提高。

Figure 10.  Comparison of total resistance and experimental values of free model

自由模型的摩擦阻力R_f和压差阻力R_p随航速变化的曲线如图11所示。

Figure 11.  The frictional resistance and pressure resistance curves of free model

由图11可见，在Fr > 0.65的高速段，压差阻力并没有随着航速的增大而迅速增加，反而呈现出非常平缓的增长趋势。这主要与船体自身的型线特征有关，说明该船型在高航速下具备优良的兴波阻力性能。

为了定量分析摩擦阻力和压差阻力在自由模型总阻力中所占的比例，绘制了两种阻力成分占总阻力的比例曲线，如图12所示。由图可见，当0.3 < Fr < 0.65时，压差阻力明显大于摩擦阻力；而在Fr > 0.65的高速段，摩擦阻力成为总阻力的主要成分，并且其比例随航速的增加而迅速增大，其中在Fr = 1.0时，摩擦阻力是压差阻力的2倍。

Figure 12.  The percentage curves of frictional resistance and pressure resistance of free model

为便于比较，定义了以下参数：δR_t，δR_f，δR_p，ε_t，ε_f，ε_p。其中：δR_t=R_t自由–R_t约束，δR_f=R_f自由–R_f约束，δR_p=R_p自由–R_p约束，δR_t，δR_f，δR_p分别为自由模型的总阻力、摩擦阻力和压差阻力变化量，R_t自由，R_t约束，R_f自由，R_f约束，R_p自由，R_p约束分别为自由模型和约束模型的总阻力、摩擦阻力和压差阻力；ε_t =δR_t/R_t约束，ε_f =δR_f/R_f约束，ε_p =δR_p/R_p约束，ε_t，ε_f，ε_p分别为自由模型的总阻力变化率、摩擦阻力变化率和压差阻力变化率。

自由模型与约束模型的数值计算总阻力曲线如图13所示。由图可见，2种模型的总阻力曲线在0.3 < Fr < 0.88时差别较明显，且自由模型的总阻力明显高于约束模型；在Fr > 0.88的甚高速段，自由模型的总阻力略低于约束模型。

Figure 13.  The total resistance curves of free model and constrained model

为定量分析不同航速下航行姿态对总阻力的影响幅度，绘制了自由模型的总阻力变化率ε_t曲线，如图14所示。

Figure 14.  The curve of total resistance change rate of free model

由图14可见，在整个航速范围内，自由模型的总阻力变化率ε_t变化显著。其中，0.3 < Fr < 0.85时的总阻力变化率ε_t明显高于其他航速段，并在Fr = 0.5，0.7时出现了2个峰值，总阻力变化率分别为24.8%和14.4%；在Fr>0.88的高速段，自由模型的总阻力反而略低于约束模型，总阻力变化率约为−2%。

由此可见，在半滑行高速段，航行过程中的航行姿态变化对三体船总阻力的影响非常显著。

自由模型和约束模型的摩擦阻力比较曲线如图15所示。由图可见，当Fr < 0.55时，两者相差较很小；而当Fr > 0.8时，两者的差距较为显著。

Figure 15.  The frictional resistance curves of free model and constrained model

自由模型的摩擦阻力变化率ε_f曲线如图16所示。由图可见，当0.3 < Fr < 0.5时，自由模型的摩擦阻力变化率约为2%~3%；在Fr > 0.55的高速段，自由模型的摩擦阻力变化率为负值，其中在Fr = 0.9时达到最低值，约为−8.4%。

Figure 16.  The curve of frictional resistance change rate of free model

与约束模型相比，自由模型的摩擦阻力在Fr < 0.55时是增加的，但在Fr > 0.55时不增反降。

为分析主体和双侧体对三体船摩擦阻力改变量的影响效果，绘制了主体、双侧体以及三体船的摩擦阻力改变量曲线，如图17所示。

Figure 17.  The curves of frictional resistance change of free model

由图17可见，当Fr > 0.4时，双侧体的摩擦阻力随Fr的增大而不断减小。这主要是由于自Fr = 0.4之后，三体船出现了非常显著的升沉和纵倾现象，双侧体的部分船体抬出水面，从而导致湿表面积不断减小，进而引起摩擦阻力的减小。主体的摩擦阻力增加量在Fr < 0.85时始终为正，而在Fr > 0.85时为负。这是因为主体艏部尖瘦而艉部宽肥，纵倾角较小时的艉倾现象能够带来湿表面积的增大，从而导致摩擦阻力增加；而当速度更高时，方艉船型会使得船体艉部小幅抬升，艉部横剖面不再接触水面，并出现“虚长度”，此时湿表面积减小，从而导致摩擦阻力减小。图18所示的图片中也显示出了这一现象。

Figure 18.  The experimental picture and computational contours at high speed（Fr=0.9）

自由模型和约束模型的压差阻力曲线比较如图19所示。由图可见，自由模型和约束模型的压差阻力曲线在Fr > 0.3时差别非常显著。相应的压差阻力变化率ε_p曲线如图20所示。

Figure 19.  The pressure resistance curves of free model and constrained model

Figure 20.  The curve of pressure resistance change rate of free model

由图20可见，在整个航速区间内，自由模型的压差阻力变化率均为正值。当Fr > 0.3时，压差阻力变化率均大于10%，并出现了2个波峰，分别对应Fr = 0.5和0.7，变化率分别为48.7%和40.7%。对比自由模型的纵倾曲线可以看到，在整个航速范围内，纵倾值始终为负，即船舶始终处于艉倾状态航行。由于该船型在船舯之后的型线偏肥大，因此这样的航行姿态会导致更大的兴波阻力和黏压阻力。

为了分析自由模型状态下不同阻力成分对总阻力变化的贡献，绘制了如图21所示的曲线，纵坐标为摩擦阻力和压差阻力的变化量。

Figure 21.  The change of friction resistance and pressure resistance of free model

由图21可见，在整个航速区间内，压差阻力的变化量始终为正；摩擦阻力的变化量在Fr < 0.55时为正，在Fr > 0.55时为负；压差阻力的变化量曲线明显高于摩擦阻力。

结合总阻力变化率分析（图14）可知，在Fr < 0.88范围内，总阻力的变化量为正值。当Fr=0.5时，压差阻力变化量δR_p =15.7 N，摩擦阻力变化量δR_f =0.5 N，总阻力增量的97%来自于压差阻力。当Fr =0.7时，压差阻力变化量δR_p =16.76 N，而摩擦阻力变化量δR_f = −1.17 N。因航行姿态变化引起的总阻力增量90%以上来自于压差阻力的变化。

本文选用标准$k - \varepsilon $湍流模型，采用静水VOF法，使用重叠网格技术，针对某前三体船型的自由模型和约束模型进行了静水阻力计算，分析了该三体船型的阻力特性，并重点分析了该三体船航行姿态对阻力性能的影响。主要得到以下结论：

1） 当弗劳德数（Fr < 0.3）较低时，航行姿态对阻力的影响可以忽略不计，这也是排水型船进行阻力预报时不考虑航行姿态影响的原因。一旦进入半滑行阶段，将发生显著的升沉和纵倾现象（纵倾可高达${1.44^ \circ }$），进而阻力也会发生显著变化。

2） 在半滑行状态，航行姿态对摩擦阻力和压差阻力这2个阻力成分的影响不同，其对压差阻力的影响可高达近50%，而对摩擦阻力的影响则不到10%。在半滑行状态，由航行姿态变化引起的阻力增加90%以上来自压差阻力的增加。