0 引言

空化机理研究，最早可以追溯到1917年Rayleigh^[1]发展的气泡动力学理论。在理想无界流中，Rayleigh对单个球形气泡运动进行了分析，假设气泡内为饱和蒸汽压，建立了气泡运动的基本方程。Plesset^[2]考虑了非凝结气体、表面张力、流体黏性的影响，得到了更全面的单泡运动方程，称为Rayleigh- Plesset方程。在理论模型中，气泡假设为球形，从而可以比较方便的描述气泡半径变化。然而在很多实际情况下，如空泡近壁溃灭、多泡相互耦合运动中，空泡会发生非球形变形，理论模型就会出现很大偏差。

空泡实验研究可以弥补理论研究的不足，同时也能揭示一些更丰富的现象。Lauterborn^[3]通过实验研究了近壁单泡的溃灭过程，发现空泡溃灭会诱导冲击波现象。Vogel等^[4]测量了空泡溃灭过程中的流场分布，指出空泡近壁溃灭会诱导微射流，射流强度与近壁无量纲距离具有紧密关系。Gonzalez-Avila等^[5]测量了溃灭过程中的冲击压强，发现瞬间压强最高可达1000个大气压的量级。与此同时，空泡运动的数值模拟方法也在不断发展和完善。鲁传敬^[6]采用边界元积分方法，模拟了三维空泡的输运和溃灭过程，获得了单泡溃灭的辐射噪声。胡影影等^[7]采用界面捕捉方法，数值分析了近壁单泡的溃灭形态。Zhang等^[8]在数值模型中考虑了液体介质的压缩性，分析了空泡溃灭诱导压力波的传播过程。

多泡间的耦合作用也是影响空泡运动的重要因素。Bremond等^[9]、Ida^[10]通过实验研究了多个泡的生长和溃灭过程，发现外侧汽泡的运动领先于内侧汽泡，泡间距的减小或泡数目的增多都会抑制空泡的运动。Bui等^[11]对37个空化泡的近壁溃灭过程进行了数值模拟，结果显示泡群总体由外向内溃灭，他们将这种由外向内的现象称为泡群的遮蔽效应。张凌新等^[12]采用直接数值模拟方法对多泡的溃灭过程进行了模拟，研究发现内层泡的溃灭运动滞后于外层泡，内层泡诱导的溃灭压力要远高于单泡溃灭压力。Quinto等^[13]利用激光多点聚焦激发空泡，研究了同一平面上25个空泡的生长和演化过程，他们同样发现了遮蔽效应的存在。

空泡溃灭时，会对附近的壁面形成冲击作用。关于冲击作用的来源，学术界主要有两种看法，一是微射流冲击，二是压力波冲击。Kornfeld等^[14]的研究中首次提出了射流冲击理论，他们认为空泡溃灭的瞬间会产生高速射流，高速射流是空蚀的主要成因。压力波冲击理论中一个重要的基础是空泡势能理论，势能理论可追溯到1963年，根据Hammitt^[15]的研究，空泡势能与空泡体积以及空泡内外的压力差值有关。Fortes等^[16]针对单个空泡，发展出了一套基于空泡势能的近壁冲击评估方法，他们首先求解空泡溃灭释放的球形压力波信号，并进一步计算得到压力波所含能量以及波能转换率，不过他们的方法仅限于球形空泡运动分析。近年来，学者们对冲击波现象^[17]和回射流现象^[18-19]比较关注，但是研究多数还集中在单泡近壁溃灭上，对多泡近壁溃灭仍然还缺乏研究。

本文将在前人研究的基础上，对单泡以及多泡的溃灭过程进行数值模拟，并采用能量分析方法探讨溃灭过程中的压力波能量变化特征。首先，在数值建模中考虑液体的压缩性，建立空泡运动的数值模拟方法。在此基础上，数值模拟空泡溃灭过程，监测压力波的传播和流场中的能量变化。最终，基于能量演化特征，分析能量的传递路径，并给出压力波能的转化率。

1 数值方法及验证

