0 引　言

当船舶高速航行时，由于兴波阻力与速度的高次方成正比，兴波阻力急剧增加，所以传统排水型船舶速度的提高受到制约。水翼艇作为高性能船型之一，航速足够高时，由水翼产生的升力会将船体托升到水面以上，从而减小艇体湿表面积，达到减小船舶总阻力的目的，在客运、军事、救援、科考等方面具有广泛的应用前景。早在1861年，英国就开始尝试建造水翼艇；直到1906年，世界上第一艘水翼艇才由意大利人恩里科·弗兰尼尼试验成功，该艇在试验期间，最高航速达到了37 kn。中国建造的第一艘水翼艇是一艘40座的水翼客艇——“水翼−”I号，由中国船舶科学研究中心设计，并于1959年下水在长江上开始运营^[1]。

水翼艇的阻力性能预报是备受关注的问题，国内外学者在这方面也进行了相关研究。早期高速艇阻力性能预报的方法有：1）半理论—半经验方法；2）根据系列试验资料整理得到的图谱方法；3）模型试验方法^[2]。杨松林等^[3]针对翼滑艇（滑行艇首部加装水翼）开展了模型阻力试验，通过对2个不同水翼安装角度的浅V滑行面翼滑艇进行试验，发现与主尺度相近的滑行艇相比，翼滑艇在同等航速下阻力降低约5% ~ 15%。唐建飞等^[4]对滑行艇、单槽道滑行艇及排水型高速船3种不同类型的高速艇加装艏部水翼的方案进行了试验研究，试验表明加装水翼能够在高航速下减少30%以上的阻力。周俊麟等采用模型试验方法研究了水翼翼型、安装高度与攻角等参数对水翼艇阻力性能的影响，为水翼艇研究设计提供了参考。

近年来，随着计算机水平的发展及CFD方法预报精度的不断提高，目前CFD方法已广泛应用于船舶的水动力性能研究。Brizzolara等^[5]利用RANSE求解器进行了滑行艇阻力的数值预报，数值结果与试验结果及半经验公式结果的对比表明，CFD方法对于总阻力的预报精度约为10%，精度要高于Savitsky等半经验公式。孙华伟^[6]利用商业软件STAR-CCM+研究了网格方案、湍流模型、时间步长对滑行艇静水阻力数值预报的影响。De Marco等^[7]采用重叠网格和变形网格方法进行了滑行艇水动力分析，认为在进行高速滑行艇水动力分析时采用重叠网格方法的结果优于变形网格方法。

综上所述，现有研究主要聚焦于特定航速或单一水翼攻角下的水翼艇阻力和航态分析，而对水翼攻角在不同航行状态下的影响规律研究较少，而选择合适的攻角对于降低阻力，提升水翼艇经济性非常重要。本文基于CFD方法，以深浸式T型水翼的双水翼艇为对象，采用重叠网格技术对水翼艇在中速过渡阶段以及高速翼航阶段进行数值模拟，研究重点分析艏部水翼攻角对不同航行状态下船舶阻力性能和航态的影响，并揭示其变化规律。

1 数值方法 1.1 控制方程

基于雷诺平均纳维-斯托克斯方程（RANS方程）数值求解的方法是目前使用最为广泛的湍流数值模拟方法。本文采用RANS方程数值求解的方法进行水翼艇的水动力性能研究。在笛卡尔坐标系下，RANS方程的表达式为：

$ \frac{\partial\overline{u_i}}{\partial x_i}=0，$

(1)

$ \rho\frac{\partial u_i}{\partial t}+\rho\frac{\partial\overline{u_i}\overline{u_j}}{\partial x_j}=-\frac{\partial p}{\partial x_i}+\mu\frac{\partial}{\partial x_j}\left(\frac{\partial\overline{u_i}}{\partial x_j}-\rho\overline{u_i^{'}u_j^{'}}\right)。$

(2)

式中：$ \rho $为密度；$ p $为流体微元体上的压力；$ t $为时间；$ \mu $为动力黏度；$ \overline {{u_i}} $和$ \overline {{u_j}} $均为速度分量时均值；$ {x_i} $和$ {x_j} $分别为i方向和j方向的坐标，i、j=1，2，3；$ \mu $为流体动力粘度；$ - \rho \overline {u_i^{'} u_j^{'} } $为雷诺应力项。

1.2 湍流模型

由于在RANS方程中出现了雷诺应力项，导致原来的控制方程不再封闭，所以需要建立新的方程对RANS方程进行封闭求解。根据对雷诺应力的假定或处理方式的不同，常用的2种湍流模型为雷诺应力模型和涡粘模型，其中涡粘模型通过引入涡粘系数来建立雷诺应力与平均速度梯度的关系，而非直接处理雷诺应力输运方程，因此在实际工程中应用更为广泛，其控制方程为：

