基于Verification & Validation的断阶双M船减阻机理分析

余泽爽，毛筱菲,沈小红

(1.中国船级社武汉分社，武汉 430022;2.武汉理工大学a.船海与能源动力工程学院；b.高性能船舶技术教育部重点实验室，武汉 430063;3.中国舰船研究设计中心，武汉 430064)

M船是一种新型的高性能滑行艇，因其艉部剖面形状而得名，双M船与M船类似，由2个滑行面和4个槽道组成，相对于常规滑行艇，M船具备更优的阻力性能、稳定性、隐身性，以及耐波性等，由于特殊的船型结构，在高速航行时，能够产生高效的滑行升力、槽道中的水混合物则能形成水气两相动升力、刚性侧裙则能有效屏蔽主艇体兴波[1-2]。

为进一步改善M船的阻力性能，借鉴常规滑行艇的设计经验，在艇体底部采用断阶结构。对断阶M船(或三体滑行艇、槽道滑行艇等)的研究方法经历了时代的更迭，如基于船模试验研究船型参数对助力和航态的影响以及纵向运动稳定性的影响[3-4]、基于船模试验和实船试航发现了槽道水翼和尾压浪板对航态与助力的有优化作用、基于Savitsky公式通过迭代求解的数值方法求解阻力与航态、基于SIT发和BP神经网格研究了压浪条和断阶对阻力的影响及其机理[5]。但是断阶对于阻力与航态的影响，不同研究得到了不同的结论。有学者通过对有、无断阶以及不同断阶数量的三体滑行艇进行系列试验，认为断阶会导致槽道中流场紊乱，从而对阻力不利，而且断阶越多阻力越大[6-7]。基于CFD技术，有学者通过对双断阶常规滑行艇的绕流场进行分析，认为某种情况下第一断阶会导致第二断阶处产生过大浸湿面积，是阻力恶化的原因，因此需要避开该情况[8]，同时讨论了断阶对三体滑行艇阻力和航态性能的影响，认为断阶能够减少主体滑行面的浸湿面积、提升升力，从而起到减阻的效果，此外还提出了断阶在艏部水线面的驻点处时，会导致有害喷溅，从而影响阻力性能[9]。有学者基于计算流体力学方法(CFD)和人工神经网络预测方法(ANN)对某滑行艇在具有无、单、双断阶时，不同重心纵向位置和断阶间距下的助力性能进行了研究，发现断阶使低航速时助力增加、高航速时助力降低，并且增加双断阶间距能降低阻力和提高纵向运动稳定性，即有效避免海豚运动[10]。

综上，断阶对于双M船阻力性能的影响，是一个比较复杂的问题。为分析断阶对双M船的减阻机理，选用一标准数值方程的双M船，在不同位置设置断阶，并保持相同排水量和重心位置，以保证计算结果可比性；同时为对数值方法进行论证，本研究首先对数值方法(主要针对网格)进行了V&V验证，并通过多船型、多工况的数值计算与试验验证，最终证明该网格方案的可靠性。

1 对象

选用的母型船来自文献[11]中提出的M船横剖线的数值控制方程。图1所示为分别在距船艉0.5L和0.25L处设置断阶后形成的四艘船分别记为DM1、DM2、DM3及DM4，主要参数见表1。

图1 不同断阶结构的双M船

表1 主要参数

2 数值方法与V&V

2.1 数值方法

基于黏性CFD方法，通过求解RANS方程来实现对N-S方程(式(1)为连续性方程，式(2)为动量守恒方程)的解析，湍流模型选用Rearlizablek-ε，采用VOF法和QUICK离散格式捕捉自由液面和水气两相，使用重叠网格技术和DFBI(Dynamic Fluid Body-Interaction,动态流-固交互作用)模块模拟刚体的六自由度运动。

(1)

(2)

式中：p为压力；ui为速度分量；xi为位置分量；v为速度矢量；F质量力矢量；ν为运动黏性系数。