考虑两相介质的流场，一相为液体水，一相为蒸汽。空泡界面通过VOF方法来捕捉，界面处施加表面张力。两相流的控制方程如下所示：

$ \frac{{{\rm{d}}\rho }}{{{\rm{d}}t}} + \rho \nabla \cdot \mathit{\boldsymbol{U}} = 0 $

(1)

$ \rho \frac{{\partial \mathit{\boldsymbol{U}}}}{{\partial t}} + \rho \mathit{\boldsymbol{U}} \cdot \nabla \mathit{\boldsymbol{U}} = - \nabla p + \nabla \cdot {\bar{\bar{\tau }}} + \sigma k\mathit{\boldsymbol{N}} $

(2)

$ \frac{{\partial \alpha }}{{\partial t}} + \mathit{\boldsymbol{U}} \cdot \nabla \alpha = 0 $

(3)

其中，U和p分别为速度矢量和压强，${\bar{\bar{\tau }}}$为黏性应力张量，k是表面曲率，σ是表面张力系数，N是表面的单位法向矢量，α是液相体积分数。流体的密度ρ由两相物质密度加权计算：

$ \rho = \alpha {\rho _1} + (1 - \alpha ){\rho _2} $

(4)

这里，ρ₁为液相的密度，ρ₂为汽相的密度。

$ {\rho _1} = {\rho _{10}} + \frac{p}{{c_1^2}} $

(5)

$ {\rho _2} = {\rho _{20}} + \frac{p}{{c_2^2}} $

(6)

这里，ρ₁₀=998和ρ₂₀=0分别为液相和汽相的密度常数，c₁和c₂为液相和汽相的声速，液相为固定声速，气相声速在等熵流假设下计算给出。为了模拟蒸汽泡的行为，设泡内压强为恒定的饱和蒸汽压p₀。

另外，为了研究溃灭过程中的能量变化，定义三种能量，分别为空泡势能、流场动能以及压力波能^[16]，它们的定义公式如下所示：

$ {E_{{\rm{pot}}}} = \frac{{4{\rm{ \mathsf{ π} }}}}{3}{R^3}\left( {{p_\infty } - {p_0}} \right) $

(7)

$ {E_k} = \int {\frac{1}{2}} \rho {U^2}{\rm{d}}V $

(8)

$ {E_{{\rm{wave }}}} = \frac{{4{\rm{ \mathsf{ π} }}{r^2}}}{{{\rho _1}{c_1}}}\int {{{\left( {p - {p_\infty }} \right)}^2}} {\rm{d}}t = \int {\frac{{{{\left( {p - {p_\infty }} \right)}^2}}}{{{\rho _1}c_1^2}}} {\rm{d}}V $

(9)

其中，p_∞为环境压强，R为空泡半径，其初始值为R₀，对应的初始空泡势能为E_{pot_0}。为了衡量压力波能的转化能力，进一步可以定义波能转化率为：

$ \eta = \frac{{{E_{{\rm{wave }}}}}}{{{E_{{\rm{po}}{{\rm{t}}_{\rm{0}}}}}}} $

(10)

控制方程组采用有限体积法离散，采用PISO算法进行变量耦合求解，计算方法在OpenFoam平台下实现，具体的步骤可以参考文献[8, 12]。为了模拟单泡以及多泡的工况，设计了一个圆柱形的计算域，如图 1所示。计算域四周距离空泡足够远，顶面距底面高为25倍的空泡直径，圆柱域直径为20倍空泡直径，底面为壁面边界，四周及上顶面为压力边界，压强为p_∞=1.0×10⁵ Pa，速度为零梯度条件。初始泡静止，初始半径为R₀=2 mm，泡内压强为p₀=3540 Pa。计算步长通过Courant数控制，Courant数取为0.2。空泡中心位置距下底面的距离为L，当L远大于空泡半径时，壁面对空泡的影响很小，可以作为远场运动工况，当L较小时，为近壁运动工况。定义无量纲距离参数γ：

$ \gamma = L/{R_0} $

(11)

