0 引　言

作为现代海军作战力量的重要组成部分，潜艇在空间立体作战及战略威慑中均发挥着无可替代的作用^[1]。当前按照动力来源，潜艇可分为核动力及常规动力两大类，前者多在深海中隐蔽航行，其续航力及水下潜航时间均远优于常规动力潜艇，然而对于近海岸侦察、打击及巡逻等军事任务，各国则更多地应用排水量相对较小的常规动力潜艇。受于作战环境及战术限制，近岸航行的常规动力潜艇往往需要完成近水面换气发电、救援、发射导弹及潜望镜深度瞭望等战术动作，较于深水潜航状态，近水面航行的潜艇受到文丘里效应的影响，潜艇水动力随潜深、航速均呈现较为明显的变化。

在早期潜艇操纵过程中，国内外学者已认识到潜艇近水面航行时水动力性能的变化^[2]，并认为当潜深 $\zeta \text{≤} 3D$ （ $D$ 为潜艇直径）时，艇体受到较为明显的吸力和力矩，且吸力随下潜深度的增加呈指数规律衰减；随着潜艇设计能力与技术的不断提升，在较长的时间内，很多学者采用切片理论及势流方法来求解潜艇近水面状态下的水动力性能^[3-4]，目前仍有较多相关研究单位基于Havelock源^[5]、Frank源^[6]等的势流方法及Ship-flow等势流软件^[7]开展考虑波浪影响等复杂工况下的潜艇近水面水动力特性研究，但上述方法忽略了流体粘性的影响，适用范围往往受到工况条件与计算经验的限制。

近年来，随着粘流CFD技术的快速发展，基于Fluent，Star-CCM+及CFX等商用软件平台可快速获取不同航态下潜艇的水动力性能，极大地提高了工程效率和研发周期，张楠等^[8-9]较早地应用粘流计算软件开展了近水面航行潜艇数值计算，其重点分析了潜艇周围流场特性及近水面对螺旋桨伴流的影响，并未对潜深及航速等因素展开分析；ZHANG^[10]通过试验及数值方法分析了AUV近水面航行时下潜困难的现象，分析了影响AUV下潜的主要因素；Wilson H S^[11]则应用CFX软件以Suboff标模裸艇体为研究对象，开展了不同潜深及航速下潜艇水动力性能分析，但其选用模型长仅为1.556 m，其结论对于带有附体及更大模型尺度潜艇是否适用还有待进一步研究。

1 数值方法 1.1 控制方程及湍流模型

文中以RANS方程作为求解潜艇近水面航行的基本方程，其具体形式如下：

$ \begin{split} \frac{\partial }{{\partial t}}\left( {\rho {u_i}} \right) + &\frac{\partial }{{\partial {x_j}}}\left( {\rho {u_i}{u_j}} \right) = - \frac{{\partial p}}{{\partial {x_i}}} + \rho {f_i} + \\ &{\frac{\partial }{{\partial {x_j}}}\left[ {{\mu _0}\left( {\frac{{\partial {u_i}}}{{\partial {x_j}}} + \frac{{\partial {u_j}}}{{\partial {x_i}}}} \right) - \frac{2}{3}{\mu _0}\frac{{\partial {u_l}}}{{\partial {x_l}}}{\delta _{ij}}} \right] + \frac{\partial }{{\partial {x_j}}}\left( { - \rho \overline {{u_i^\prime}u_j^\prime } } \right)}\text{，} \end{split} $

(1)

式中： $\rho $ 为流体密度； $\mu $ 为流体粘度； $p$ 为静压； $f_i$ 为单位质量的质量力； ${{u_i}}$ 和 ${{u_j}}$ 为速度分量。

湍流模式选取SST $k - \omega $ 模型，详细推导过程参见文献[12]，下面仅给出该模型的数学表达式：

$ {\frac{\partial }{{\partial t}}(\rho k) + \frac{\partial }{{\partial {x_i}}}\left( {\rho k{u_i}} \right) = \frac{\partial }{{\partial {x_j}}}\left( {{\Gamma _k}\frac{{\partial k}}{{\partial {x_j}}}} \right) + {G_k} - {Y_k} + {S_k}}\text{，} $

(2)

$ {\frac{\partial }{{\partial t}}(\rho \omega ) + \frac{\partial }{{\partial {x_i}}}\left( {\rho \omega {u_i}} \right) = \frac{\partial }{{\partial {x_j}}}\left( {{\Gamma _\omega }\frac{{\partial \omega }}{{\partial {x_j}}}} \right) + \!{G_\omega } - \!{Y_\omega } +\! {S_\omega }}\text{。} $

(3)

式中： ${{\Gamma _k}}$ 和 ${{\Gamma _\omega }}$ 表示 $k$ 和 $\omega $ 的有效扩散率； ${{G_k}}$ 表示由于平均速度梯度产生的湍流动能； ${{G_\omega }}$ 表示特殊湍流动能耗散 $\omega $ 的产生； ${{Y_k}}$ 和 ${{Y_\omega }}$ 表示由于湍流 $k$ 和 $\omega$ 的耗散； ${{S_k}}$ 和 ${{S_\omega }}$ 为用户自定义源项。