采用的数值算法为SIMPLE(semi-implicit method of pressure-linked equation，压力耦合方程的半隐式方法)算法，SIMPLE算法求解过程见图2。计算域与边界条件见图3。

图2 SIMPLE算法求解简图

图3 计算域与边界条件

2.2 数值计算

为了对数值方法(主要针对网格)进行收敛性及不确定度分析，选取某双断阶双M船(记为M-1)作为计算对象，通过改变船体表面网格尺寸，并保持网格投影、边界层过渡系数等其他参数不变，来控制网格密度。4种不同密度的网格方案见表2，断阶处网格示意于图4。

表2 不同网格方案

图4 断阶处网格示意

采用上述4套网格方案对M-1从Fr▽=0.820到Fr▽=4.935的静水直航2-DOF运动进行数值模拟，阻力的无量纲系数(阻升比R/Δg)与航态(纵倾角φ)的计算结果以及试验结果见图5、6。

图5 阻力计算结果与试验结果对比

图6 纵倾角计算结果与试验结果对比

比较可见，纵倾角的计算结果受网格数量的影响并不大，这体现了重叠网格技术在捕捉动航态问题时的优势，因此这里主要对阻力进行分析。

为了保证数值计算的可靠性，采用ITTC的V&V法进行网格独立性分析、不确定度分析以及试验验证，从标准化的规程中，确定可信度在95%的数值计算方案，为后续研究提供理论基础。

2.3 Verification

数值计算误差δS为数值计算结果S与真实值T之差，包括建模误差δSM与数值误差δSN，这里主要考虑后者，数值误差及其不确定度的计算如式(3)和(4)所示。

δSN=δI+δG+δT+δP

(3)

(4)

式中:I为迭代;G为网格；T为时间步长；P为其他参数。

本研究中主要针对网格。网格独立性研究是为了判断计算结果是否随网格数量增加而收敛[12]，收敛因子RK的计算结果见图7，可见计算结果均满足单调收敛，故网格收敛性得到满足。

图7 收敛因子计算结果

网格不确定度是否需要修正，取决于网格修正系数CK与1的关系，通过式(5)和(6)进行修正，修正系数计算结果见图8。

图8 修正系数计算结果

可见在Fr▽<2.5时，两者的修正系数均远小于1；在Fr▽>2.5时，Grid2-3-4的修正系数开始接近甚至大于1，因此其网格不确定度需要进行修正。

当CK1或CK1时，

(5)

当CK→1时，

(6)

网格不确定度的计算结果见图9，可见Grid2-3-4的不确定度均比Grid1-2-3大，并且经过修正后的不确定度均得到减小，尤其是Fr▽>2.5时，说明该航行阶段数值计算的复杂性。

图9 网格不确定度计算结果

2.4 Validation

Validation即为验证试验误差δE(试验结果D与计算结果S的误差)是否小于试验验证的不确定度UVAL，如式(7)和(8)所示，其中读数误差δSPD与模型误差δSMA难以确定，故采用保守计算，即如式(9)所示。

(7)

(8)

(9)

式中：δD为试验误差，UD为试验不确定度，按试验结果的2.5%给出。

试验验证不确定度与试验验证误差的计算结果见图10，可见，在Fr▽<2.5时，均通过了试验验证，在2.54时，由于不确定度过小，所以均未能通过试验验证。

图10 试验验证计算结果

2.5 小结

1)在Fr▽<2.5时，Grid2和Grid3的计算结果通过了V&V验证，选择Grid2的方案较为合适。

2)在2.5

3)在Fr▽>4时，尽管Grid3的计算结果与试验误差均在10%以内，但未能通过Validation，故在Grid2和Grid3选择中间值，既能保证较小的不确定度，又能得到更精确的结果。

3 数值计算与结果分析

3.1 阻力与纵倾角的比较

基于上述V&V的数值计算与验证过程，采用推荐的网格方案，对4艘双M船静水直航的2-DOF运动进行数值模拟。阻力R与纵倾角φ的计算结果见图11、12，其中纵倾角φ>0为尾倾，计算过程中，保证体积傅汝德数Fr▽一致。