首先，对单泡远场溃灭和近壁溃灭过程进行了数值模拟，并与相应的理论和实验结果进行了对比。对于远场溃灭，空泡基本保持球形，可以由Rayleigh方程给出其解析解，根据理论结果，空泡溃灭的时间为：

$ {T_c} = 0.915{R_0}\sqrt {\frac{{{\rho _1}}}{{{p_\infty } - {p_0}}}} $

(12)

将空泡半径和演化时间进行无量纲化：

$ {R^*} = R/{R_0} $

(13)

$ {T^*} = T/{T_c} $

(14)

数值模拟的结果和解析解的结果如图 2所示。在数值模拟中，检验了不同网格精度对结果的影响，以一个泡直径方向覆盖的网格数为量化指标，当网格数为20时，在空泡溃灭后期的结果与解析解偏差较大，当网格数达到30以上时，数值模拟能给出较为精确的结果。在图 3中，给出了近壁空泡的溃灭形态。当空泡距离壁面较近时，由于上下压差的不对称，造成空泡上部运动要快于空泡下部，从而形成指向壁面的射流，空泡形态也逐渐发展为凹形塌陷。图 3结果中，给出了空泡形态在不同时刻的数值结果，并与Kling等^[20]的实验观测结果进行了比较，无论从趋势以及细节上看，数值结果都达到了较好的吻合度。

2 结果与讨论 2.1 单泡溃灭

在单泡研究中，首先对远场溃灭过程进行了模拟，随后研究了不同近壁距离下的空泡溃灭。模拟在同一计算域中完成，通过改变不同近壁距离来实现工况的改变。远场工况下，泡离壁面的距离足够远，γ=12.5。近壁工况下，分别研究了γ=4.0、3.0、2.0、1.85、1.65、1.5的运动。

图 4给出了远场空泡溃灭过程中的压力云图和流线图。可以看出，在伴随空泡溃灭的过程中，全场压强在发生剧烈变化。首先，在空泡溃灭完成前，液相主要是指向空泡中心的径向运动，由于半径小的地方速度更快，所以空泡周围压强比较高。在空泡溃灭瞬间，径向流动在空泡中心对撞，产生冲击波，并在随后的过程中向四周呈球状传播。以泡中心为原点，在半径位置为2~8 mm处，设置了7个观察点来监测压力波传播，其结果如图 5所示。可以看出，压力波依距离远近，先后通过监测点，压力波呈现为一种近似对称的分布。距离近的监测点，压力峰值高，距离远的监测点，压力峰值低。根据压力峰值的时间和位置，可以计算出压力波的波速，其波速约为1500 m/s，与液相声速接近。这些是压力波的一些基本特征。

随后，我们对全场中的能量变化进行了分析，其结果如图 6所示。从趋势上看，随着空泡的收缩溃灭，其蕴含的势能逐渐变小。空泡的收缩造成了液体径向流动，全场动能逐渐升高。势能与动能间存在此消彼长的关系，当空泡溃灭瞬间，势能消失为0，动能达到其最大值，其值基本相当于初始时刻的空泡总势能。空泡溃灭后，动能陡然下降，与此同时，压力波释放，压力波能快速上升，其最大值约为初始空泡总势能的30%。在随后的时间里，动能下降，压力波能下降。在这里，我们基本可以判断出能量的传递路径。首先，在空泡溃灭过程中，空泡势能近乎全部转化为流场动能，空泡溃灭瞬间，动能中的一部分转化为压力波能，其余部分在随后的运动中耗散掉。

进一步，我们研究了不同近壁距离下的空泡溃灭过程。近壁距离γ=1.5下的能量变化在图 7中给出，其整体变化趋势与图 6类似。有所不同的是，该工况下产生的压力波能峰值较小，大约为初始空泡势能的10%左右。同时，在压力波能达到峰值时，动能的水平高于压力波能。分析其原因，在近壁溃灭中，会诱导微射流，因此在溃灭瞬间，全场动能有一部分是保存在射流中的，这部分射流动能并不会贡献给压力波能，而是沿自身的路径随后与壁面发生作用，因此近壁溃灭中压力波能小，动能留存时间长。其他近壁工况下的压力波能及其转化率见表 1，可以看出基本的转化趋势，远场转化率最高，达到30%，随着近壁距离减小，转化率逐渐降低，在γ=1.5时，其转化率约为10%。

