航行体沾湿区域对空泡尾流结构的影响

王威，王聪，李聪慧，宋武超

(哈尔滨工业大学 航天学院，黑龙江 哈尔滨 150001)

0 引言

水下武器借助通气超空泡的包裹，可以大幅度地减小水下阻力、提高航行速度。超空泡多相流动特性研究始终是各国学者关注的热点问题，其中空泡尾流特性是通气超空泡多相流动研究的重要方面，空泡尾部的双涡管泄气现象多在低弗劳德数Fr条件下产生，水洞实验也易于观察，因此针对双涡管尾流特性的研究备受各国学者青睐。

从工程应用的角度看，带有后体的超空泡航行体多相流特性也是需要考虑的问题，袁绪龙等[7]和陈伟政等[8]通过水洞实验分析了沾湿区域及重力对空泡形态的非对称性影响。航行体后体与流动介质相接触[9]，沾湿区域是维持航行体受力平衡的关键。沾湿区域从空泡下部穿出，必然会影响空泡多相尾流结构，这与无后体的自由闭合空泡尾流存在较大差异。Yu等[10]对航行体后体沾湿区域的升力进行了研究，分析了升力突变的原因。Wang等[11]研究了超空泡航行体的尾流结构，采用定常计算方法定性分析了空泡尾部流场压力和涡量分布特性，指出了沾湿区域对涡量方向的影响，并给出了涡管内部流动参数的变化规律。

从相关文献报道可以看出，对圆盘空化器模型的研究多侧重于空泡尾流结构，对带有后体的航行体模型研究多侧重于流体动力，但研究后体沾湿区域对航行体空泡尾流结构的影响并不多见。数值模拟方法是研究超空泡尾流结构的重要方法，其优点是易于捕捉实验方法较难获得的流场参数变化特性，分析流场参数非定常变化。

本文利用数值模拟方法，在文献[11]的基础上，采用非定常计算方法，进一步定量分析低Fr条件下航行体攻角逐渐变化时沾湿区域形成过程中空泡尾流场压力及涡量分布的演化特性，揭示沾湿区域对空泡尾流结构的影响机理，给出沾湿区域比例与空泡尾流结构的转变关系。

1 数学模型及验证

1.1 航行体模型

