0 引　言

随着技术的发展，以空投鱼雷入水、水上飞机降落、返回舱回收为代表的实际应用问题不断增加，这推动着人们对入水问题进一步开展研究。近年来，由于超空泡技术的发展，跨介质航行体的入水过程得到了极大的关注。在高速入水过程中，由于空化现象的产生形成了包裹运动体的空泡，可实现对航行体大幅减阻的效果。因此对高速入水过程的研究有助于掌握航行体的运动规律和空泡演化过程，对水中兵器的发展具有重要意义。

自从Von Karman^[1]和Wangner^[2]在20世纪对物体入水过程进行研究以来，人们对入水现象进行了深入的研究。针对空泡理论研究的情况，Logvinovitch^[3-4]根据理论分析和实验结果提出了空泡截面独立膨胀原理，获得了空泡轮廓计算公式。Savchenko^[5]结合大量水中实验数据提出了对高速航行体产生的空泡计算公式，可应用于小空化数下的空泡轮廓计算。Serebryakov等^[6-7]将前人对高速航行体的研究进行总结，并推导了亚音速、跨音速和超音速等情况下的空泡计算公式。对于入水过程的实验研究方面，Worthington^[8]研究了小球入水的实验过程，并记录了空泡的演化过程。Hrbues^[9]应用火炮发射超空泡射弹达到了1540 m/s的量级，并对实验结果进行了分析。Truscott^[10]总结了针对入水过程的理论、实验和数值仿真等方面的研究成果。Schaffar^[11]进行了高速入水实验，并与仿真结果进行对比。国内郭子涛等^[12]进行了不同头型弹丸的水平入水实验并与仿真结果进行了对比。蒋运华等^[13]对运动体在小扰动下的入水过程开展了实验研究，对比了空泡在有无扰动下的形态演化规律。施红辉等^[14]对细长体倾斜入水过程开展了实验研究，分析了细长体的入水过程及其与自由面的相互作用。高建国等^[15]开展了弹体高速斜入水的研究，对比了不同尺寸的弹丸入水过程，得到了面闭合时间与弹体入水速度的关系。孙玉松等^[16]总结了国内外对结构体入水冲击问题的研究概况，从理论研究、数值计算和试验研究3个方面进行了梳理分析。

本文采用VOF模型和重叠网格技术对圆柱体高速入水过程进行仿真模拟，得到了圆柱体入水过程中的空泡形态演变、入水压力和姿态变化曲线，对圆柱体入水所受到的流体动力及其运动情况进行分析，获得了圆柱体入水过程的运动规律。

1 仿真模型 1.1 建立计算模型

建立包含背景网格和重叠网格的圆柱形求解域，如图1和图2所示。求解域侧面为滑移壁面边界，底部为压力出口边界，顶部为压力入口边界。高速入水圆柱体外形长度为L=25.4 mm，直径D=12.66 mm。求解域直径为8L，空气域高度为5L，水域高度为9L。

1.2 流体控制方程

采用VOF方法对入水过程中的流动进行模拟，同时捕捉空泡轮廓和自由表面的演化过程。VOF方法采用均质平衡流理论将计算域中的气、汽、液三相流体当作单一混合流体处理，通过求解各相在流体域中的体积分数从而得到各相在流动中的分布情况。本文描述流体运动的控制方程为连续性方程和动量方程组，即

$ \frac{{\partial \rho }}{{\partial t}} + \frac{{\partial (\rho {u_i})}}{{\partial {x_i}}} = 0, $

(1)

$ \begin{split}\frac{{\partial (\rho {u_i})}}{{\partial t}} +& \frac{{\partial (\rho {u_i}{u_j})}}{{\partial {x_j}}} = - \frac{{\partial p}}{{\partial {x_i}}} + \frac{\partial }{{\partial {x_j}}}[\mu + {\mu _t}]\times \\ & \left(\frac{{\partial {u_i}}}{{\partial {x_j}}} + \frac{{\partial {u_j}}}{{\partial {x_i}}} - \frac{2}{3}\frac{{\partial {u_k}}}{{\partial {x_k}}}{\delta _{ij}}\right) 。\end{split}$

(2)

其中：ρ为混合流体密度；u_i为各方向上的流体速度分量；p为压强；μ和μ_t分别表示层流粘度和湍流粘度。

流体密度和层流粘度通过各相体积分数来确定，即

$\begin{split} &\rho = {\alpha _g}{\rho _g} + {\alpha _v}{\rho _v} + {\alpha _l}{\rho _l}，\\ &\mu = {\alpha _g}{\mu _g} + {\alpha _v}{\mu _v} + {\alpha _l}{\mu _l}。\end{split}$