2.2 多泡溃灭

多泡运动要更为复杂一些，泡间相互作用不仅会影响相邻空泡运动，而且会对整体动力学过程产生重大影响。在本文的工作中，简化了部分物理参数，将泡间距设定为2.5 mm，空泡呈规则立方体排列，分别研究了8个、27个以及64个空泡运动。在多泡情况下，定义最下层空泡中心位置距壁面的距离为近壁距离。

图 8给出了27个泡远场溃灭的过程，可以明显看出泡与泡之间的相互作用，中间的泡在开始阶段运动缓慢，外层的泡首先溃灭，在溃灭过程中，外层空泡也会产生凹形变形，伴随着局部微射流和局部高压。中心空泡最后溃灭，类似于单泡远场溃灭，27个泡溃灭也产生了球形的冲击波。当然，此时的总体能量要远高于单泡溃灭。如图 9所示，远场溃灭中，其动能峰值约为单泡动能峰值的27倍，压力波能约为单泡波能的14倍。近壁情况下，整体运动时间有所延长，动能峰值几乎不变，波能峰值有所下降，约为远场工况的一半左右。因此，从总体趋势上看，无论是单泡还是多泡，近壁溃灭产生的波能要小于远场溃灭。

压力波是由于高速流体撞击而形成。在单泡远场溃灭中，流动可以看作球对称，空泡溃灭诱导的射流均指向空泡中心位置，空泡四周的射流同步产生、同步在中心点对撞，所以其压力波能转化率最高。而在多泡远场溃灭中，外层泡在溃灭的同时，诱发了局部微射流，这些射流既有指向中心位置的，也有偏离径向方向的，有一部分动能并没有参与到中心空泡溃灭的流动撞击中去，所以多泡远场溃灭的波能转化率比单泡溃灭要低。而在多泡近壁溃灭中，除了外层空泡溃灭诱导的局部射流，还诱发了指向壁面的主射流，射流会形成自身的壁面冲击，但是对压力波能的转化贡献就弱了一些。所以从总体上说，近壁溃灭强化了射流，而抑制了压力波能的释放。

最后，给出三种多泡工况下的压力波能演化过程，如图 10所示，其中纵坐标已用各自的空泡起始势能无量纲化，分别对应8泡、27泡、64泡。随着泡的层数增多，其溃灭时间也逐渐延长，符合遮蔽效应特征。从波能转化率值来看，几种工况下的转化率均稳定在10%~13%左右。联系到单泡在近壁γ=1.5溃灭中的波能转化率也在10%左右，因此波能转化率可能受近壁距离的影响较大，但与其他因素相关性较小。当然，现在研究的空泡分布，其长宽高比较接近，对于极其狭长的空间分布也许会有不同的结果。目前的结果展现了一种可能性，即如果波能转化率是一个有规律参量，那么就可能以波能转化率来评估更大规模的空泡云团，进而建立一种服务于工程的压力波能评估方法。

3 小结

本文通过直接数值模拟，研究了单泡以及多泡的溃灭过程，在监测能量变化的基础上，分析了溃灭过程中的能量传递机制，并比较了远场溃灭和近壁溃灭的区别。

研究发现，无论是远场溃灭还是近壁溃灭，其能量演化存在着比较接近的传递路径，即首先在空泡溃灭完成前，空泡势能近乎全部转变为流场的动能，在溃灭瞬间，部分动能转变为压力波能，部分动能最终耗散掉。所以压力波能并不直接与空泡势能相关，而是与流场动能相关。

近壁溃灭中，总体上会强化射流，抑制压力波能的释放。所以近壁溃灭中，其压力波能转化率相对较低，而流场动能留存时间更长。本文对单泡以及立方体布置的8泡、27泡、64泡的近壁溃灭过程进行了研究，发现在近壁距离γ=1.5下，它们的波能转化率约在10%的量级。