$ - \rho \overline {u_i^\prime u_j^\prime } = {\mu _t}\left( {\frac{{\partial {u_i}}}{{\partial {x_j}}}} \right) - \frac{2}{3}\left( {\rho k + {\mu _t}\frac{{\partial {u_i}}}{{\partial {x_j}}}} \right){\delta _{ij}} 。$

(3)

式中：$ k $为湍动能；$ {u_i} $为时均速度；$ {\mu _t} $为湍流粘度；$ {\delta _{ij}} $为克罗内克尔（Kronecker）数。计算选用基于式（3）的Realizable k-ε湍流模型封闭控制方程，对于不可压缩粘性流体，Realizable k-ε模型中关于k和ε的输运方程为：

$ \rho\left[\frac{\partial k}{\partial t}+\frac{\partial\left(ku_i\right)}{\partial x_i}\right]=\frac{\partial}{\partial x_j}\left[\left(\mu+\frac{\mu_t}{\sigma_k}\right)\frac{\partial k}{\partial x_j}\right]+G_k-\rho\varepsilon，$

(4)

$ {\rho \left[ {\dfrac{{\partial \varepsilon }}{{\partial t}} + \dfrac{{\partial \left( {\varepsilon {u_i}} \right)}}{{\partial {x_i}}}} \right] = \dfrac{\partial }{{\partial {x_j}}}\left[ {\left( {\mu + \dfrac{{{\mu _t}}}{{{\sigma _\varepsilon }}}} \right)\dfrac{{\partial \varepsilon }}{{\partial {x_j}}}} \right] + \rho {C_1}E\varepsilon - \rho {C_2}\dfrac{{{\varepsilon ^2}}}{{k + \sqrt {v\varepsilon } }}} 。$

(5)

式中：各参数的具体含义见文献[8 - 11]。该模型能够处理时均应变率较大的流动问题，对于有旋度、带分离和应变率相对较大的湍流问题计算精度较高，故本文的数值模拟采用Realizable k-ε模型。

1.3 自由液面处理方法

自由面的处理方法主要有界面追踪法和界面捕捉法2种。在处理船舶水动力两相流问题方面，流体体积（Volume of Fluid，VOF）方法具有计算效率高，存储占用小的优点，应用比较广泛。其控制方程为：

$ \frac{\partial\alpha_1}{\partial t}+\nabla\left(\alpha_1\vec{U}\right)=0。$

(6)

式中：$ {\alpha _1} $为气体体积分数；$ \vec U $为流场速度矢量场。

1.4 消波方法

由于计算域的尺度有限，出口边界的波浪反射可能影响数值模拟的精度。为降低反射影响，本文在计算域出口区域设置阻尼层，以实现有效的消波处理。消波区动量方程为：

$ \frac{\partial u}{\partial t} + u\frac{\partial u}{\partial x} + v\frac{\partial u}{\partial y} + w\frac{\partial u}{\partial z} = \mu\left( \frac{\partial^2u}{\partial x^2} + \frac{\partial^2u}{\partial y^2} + \frac{\partial^2u}{\partial z^2} \right) - \frac{\partial p}{\partial x} + C(x)u，$

(7)

$ \begin{split}&\frac{{\partial w}}{{\partial t}} + u\frac{{\partial w}}{{\partial x}} + v\frac{{\partial w}}{{\partial y}} + w\frac{{\partial w}}{{\partial x}} =\\ &g + \mu \left( {\frac{{{\partial ^2}u}}{{\partial {x^2}}} + \frac{{{\partial ^2}u}}{{\partial {y^2}}} + \frac{{{\partial ^2}u}}{{\partial {z^2}}}} \right) - \frac{1}{\rho }\frac{{\partial p}}{{\partial x}} + C(x)w \end{split}。$

(8)

式中：$ C(x) $为阻尼衰减系数；$ \mu $为流体动力黏度；$ g $为重力加速度；$ p $为压强。

$ C(x) = c{\left( {\frac{{x - {x_0}}}{{{x_{\text{{{L}}}}} - {x_0}}}} \right)^2} 。$

(9)

式中：$c$为经验系数；${x_0}$为消波区左边界的x坐标值；$ x\mathit{_{{L}}} $为消波区右边界的x坐标值。

2 数值模型及离散网格

本文基于STAR-CCM+软件，以深浸式T型水翼的双水翼艇为对象，基于RANS方程的数值求解，对不同艏翼攻角下的船舶绕流场进行数值模拟以及船舶阻力、航态的数值预报。

2.1 几何模型