1.2 网格离散方法

Star CCM+商用软件中提供了多种网格划分方式，由于潜艇近水面航行计算的关键在于自由液面兴波的准确捕捉，因此本文中选用六面体切割体网格划分策略对计算域进行网格离散。

切割体网格技术的优势在于能够确保流体域内绝大部分网格为排列规则的六面体网格，同时通过对近艇体壁面区域进行网格加密设置可对潜艇周围的湍流区域进行准确模拟，且对于过于复杂的几何结构，边界层内网格自适应地采用其他形状的网格结构形式，以满足物面形状及网格质量的要求。其生成原理及步骤如下：

步骤1　首先根据艇周围加密区设置及域边界确定网格的生成边界，如图1（a）所示，边界层外边界与域边界圈定的区域为中心网格区，艇体壁面与边界层外边界间的区域为边界层区；

步骤2　在整个计算域及周围生成初始体网格（见图1（a）），用第1步确定的边界进行裁剪，裁剪后的网格如图1（b）所示；

步骤3　根据网格密度设置对中心网格区内网格疏密程度进行控制，并进一步提高网格质量；

步骤4　在边界层区域内，根据物面形状自适应地生成网格。

2 研究对象及计算工况

以美国国防高级研究计划局（DARPA）公布的Suboff全附体模型为研究对象，该模型由主艇体、指挥台围壳及十字型尾舵组成，模型主参数如图2所示。将潜艇中心线距自由液面的距离与艇体直径的比值定义为无因次潜深（见图3），即 ${H^ * } = H/D$ ；无因次航速用艇长傅氏数表示，即 $Fr = U/\sqrt {gL} $ ，式中 $U$ 为潜艇航速， $L$ 为艇长， $g$ 为重力加速度，取值为9.81 m/s⁻²。

本文开展傅氏数 $Fr = $ 0.1～0.5、潜深 ${H^ * } = $ 1.0～2.9D及无限潜深工况条件下潜艇近水面航行水动力数值计算。

3 计算域设置及网格划分

由于Suboff潜艇模型及流场均关于艇体中纵剖面对称，为减小计算量取一半潜艇开展计算，如图4所示。计算域大取为长宽高为5.5L×1.5L×3.0L的方形域，边界条件设置如下：

1）距离艇首1.5倍艇长，边界条件设置为速度入口；

2）距离艇尾4倍船长，边界条件设置为压力出口，压力为未扰动时边界压力；

3）流域的上、下及右边界设置无滑移、不可穿透壁面；

4）流域左边界设置为对称面；

5）艇体表面定义为无滑移、不可穿透边界条件，边界条件设置为无滑移壁面。

在网格离散过程中，为准确地捕捉两相流的水气分界面，对自由液面附近区域进行网格加密。此外，在艇体近壁面设定了局部加密区，确保壁面y+～1，并对指挥台及尾舵处的网格进行了细化处理，船体表面及域剖面网格如图5所示，总网格数为215万。

4 计算方法可行性分析

为验证上述计算方法的可行性，以Suboff模型裸艇体为计算对象，模型长1.556 m，试验数据参见文献[11]。计算域设置及模型无因次潜深的定义与前述模型相同。图6为该模型无因次潜深 ${H^ * } = 1.1$ ，Fr = 0.23～0.64工况下艇体阻力的计算结果与试验值。对比可知，本文方法可较为准确地计算潜艇近自由液面的阻力性能。

由图6中近水面裸艇体阻力随傅氏数的变化可知，随着航速的增大，潜艇近水面航行时艇体阻力呈现先增大后减小再增大的变化趋势，这与深水航行存在较大差异。

5 计算工况与结果分析 5.1 力、力矩定义

为便于计算与分析，设定潜艇重心与浮心纵向位置相同，即距艇首2.0095 m，并定义潜艇近水面航行时所受阻力、垂向力及对重心的纵倾力矩的正方向如图7所示。其中，潜艇垂向力中不包含艇体所受浮力。力及力矩系数定义如下：

$R' = R/(0.5 \times \rho \times {L^2} \times {U^2})\text{，}$

$Z' = Z/(0.5 \times \rho \times {L^2} \times {U^2})\text{，}$

$ M' = M/(0.5 \times \rho \times {L^3} \times {U^2})\text{。}$

式中： $ R'$ ， ${\rm{Z}}'$ 及 $M'$ 分别为阻力系数、纵倾力矩系数， $\rho $ 为流体密度， $L$ 为艇长， $U$ 为潜艇航速。

5.2 航速的影响

图8～图10分别给出了不同潜深下潜艇阻力系数、垂向力系数及纵倾力矩系数随傅氏数Fr的变化规律，分析可得：