图11 阻力计算结果对比

图12 纵倾角计算结果对比

计算结果显示，DM3在Fr▽=5.260时发生了海豚现象，其阻力与纵倾角均取自时均值，其他船型均未发生失稳现象。

由上述比较可见，在排水和过渡阶段(Fr▽<3)时，阻力和纵倾角：DM4>DM2>DM3>DM1；DM在进入滑行阶段 (Fr▽→3)时，DM2的阻力开始小于DM1；在Fr▽≈3.5时，DM45)时，阻力与纵倾角与排水和过渡阶段完全相反：DM4

3.2 艇底压力分布

DM1和DM3艇底压力载荷中心位置随航速的变化曲线见图13，其中载荷中心以船长L的比例给出，其含义为船体所受载荷的合力位置，代表动负载荷作用于滑行面的情况。

图13 艇底载荷中心位置

Fr▽=3.682和4.734时，艇底压力分布情况见图14、15，结合图13，可见在排水和过渡航行阶段时，艇底载荷中心位置随航速增加而前移；在滑行阶段时，该位置开始后移，当其刚好处于或接近断阶时，艇底压力被重新分配，因此阻力得以改善，如DM2

图14 艇底压力分布(Fr▽=3.682)

图15 艇底压力分布(Fr▽=4.734)

同时对比DM1和DM3在Fr▽=5.260时的艇底载荷中心位置，说明海豚运动的发生与载荷位置相关，而DM3由于断阶的缘故，艉部的滑行面较宽，故而其载荷位置较为靠后。

3.3 艇底水气分布

DM1和DM2在航行过程中，水线面与艇底相交处的位置(即驻点)随航速的变化曲线见图16，其位置表征了艇体的航态，与纵倾和升沉模态有关。

图16 水线面与艇底的驻点位置

Fr▽=3.682和4.734时，艇底水气分布情况见图17、18，可见在Fr▽≈1.5时，驻点位置出现一个极小值，这是由于航态使然；在滑行阶段时，驻点位置逐渐后移，当其刚好位于断阶前不远处时(如图17的DM1与DM3)，会产生反向的不利喷溅，反而使阻力恶化；但另一方面，由于断阶使DM3滑行面与槽道均有所增宽，从而提高了其滑行效率与槽道中的水气通畅，故而又使阻力得以改善，综合作用下，DM3与DM1的阻力相差不大，且变化趋势一致。

图17 艇底水气分布(Fr▽=3.682)

图18 艇底水气分布(Fr▽=4.734)

同时对比DM1与DM3以及DM2与DM4可知，断阶形成的有利影响是使断阶后的滑行面上形成水气混合现象，并增加了槽道宽度；而不利影响则是当断阶与驻点处于不合适的相对位置时，易形成不利喷溅，从而增加阻力。

4 结论

1)在Fr▽<2时，Grid2和Grid3均通过V&V验证，选择Grid2较为合适；在24时，由于Grid3网格不确定度过小。因此，应在Grid2和Grid3之间进行网格调整，以保证不确定度，并提高计算精度。

2)在排水和过渡航行阶段(Fr▽<3)时，阻力与纵倾角均满足：无断阶(DM1)<单断阶(DM35)时，阻力与纵倾角有：无断阶(DM1)>单断阶(DM3双断阶(DM4)，这是由于断阶改善了艇底的压力与水气分布情况。

3)在滑行航行阶段，当断阶(比较靠后时)位于艇底载荷中心附近时，艇底压力重新分配，引起滑行面上载荷情况得以改善。

4)在滑行行阶段，当断阶(比较靠前时)位于驻点后不远处时，断阶会引起反向的不利喷溅，从而使阻力恶化；当断阶远离驻点时，槽道中的空气会由于断阶处的负压，而被吸入断阶，并进入滑行面，从而改善滑行面上的水气分布情况。

5)在超滑行航行阶段，其阻力与纵倾角只与艇体后半部的形状相关，由于断阶引起后半部的槽道宽度得以增加，因此对阻力和纵倾角均有利。