式中：ρ_g，ρ_v，ρ_l和μ_g，μ_v，μ_l分别为气、汽、液三相介质的密度和粘度；α_g，α_v，α_l分别为三相介质的体积分数。

1.3 湍流模型

入水过程是一个复杂的流体和运动体交互作用的多相流动过程，其表现为流体运动的强非线性，因此需要采用湍流模型对方程组进行封闭。本文对入水过程中湍流流动采用可实现的壳k-ε模型对流体方程组进行封闭。湍动能k方程和湍流耗散率ε方程为：

$\begin{split} \frac{{\partial (\rho k)}}{{\partial t}} + \frac{{\partial (\rho k{u_i})}}{{\partial {x_i}}} = &\frac{\partial }{{\partial {x_j}}}\left[ {\left( {\mu + \frac{{{\mu _t}}}{{{\sigma _k}}}} \right)\frac{{\partial k}}{{\partial {x_j}}}} \right] + {G_k} + {G_b}- \\ &\rho \varepsilon - {Y_m} + {S_k}，\end{split}$

(3)

$\begin{split}\frac{{\partial (\rho \varepsilon )}}{{\partial t}} + \frac{{\partial (\rho \varepsilon {u_i})}}{{\partial {x_i}}} = &\frac{\partial }{{\partial {x_j}}}\left[ {\left( {\mu + \frac{{{\mu _t}}}{{{\sigma _\varepsilon }}}} \right)\frac{{\partial \varepsilon }}{{\partial {x_j}}}} \right] + \rho {C_1}S\varepsilon- \\ &\frac{{\rho {C_2}{\varepsilon ^2}}}{{k + \sqrt {\nu \varepsilon } }} + {C_{1\varepsilon }}\frac{\varepsilon }{k}{C_{3\varepsilon }}{G_b} + {S_\varepsilon }。\end{split}$

(4)

其中：μ_t为湍流粘性系数；σ_k和σ_ε为湍流普朗特数；G_k为平均速度梯度引起的湍动能生成项；G_b为浮力产生的湍动能项；Y_m为可压缩湍流脉动的影响； ${C_1} = $ $ \max \left[ {0.43,{\eta \mathord{\left/ {\vphantom {\eta {\left( {\eta + 5} \right)}}} \right. } {\left( {\eta + 5} \right)}}} \right]$ ， $ \eta = {{Sk} \mathord{\left/ {\vphantom {{Sk} \varepsilon }} \right. } \varepsilon } $ ；C₂，C_1ε和C_3ε为常数；S_k和S_ε为源项。

1.4 空化模型

物体高速入水时，水在物体表面高速流动，从而引起局部压强降低，当达到水的饱和蒸汽压时就产生了空化现象。本文采用Schnerr-Sauer空化模型^[17]对空化现象进行建模，其质量输运方程及源项为：

$ \frac{{\partial ({\rho _v}{\alpha _v})}}{{\partial t}} + \frac{{\partial ({\rho _v}{\alpha _v}{u_i})}}{{\partial {x_i}}} = {\dot m^ + } - {\dot m^ - }，$

(5)

$ {\dot m^ + } = \frac{{{\rho _v}{\rho _l}}}{\rho }{\alpha _v}(1 - {\alpha _v})\frac{3}{{{R_B}}}\sqrt {\frac{2}{3}\frac{{({p_v} - p)}}{{{\rho _l}}}} \text{，} p \leqslant {p_v}，$

(6)

$ {\dot m^ - } = \frac{{{\rho _v}{\rho _l}}}{\rho }{\alpha _v}(1 - {\alpha _v})\frac{3}{{{R_B}}}\sqrt {\frac{2}{3}\frac{{(p - {p_v})}}{{{\rho _l}}}} \text{，} p > {p_v} 。$

(7)

其中：ρ_v为水蒸气密度；ρ_l为液体密度；ρ为按体积加权得到的混合密度；α_v为蒸汽相体积分数；p_v为蒸汽压强；p为当地压强；R_B为气泡半径。

1.5 运动方程

根据牛顿定律，描述圆柱体在三维空间中运动的动力学方程为：

$ m\frac{{{\rm{d}}{\boldsymbol{v}}}}{{{\rm{d}}t}} - F = 0，$

(8)

$ J\frac{{{\rm{d}}{\mathbf{\omega }}}}{{{\rm{d}}t}} - M = 0 。$

(9)

其中：v和ω分别为物体运动的速度、角速度；m和J为质量和转动惯量矩阵；F和M分别为物体受到的合外力和合外力矩。

1.6 离散方法

本文采用有限体积法仿真圆柱体垂直入水过程，利用SIMPLE算法对压力和速度耦合求解。其中对流项采用二阶离散格式，扩散项采用二阶中心差分格式，时域采用隐式离散方法进行求解。