1）潜艇无因次潜深 ${H^ * }\text{≤} 1.6$ 时，阻力系数随 $Fr$ 呈先增大后减小再增大的趋势，而当潜深 ${H^ * }\text{≥} 1.9$ 时该变化规律并不明显，在无限潜深下，潜艇阻力系数随傅氏数逐渐减小并最终趋于稳定值。

2）不同潜深下，潜艇垂向力系数随傅氏数 $Fr$ 的变化规律基本一致，除无限水深工况外，均在 $Fr = 0.4$ 时出现峰值状态。由图也可看出，无因次潜深 ${H^ * } = 2.9$ 时，自由液面对潜艇向上的吸力作用较弱，且随航速变化垂向力系数变化幅度较小。

3）在无限水深下，潜艇纵倾力矩系数均远小于近水面航行工况，其值趋近于0；而在近水面工况下（ ${H^ * }\text{≤} 2.9$ ），当傅氏数 $Fr　\text{≤} 0.3$ 时，潜艇所受纵倾力矩系数 $M'$ 随航速的变化规律并不明显，而当 $Fr\text{≥} 0.3$ 时，纵倾力矩系数绝对值随航速逐渐增大。

5.3 潜深的影响

通过对比分析不同潜深下潜艇水动力的差异可直观分析潜艇水动力随潜深的变化趋势，潜艇阻力系数、垂向力系数及纵倾力矩系数随无因次潜深的变化曲线如图11～图13所示。

由图可知：

1）在低航速下（ $Fr \text{≤} 0.3$ ），潜艇阻力系数随潜深波动变化，而在 $Fr = 0.4$ 与 $Fr = 0.5$ 工况下，潜艇阻力系数随潜深增大而减小；在研究航速范围内，当潜深 ${H^ * } = 2.9$ 时，潜艇阻力系数与无限潜深下潜艇的阻力系数较为相近，即可认为当潜深 ${H^ * }\text{≥} 2.9$ 时，潜艇阻力性能受自由液面的影响较小。

2）相较于阻力性能，潜艇近水面航行时受到向上的垂向吸力，且一定航速下，随着潜深的不断增大，自由液面的吸力作用逐渐减弱，并在 ${H^ * }\text{≥} 2.9$ 后的较深水域趋于0。

3）由图13可知，在研究航速及潜深范围内，潜艇均受到埋首的纵倾力矩，而当 $Fr\text{≤} 0.3$ 时，潜艇所受纵倾力矩值较小，且随潜深的变化规律并不明显；在 $Fr = 0.5$ 工况下，潜艇纵倾力矩系数绝对值随潜深逐渐下降，潜深 ${H^ * } = 2.9$ 工况下潜艇纵倾力矩值约为 ${H^ * } = 1.0$ 时的3倍。

5.4 流场分析

由潜艇所受力及力矩分析可知，不同潜深下潜艇近水面航行时的水动力随航速的变化规律基本趋于一致，潜深越小，潜艇水动力受自由液面的影响也越大，因此，取定潜深 ${H^ * } = 1.0$ 下，不同航速下潜艇表面压力分布进行对比分析，如图14所示。在近自由液面工况下，低航速时（ $Fr = 0.1$ ）艇体表面压力主要受到静水压力的影响，而随着航速的增大，艇体上表面压力沿船长呈先增大后减小的趋势，在 $Fr = 0.5$ 航速下，艇体尾部上表面低压区较为明显。艇体前后的压力差使得潜艇在近水面高速航行时受到较大的埋首力矩，此现象与前述数据分析所得结论较为吻合。

同等航速下，潜深是影响潜艇水动力的重要因素，图15给出了航速 $Fr = 0.5$ 工况下，不同潜深时自由液面波高图。分析可知，在 ${H^ * } = 1.0$ 工况下，潜艇兴波较为明显，艇首波高为正，艇尾尾舵处存在较明显的负波高区域，而在艇后又会出现扇形波。随着潜深的不断增大，艇体兴波衰减较快，至 ${H^ * } = 2.9$ 潜深下，兴波只呈现较小的波高增值区域。

6 结　语

通过开展潜艇近水面水动力性能数值分析，得出如下结论：

1）基于切割体网格技术及RANS方法能够较为准确地模拟潜艇近自由液面水动力性能；

2）潜深直接影响潜艇水动力，当无因次潜深 ${H^ * } \text{≤}1.3$ 时，自由液面兴波较为显著，并直接影响艇体表面压力分布，使得潜艇水动力性能较深水状态存在较大差异，而当潜深 ${H^ * }\text{≥} 2.9$ 时，可认为潜艇水动力性能与深水状态无显著差异；

3）近自由液面条件下（ ${H^ * }\text{≤} 1.3$ ），潜艇所受阻力及垂向力随航速的增大均呈现明显的波动变化规律，此现象主要是由自由液面兴波使得艇体表面压力的变化所引起。