水翼艇的船型参数如表1所示，数值模拟对象为缩尺比2的船模。前后水翼剖面均为NACA6409翼型，展弦比均为8；其中艏部水翼采用三段式翼板外形，上反角为30°；艉部水翼采用一段式翼板；艏部水翼的安装位置在距离船首42%L_wl处，深度为73%D；艉部水翼的安装位置在距离船尾6.24%L_wl处，深度为150%D。水翼艇的三维模型如图1所示。

2.2 计算域及网格划分

在水翼艇绕流场数值模拟中采用重叠网格方法，针对二自由度船舶绕流场的数值模拟建立2个不同的区域，即背景域和重叠域。背景域和重叠域的交界面设置为重叠网格，2个区域通过此交界面进行流场数据的传递。流场背景域入口距离船艏1.5L（L为船长），出口距离船尾4L；船宽方向由于船体沿中纵剖面对称，为了减小计算量，取1/2船体（左舷）进行研究，计算域的侧面距离船体中纵剖面2L；计算域的顶部距离船体1L，底面至船底的距离为2.5L。计算域边界设置如图2所示。

计算域采用切割体网格与棱柱层网格相结合的形式进行网格划分，同时在部分流场参数变化较快的位置（如自由液面、船首船尾处及水翼安装处等）进行网格加密。流体域自由液面处划分情况如图3所示。

2.3 网格收敛性分析

为验证数值模拟方法的有效性，对不同网格密度下的数值模型进行网格收敛性研究。选取F_rσ=3.455下，T型水翼攻角为5°的工况为对象建立3个不同网格密度的计算模型，相邻网格之间基础网格尺寸比r_G=$ \sqrt{2} $，根据网格数量由小到大分别为M₁（低密度网格）、M₂（中密度网格）和M₃（高密度网格），不同网格密度的模型网格数量及总阻力R_t的数值计算结果见表2。

根据公式$ {{{R}}_{{G}}} = \displaystyle\frac{{{{{R}}_{{{tm2}}}}{{ - }}{{{R}}_{{{tm3}}}}}}{{{{{R}}_{{{tm}}1}}{{ - }}{{{R}}_{{{tm}}2}}}} $，求得$ {{{R}}_{{G}}} $=0.353，可以得到收敛因子0<$ {{{R}}_{{G}}} $<1，表明以上3套网格对水翼艇总阻力满足单调收敛。此外正可以看出，高密度网格模型总阻力数值结果相对于中密度网格模型变化相对有限，为了平衡计算精度以及计算效率，后续计算将采用中等密度网格进行。

3 不同攻角下水翼艇数值模拟结果

基于第2节中所述的数值方法，采取中密度网格方案对水翼艇在F_rσ=2.488（过渡状态）及F_rσ=3.455（翼航状态）下，无水翼工况与艏部水翼攻角在1°～10°范围内10种攻角变化的水翼艇绕流场进行数值模拟，计算采用的水翼艇的几何参数见表1。不同攻角下水翼艇过渡状态及翼航状态的阻力及航态变化如表3、表4及图4～图6所示，其中纵倾角艏倾为正，升沉值上浮为正。

结合图表可以看出：

1）在过渡状态下，水翼艇的总阻力随艏部水翼攻角的增大呈现先降低后上升的趋势，其中攻角5°时总阻力最低，相较于无水翼工况，该攻角下的减阻效果达到32.31%；翼航状态下，攻角4°时阻力性能最优，相比于无水翼的工况减阻效果能够达到35.38%；另外从阻力成分的曲线对比中可见：过渡状态及翼航状态下，在过渡状态及翼航状态下，水翼艇的总阻力主要由剩余阻力（包括兴波阻力和诱导阻力）构成，其占比高于摩擦阻力。

2）在过渡状态下，水翼艇整体表现出艉倾趋势，且艉倾角随着攻角的增大而增加，在攻角10° 工况下，艉倾角达到3.021°；在翼航状态下，随攻角增大水翼艇的航态由艏倾转变为艉倾，在攻角为10°的工况下船舶艉倾值为1.304°，在攻角为5°的工况下，艇体基本处于正浮状态。

3）在过渡状态时尽管未到达翼航航速，但得益于水翼的升力，艇体相较于无水翼状态仍有明显的抬升，艇体湿表面积逐渐减少，而在水翼攻角逐渐增大超过最佳攻角的过程中，艉倾加剧引起剩余阻力增加，造成了总阻力的增大；而当水翼艇进入翼航状态时，艇体被完全抬升出水面，与水面几乎没有接触，湿表面积大幅减小，导致摩擦阻力减少，此时由于水翼攻角变化导致的水翼阻力变化是水翼艇总阻力变化的主要因素。

升阻比是衡量水翼水动力性能的关键指标之一，较高的升阻比通常意味着更高效的升力与较低的阻力。不同攻角下艏部水翼的升阻比变化如图6所示。

升阻比计算式为：