1.7 模型验证

为了验证建立模型的正确性，采用文献[12]中的圆柱体入水实验数据与本文通过仿真模型得到的计算结果进行对比，速度-时间曲线、水中位移-时间曲线的仿真结果与文献试验结果如图3所示。可以看出，仿真得到的速度和位移变化情况与文献保持一致且误差较小，从而验证了计算方法的可行性和仿真模型的正确性。

2 结果及分析 2.1 空泡形态分析

图4为入水过程中圆柱体在不同位置处气相体积分数分布的演化过程。可以看出，圆柱体距离液面较远时，自由表面无变化。当圆柱体与液面接触后，柱体周围附近液面有所上升。随着圆柱进入液体继续向下运动，空泡开始形成并随着直径扩张，空泡长度也在逐渐增长。图中所示从左至右依次为t = 0，72.5 μs，183.07 μs，297.12 μs，535.24 μs，660.289 μs各时刻的分布情况。

Savchenko在文献[5]中总结了Logvinovich^[18]的研究成果，根据大量试验结果提出了预测空泡形态的经验公式：

$ \bar R = \left\{{\begin{array}{*{20}{c}} {{{(1 + 3\bar x)}^{{1 \mathord{\left/ {\vphantom {1 3}} \right. } 3}}},}&{\bar x < 2}，\\ \begin{gathered} (3.659 + 0.847(\bar x - 2) \hfill - 0.236\sigma {(\bar x - 2)^2}{)^{{1 \mathord{\left/ {\vphantom {1 2}} \right. } 2}}}, \hfill \\ \end{gathered} &{\bar x \geqslant 2} 。\end{array}} \right. $

(10)

式中： $ \bar R = {{R(x)} \mathord{\left/ {\vphantom {{R(x)} {{R_n}}}} \right. } {{R_n}}} $ ，R(x)为空泡半径，R_n为空化器半径； $ \bar x = {x \mathord{\left/ {\vphantom {x {{R_n}}}} \right. } {{R_n}}} $ 为归一化轴向坐标；σ为空化数。

经验公式得到的空泡轮廓与仿真结果如图5所示。可以看出，圆柱体运动至0.5L，1.5 L和2.5 L水深处的仿真空泡轮廓与经验公式基本吻合，因而再次说明了计算方法和模型的正确性。

2.2 流体动力分析

图6为圆柱体在入水过程中的阻力系数变化过程。在空气中运动时，圆柱体所受流体阻力较小。当运动至与液面接触瞬间产生了砰击现象，圆柱体阻力急剧上升达到峰值，然后又快速下降。随后圆柱体继续运动，阻力系数下降并逐渐稳定。图7为圆柱前端面中心压强随时间变化曲线。在入水前该点压强较低，撞水瞬间产生一个尖锐脉冲，入水后趋势逐渐平稳。

圆柱前端面上的压强分布如图8所示。圆柱体在空气中运动时，压强呈中心对称分布，中心处压强较高。入水较短距离后，中心处压强较低，端面压强沿径向呈先增大后减小趋势。随着圆柱体继续运动，端面中心区域内压强最大，沿径向向外逐渐减小。

2.3 速度变化分析

图9为圆柱体入水过程的速度随位移变化曲线。可知，圆柱体入水前由于空气阻力较小，速度衰减较为缓慢。当圆柱体运动至水面处时与自由液面发生砰击，运动速度急剧减小。随后圆柱体继续向水中运动呈自由减速，速度下降趋势平稳。

2.4 空化现象分析

图10为水蒸气相体积分数在计算域中的变化情况。在圆柱体运动至与液面接触后，圆柱体头部端面沾湿，端面边缘处由于流体分离出现高速射流，此处便产生了空化现象。随着圆柱体继续运动，在端面边缘空化持续进行。

3 结　语

通过采用计算流体力学方法对圆柱体高速入水过程进行仿真研究，分析了空泡演化过程、圆柱体所受流体动力及其速度变化规律，得到以下结论：

1）根据建立的仿真模型得到的仿真计算结果与试验结果一致，表明了仿真模型的正确性和计算方法的可行性。

2）在高速入水过程中，圆柱体头部受到了较大的冲击载荷。入水后，圆柱体头部载荷逐渐趋于平稳，撞水过程的持续时间在微秒量级。

3）在水面发生砰击后，圆柱体速度急剧下降，随后呈自由减速运动，速度衰减趋于稳定。

4）在圆柱体与水面接触时，空化现象首先发生在圆柱体前端面边缘射流处。随着圆柱体继续运动，在前端面边缘处持续产生空化水蒸气。