根据低Fr条件(5

航行体模型如图1所示，该模型由圆盘空化器、圆锥段、圆柱段、圆柱形尾喷管组成。航行体模型圆柱段直径为D，航行体全长L=12D，锥段长度为L/2，圆盘空化器直径Dn=0.3D，航行体质心位置在xc=6.5D处，坐标系原点建立在航行体质心。

图1 计算模型Fig.1 Computation model

1.2 计算流域边界条件

计算流域的边界如图2所示(图中p∞为远场边界压力)，计算域为圆柱形，直径为20D，长度为48D，重力方向沿y轴负向，计算流域采用六面体结构化网格划分。航行体表面采用壁面边界，在空化器后部设有通气孔，采用质量入口条件；流域入口和流域圆柱面边界都采用速度入口条件，流域出口采用水深5 m压力出口条件。

图2 计算域及边界条件Fig.2 Computational domain and boundary condition

Zhou等[12]研究了无限流域尺度对数值模拟结果的影响，指出p∞的等值线应在计算流域内部。本文利用p∞的等值面判断三维无限流域数值模拟对流域的最小尺度要求。从图2可以发现，p∞的等值面并未与计算流域边界交叉，且在计算流域内部，表明本文的流域边界尺度用于模拟无限流域情况下的通气空泡流动是合理的。

比较不同网格密度条件下航行体全湿状态时的流体动力系数随网格数量Ne的变化，结果如图3所示。图3中，CD为阻力系数，CL为升力系数。

图3 不同网格密度下航行体的流体动力系数Fig.3 Hydrodynamic coefficients of vehicle under different grid densities

从图3中可以看出：在0.6×106～0.9×106网格密度条件下，航行体的阻力系数CD和升力系数CL均随着网格数量而变化；在0.9×106～1.2×106的网格密度条件下，航行体的阻力和升力系数不再随网格数量变化。采用0.9×106以上网格数量进行数值模拟时，计算结果不会因为网格疏密程度而造成计算误差，综合考虑数值计算速度等因素，本文选取网格数量Ne=1.0×106的条件进行数值模拟研究。

1.3 模型验证

本文的通气率系数为

Cq=Q/(v∞D2)，

(1)

式中：Q为气体体积流量。

在相同的Cq和v∞条件下，将文献[10]的水洞实验结果和本文数值模拟结果对比，如图4和图5所示。图4分别给出了数值模拟和实验空泡形态发展过程，上方是数值模拟结果(本文以空气体积分数为0.5作为空泡界面)，下方是实验结果，可见空泡尾部都出现了涡管泄气的现象。

图4 数值模拟结果和实验结果对比Fig.4 Comparison of simulated and experimental results

图5 数值模拟和实验的空泡尺度对比Fig.5 Comparison of simulated and experimental cavity dimensions

空泡尺度的对比结果如图5所示。由图5可以看出，空泡长度Lc和空泡直径Dc的变化曲线具有很好的一致性，并且数值模拟结果和实验结果的空泡尺度相差在6%以内。图4和图5表明，数值模拟结果与实验结果吻合较好，选用的数学模型可以很好地预测通气空泡多相流动现象。

2 模拟结果分析与讨论

2.1 空泡尾流结构

航行体在不同工况下的空泡形态如表1所示。

表1 不同工况下航行体的空泡形态

Tab.1 Supercavity shapes under different working conditions

比较表1中工况1和工况2，发现：在无沾湿区域条件下，通气超空泡的尾流在重力环境中呈现双涡管的泄气形式；工况1仅起对比作用，表明在给定流速条件下，不计重力环境的压差因素时，空泡尾部为单一涡管的流动结构。

航行体后体沾湿区域主要为航行体提供纵平面内平衡的升力作用[13-14]。比较表1中工况2、工况3、工况4可以发现，随着攻角α的增加，航行体沾湿区域面积逐步增大，航行体超空泡尾流结构从双涡发展到三涡，表明沾湿区域对空泡尾流结构存在影响。

表1中工况5仅起与工况4的对比作用，给出了不计环境压差因素时航行体的超空泡尾部形态。比较工况4、工况5可以发现，在相同攻角(α=1.8°)条件下，两种工况中航行体的超空泡形状和沾湿区域都非常相似。二者环境压差不同，进一步表明航行体沾湿区域是影响超空泡尾部流场结构的重要因素。

航行体直线运动时空泡形态演化过程如图6所示。由图6可见，0.50 s以后，通气超空泡已经包裹整个航行体，形态尺度趋于稳定，空化数σ≈0.028，此时压力曲线平直，整个航行体表面压力分布趋于平稳，空泡内部压力差别较小[11]。

图6 航行体空泡演化Fig.6 Evolution of cavity of vehicle

本文对压力系数Cp及涡量Ω分别定义如下：

(2)

(3)

式中：p为流域压力；Ωx、Ωy、Ωz分别表示x轴、y轴、z轴方向的涡量，涡量的符号按照右手螺旋定则取值；vx、vy、vz分别表示x轴、y轴、z轴方向的流动速度。

空泡在0.30 s形成自由闭合，航行体从0.30 s开始匀速改变攻角至0.34 s结束，航行体的攻角α从0°增加到1.8°. 航行体尾流场的压力及涡量变化如图7所示。

图7 空泡尾部压力及涡量云图Fig.7 Contours of pressure and vorticity at the tail of supercavity

统计图7中空泡尾部流场最右侧截面在航行体攻角改变过程中空泡尾涡的非定常变化，结果如表2所示。表2中S0为沾湿区域面积，S为航行体表面积。

表2 空泡尾流横截面涡方向的非定常变化

Tab.2 Unsteady variation of vortex direction in cross section of cavity wake

从表2中可以看出：t取值为0.300～0.310 s时，空泡尾部的双涡管逐渐形成，流域环境压差使空泡尾部主要以双涡管闭合为主，虽然有小部分沾湿区域影响，但S0/S只有0.87%左右，空泡尾部的涡量方向在此期间保持一致，为逆时针和顺时针，如图7(b)所示；t取值为0.310～0.319 s时，空泡尾部双涡逐渐消失，随着S0/S的不断增大，空泡尾部形成四涡，该情况在圣安东尼瀑布实验室的研究中也有出现[5]，实验指出四涡情况发生在超空泡尾流从回射流到双涡管之间的转变过程中，由于数值模拟研究中选择空气体积分数不同以及数值方法在处理细微气泡的旋涡结构时还与实验结果存在一定差异，故数值模拟暂时还不能精细地反映实验中涡管细节；t取值为0.319～0.321 s时，随着沾湿区域对压差影响的增大，四涡中的上方双涡逐渐减小，下方双涡逐渐增大，涡量云图如图7(c)所示；t取值为0.321～0.322 s时，随着S0/S的进一步扩大，空泡尾流上方双涡逐渐被下方双涡取代；t取值为0.322～0.340 s时，空泡尾部上方的涡管尺度逐渐增大，形成三涡管尾流，上方涡管流域的涡量很小，下方双涡的涡量保持顺时针和逆时针方向，如图7(d)所示；t取值为0.340～0.440 s时，航行体攻角保持固定，空泡形态趋于稳定，尾涡形态如图7(e)所示。

2.2 空泡尾部压力

航行体尾部涡管附近的流线方向如图8所示，从图8(a)和图8(b)中的流线方向对比可以发现，沾湿区域致使空泡尾部横截面中涡管附近的流线方向发生了改变。

图8 空泡尾部涡管附近的压力和流线Fig.8 Pressure and streamline near vortex tube at cavity tail

图9 攻角、沾湿区域和空泡尾流高压区无量纲压力差系数随时间变化曲线Fig.9 Changes of angle-of-attack, wetted area dimensionless pressure difference in the high pressure region of wake flow of cavity with time

由图9可以看出：t取值为0.300～0.318 s时，受流域压差影响，上方高压区的压力小于下方，双涡管之间的流动方向向上，在空泡尾部形成一对中心上卷的双涡；t取值为0.318～0.322 s时，航行体的沾湿区域开始逐渐影响尾流压力分布，使上下方高压区逐渐融合，压力差减小；t取值为0.322～0.340 s时，攻角不断增加，沾湿区域绕流的诱导作用对尾流压力分布的影响进一步扩大，高压区再次分成上下两个部分，上方高压区的压力大于下方，压力差值逐渐增加，双涡管之间的介质转为向下流动，在空泡尾部形成一对中心下卷的双涡；t取值为0.340 s后，航行体攻角不在变化，沾湿区域趋于稳定。

2.3 空泡尾部涡量

为进一步分析四涡出现期间沾湿区域对尾涡方向的影响规律，考虑到流域环境压差和沾湿区域对空泡尾涡存在相反方向的影响，可以通过空泡尾部上下流线的速度分布进行深入研究[7-8]，如图10所示。空泡纵向平面上围绕整个空泡闭合轮廓线的环量Γ[15]近似表示为

(4)

式中：v为空泡轮廓线速度矢量；s为空泡轮廓线长度矢量；l为空泡轮廓线周长；vx1和vx2分别为空泡尾部上、下表面流线的x轴方向速度。

图10 空泡尾部上下流线的速度vx1和vx2Fig.10 Velocities vx1 and vx2 of upper and lower streamlines at cavity tail

由于空泡上下表面的速度差异在空泡尾部处较大，用vx1曲线和vx2曲线所围的面积可以近似地表征空泡尾部的涡量变化，如图11所示。图11(a)中空泡刚好闭合在航行体尾部，受流域环境压差作用，空泡下表面的速度略大[11]，形成双涡尾流。图11(b)中空泡尾流既受环境压差影响，又受到不断扩大的沾湿区域影响，空泡下表面速度vx2受沾湿区域阻挡逐渐小于vx1，涡量方向(顺时针区域和逆时针区域的面积之差)变化较快，空泡尾部形成了四涡流动结构。图11(c)中，空泡尾流随着沾湿区域的进一步扩大，空泡上表面的速度逐渐大于下表面，涡量方向发生变化，再次形成方向相反的双涡尾流。

图11 空泡尾部上下流线的速度比较(g=9.8 m/s2)Fig.11 Comparison of velocities of upper and lower streamlines at cavity tail(g=9.8 m/s2)

进一步研究沾湿区域和流域环境压差对空泡尾涡影响的主次关系，对比分析了不考虑环境压差条件下(如工况1，假设g=0)航行体攻角变化时，沾湿区域对空泡尾部涡量的影响。如图12所示，0.300 s时航行体攻角α=0°，空泡尾部不受沾湿区域影响，航行体空泡尾部上下速度差异较小，空泡横截面上有4个尾涡。0.340 s时航行体攻角α=1.8°，沾湿面积比例为36%左右，沾湿区域影响了航行体尾部的流速分布，空泡尾部出现对称分布的涡量。

图12 空泡尾部上下流线的速度比较(g=0 m/s2)Fig.12 Comparison of velocities of upper and lower streamlines at cavity tail(g=0 m/s2)

将vx1和vx2沿x轴方向积分，两条曲线所围的无量纲面积A为

(5)

用A的统计结果分析比较考虑沾湿区域引起的流速差异和环境压差引起的流速差异，比较二者对空泡尾部涡量方向的影响关系，如图13所示。

图13 沾湿区域对空泡尾涡的影响规律Fig.13 Influence law of wetted area on cavity wake vortex

从图13中可以看出，沾湿区域影响的曲线和流域环境压差影响的曲线相交在S0/S处于7～10%之间：当S0/S<7%时，环境压差对空泡尾部涡量方向起到主导作用，呈现一对中心上卷的双涡；当7%10%时，沾湿区域绕流的诱导作用对空泡尾部涡量方向起到主导作用，呈现一对中心下卷的双涡。另外，由于流动过程中存在多相流之间复杂的扰动影响，两条曲线交点处S0/S的数值是一个范围。

3 结论

本文研究了超空泡航行体沾湿区域对空泡尾流结构的影响，对空泡尾流结构的压力及涡量的非定常变化进行了深入分析，得到了沾湿区域对空泡尾流结构的影响规律。主要得出以下结论：

1)揭示了低弗劳德数条件下通气超空泡尾流结构特性。航行体无攻角时，空泡尾部形成一对中心上卷的双涡；航行体带有攻角的条件下，空泡尾部形成一对中心下卷的双涡。

2)给出了沾湿区域对尾涡方向的影响规律。在考虑重力环境影响的条件下，沾湿区域比例小于7%时，环境压差是空泡尾涡方向的主导因素；沾湿区域比例大于10%时，沾湿区域是空泡尾涡方向的主导因素。

3)指出了空泡的四涡尾流结构产生的原因。沾湿区域的比例处于7%～10%之间时，在环境压差和沾湿区域的双重作用下，空泡尾部流场中出现了四涡流动结构。

致谢本文撰写过程中得到了哈尔滨工业大学航天学院曹登庆教授和徐敏强教授的大力帮助与支持，在此对两位老师表示衷心的感谢！