$ K = \frac{{{C_L}}}{{{C_D}}} = \frac{{{L_{{\text{foil}}}}/\dfrac{{\text{1}}}{{\text{2}}}\rho {V^{\text{2}}}{S_{{\text{foil}}}}}}{{{D_{{\text{foil}}}}/\dfrac{{\text{1}}}{{\text{2}}}\rho {V^2}{{\text{S}}_{{\text{foil}}}}}} 。$

(10)

式中：$\rho $为流体密度；$ {{{C}_{{L}}}} $、${{\text{C}}_{{D}}}$分别为水翼的升力系数及阻力系数；$ \mathit{{L}_{\text{foil}}} $、$ \mathit{{D}_{\text{foil}}} $分别表示水翼的升力及阻力；V为航速；${{{S}}_{{\text{foil}}}}$为水翼的投影面积。

可以看出：当攻角小于4°时，随着攻角的增大，升阻比呈现上升的趋势；在攻角为4°的工况时，过渡状态及翼航状态的升阻比均出现最大值；而随着攻角的继续增大，升阻比开始减小。

图7和图8给出了水翼艇在过渡状态和翼航状态下无水翼以及水翼攻角3°、5°、8°的自由液面波形图。可以看出，当Fr_σ=2.488时，水翼艇前部两侧的船体与水接触较多，导致明显的兴波效应，其产生的船行波沿船体两侧扩散，并在船首区域形成典型的V型波；随着船速增加至Fr_σ=3.455，流体从艇体下方流过后在靠近船尾区域形成“空穴”，艇后“空穴”的长度随着航速的增加而增长。相比于有水翼的工况，无水翼的工况会在艇后出现明显的鸡尾状“水丘”，“水丘”的高度随着航速的增加而降低，且最大波峰向后移动。

压力系数$ C_p$定义为：

${C_p} = \frac{{P - \rho gh}}{{0.5\rho {V^2}}} 。$

(11)

式中：P为监测点的压力值。

图9和图10给出了水翼艇在过渡状态和翼航状态下艇底压力分布。可以看出，在没有水翼的工况下，艇体直接受到来流冲击，艇首较强的正压力分布，压力驻点随着航速的增加逐渐向前移动，在有水翼的情况下，压力主要集中于水翼下表面前缘，并沿翼型向后缘逐渐降低。水翼底部存在高压区，在同一航速下，随着攻角的增大，高压区的范围也增大，导致水翼上下两侧的压力差增大，升力增大。攻角较小时，艇底压力整体较均匀，表明此攻角下水翼对艇体的升力贡献开始显现；攻角为5°时，水翼升力进一步增强，但艇底压力分布相对均匀，有助于减少总阻力；攻角增加到8°时，艇底前部压力降低，表明水翼升力继续增加，使艇首尾部压力增大，可能导致艉倾加剧，使得阻力增加。当航速处在过渡状态时，艇底压力仍然较高，表明艇体与水面的接触较多，摩擦阻力和兴波阻力仍然较大；随着航速增加，水翼升力作用加强，将艇体抬升出水面而不与水面接触，翼航状态时艇底与空气接触多为低压分布，此时阻力降低。

4 结　语

本文通过STAR-CCM+软件，采用重叠网格技术进行了水翼艇不同艏部水翼攻角下2种航行状态的阻力及航态预报，并对网格的收敛性进行了分析。可以得到以下结论：

1） 当水翼艇达到过渡状态和翼航状态需要的航速时，水翼能够减小航行阻力。减阻效果最好的工况分别是5°攻角（过渡状态）和4°攻角（翼航状态）的工况，总体的减阻效果在20.64% ~ 35.38%。合适的攻角下水翼的存在均能起到一定的降低纵倾、改善航态的作用。但当攻角过大时，虽然升力进一步增强，但造成艇体艉倾现象加剧，导致艇尾与水体接触增加，使兴波阻力上升，不利于航行效率及安全。

2） 随着艏部水翼攻角的增大，水翼的升阻比呈现先增大后减小的趋势，升阻比的最大值出现在攻角为4°的工况；而水翼产生的升力不断增大，在其他条件不变的情况下，艏部升力增大会增加船舶的纵倾力矩，导致船舶尾倾加剧。

3） 在本文研究的工况下，过渡状态下艏部水翼的最优攻角为4°，翼航状态下的最优攻角为5°，该攻角条件下，船舶表现出最低阻力和最佳航态平衡。

4） 本研究系统地探讨了攻角对水翼艇阻力及航态的影响，并结合流场分析解释了阻力变化的机理，这一研究结果可为水翼艇的设计与操纵提供理论依据。但本研究仍有拓展的空间，未来可以开展水翼智能控制等方面的应用，使得水翼可以根据自身航态、环境载荷等能对攻角进行实时调整，以便水翼艇能够在航行时能动态地降低阻力，提高稳定性。