近壁面水下爆炸气泡运动的数值计算研究

李 健，林贤坤，荣吉利，项大林

(1. 广西科技大学 汽车与交通学院，广西 柳州 545006； 2. 北京理工大学 宇航学院，北京 100081)

当气泡在结构表面附近脉动时，在膨胀阶段被结构表面轻微地排斥开，而在坍塌阶段被结构强烈地吸引，并产生高速射流。射流的成因可以用著名的Bjerknes效应来解释。一般情况下，由于Bjerknes力作用，气泡会被水中结构吸引而向结构方向运动，在浮力、惯性力等因素的影响下，沿着气泡坍塌方向会形成高速射流，且射流会穿透气泡并冲击气泡壁的另一面，最终作用于结构造成结构的破坏。由于气泡脉动、坍塌、射流形成等一系列物理现象非常复杂，难以通过理论方法得出气泡运动过程的解析解，因此，气泡与结构相互作用主要是通过采用实验或仿真的方法进行研究。

实验方法主要是利用小药量炸药爆炸[1-2]或借助气泡发生器在水中形成气泡[3-4]，并将结构置于气泡附近，以此研究近壁面气泡的运动特性。数值仿真计算方面目前大都采用边界积分法[5-8]或基于有限体积法的非线性仿真软件[9-11]。边界元方法只在边界离散，大大降低了计算成本，但是在射流产生、气泡破碎及环形气泡形成方面有一定的局限性。有限体积法或流体体积法可以较好的模拟流体流动过程中的间断问题，其缺点是需要对每个空间位置进行网格离散，导致计算量大且介质间的交界面难以精确追踪、捕捉。

由于气泡与结构之间相互作用的现象非常复杂，目前仍有许多现象没被揭示，关于射流特征、环形气泡形成等方面的研究文献甚少，而基于有限体积法结合界面追踪、捕捉方法进行气泡特性的研究更是十分罕见。为此，为了对结构附近气泡运动数值模拟研究提供有益的补充，本文从气泡与刚性壁面相互作用的基本现象入手，系统地研究刚性壁面附近气泡的动力学特性，包括刚性壁面角度、气泡与壁面距离等参数对气泡运动和射流方向、速度的影响，旨在为相关研究提供参考。

1 理论背景

1.1 气泡运动初始条件

要模拟气泡的运动，需要先求出气泡的初始半径。目前关于气泡初始状态求解的方法较多，主要包括Rayleigh模型[12]和Frost和Harper[13-14]提出的气泡体积加速度模型，后者主要是分别建立了药包参数和流场压力以及气泡体积加速度与流场压力之间的关系，并利用压力等效关系最终确定气泡的初始半径与初始膨胀速度。此方法已经能够较准确地描述气泡的初始条件，但是考虑到从炸药点火到气泡脉动是一个连续的物理过程，为了更为准确的模拟整个水下爆炸过程，本文将基于MSC.DYTRAN软件，先采用一维模型模拟炸药的整个爆轰过程，通过将一维计算结果映射到三维流场中，用以确定三维气泡脉动的初始条件，以此实现炸药从爆轰到气泡脉动全过程进行数值计算，为水下爆炸数值计算提供一种新的解决方案。

为验证一维数值计算的正确性，将爆炸过程中流场不同距离冲击波压力峰值与经验公式进行对比。图1为装药量为55 g的TNT炸药在水下爆炸不同测点距离处冲击波的波形曲线，其中炸药采用JWL状态方程描述，各参数的选择可参见文献[1]。从图中可以看到，随着距离的增加，冲击波压力峰值衰减较快，而冲击波波形与真实情况相吻合，表现为波峰“尖锐”且后续无明显的尾随振荡。

图1 不同测点距离时冲击波波形曲线

表1为不同工况下冲击波压力峰值与经验公式的对比，通过对比可以看到，数值计算结果与经验公式具有较好的一致性，误差能够控制在7%以内，验证了数值计算的正确性。

表1 不同测点冲击波压力峰值仿真结果与经验公式的对比

1.2 气泡状态方程

由于气泡脉动过程中其内部压力远低于炸药爆轰过程中的压力，因此不太适宜采用JWL状态方程描述气泡，一般情况下，气泡脉动过程中，认为气体的压力仅和气泡的初始状态及体积有关[14]，即气泡内的压力P与气泡体积V的关系可以表示为：

(1)

(2)

其中：Pc为可冷凝气体的饱和蒸汽压，量级与大气压相当，一般可忽略；P0和V0分别为气泡形成时的初始压力和初始体积，W为装药量，γ为气体的比热比，对于TNT水中爆炸的爆轰产物而言，γ取1.25。

1.3 流场初始条件定义

在利用一维有限元模型完成了炸药爆轰与冲击波传播过程的数值计算后，可以获得气泡壁运动速度、气泡内压力等参数，若要进行三维气泡的数值模拟，还需要定义流场的初始条件(如静水压)与边界条件。为此，根据初值及边界条件特点，利用定义初始条件和边界条件的子程序接口EXINIT及EXFLOW2，开发了定义流场初值和边界条件的子程序，其思路就是根据模型中各单元节点的实际位置，将与之相对应的静水压力通过EXINIT和EXFLOW2子程序逐一赋到整个流场的单元节点上。

2 结果与讨论

2.1 模型验证

为验证三维气泡运动数值计算的正确性，将装药量为55 g的球形炸药置于水下3.5 m处，并将气泡脉动体积时程曲线与文献[1]中的实验进行对比，对比结果如图2所示。从图中可以看出，本文所采用的由一维向三维映射方法所计算的气泡体积与脉动周期较文献[1]中所采用的仿真方法与实验测试结果更为接近，验证了数值计算模型、流场初始条件与边界条件子程序开发、状态方程选取的正确性。

图2 气泡脉动半径时间历程曲线

为验证本文所建立数值计算模型在模拟刚性壁面附近三维气泡运动的有效性，根据文献[1]中的实验工况，即将55 g炸药置于水下3.5 m处，并在其水平距离0.4 m处垂直放置一刚性平板，通过高速相机捕捉各时刻气泡运动特性(实验结果如图3所示)。

根据实验测试工况，采用流-固耦合方法定义流体与结构的相互作用，建立近刚性壁面三维气泡运动模型，气泡与壁面相互作用过程的计算结果如图4所示。通过对比图3与图4发现，气泡在各时刻的运动均具有比较好的一致性，特别是在t=96 ms时，可以通过仿真结果清晰地看到气泡坍塌会发生于远离结构一侧，形成由右下指向左上的射流，射流最终击穿气泡而作用在结构上。

图3 55 g炸药在水下3.5 m距垂直刚性面40 cm处爆炸气泡运动实验结果

图4 55 g炸药在水下3.5 m距垂直刚性面40 cm处爆炸气泡运动仿真结果

2.2 特征参数对气泡运动的影响

为研究刚性壁面附件气泡的运动规律，本节将以气泡与垂直放置刚性壁面相互作用数值计算模型为基础，研究气泡与刚性壁面角度、距离等参数对气泡运动特性的影响规律。由于气泡坍塌、射流方向、射流速度等物理现象的产生与爆心初始位置、气泡脉动周期、气泡最大半径等关键参数有密切关系，因此以气泡的最大半径Rm为特征长度，中点处静水压力pamb为特征压力，则特征时间为t=Rm(ρ/pamb)0.5，气泡与刚性壁面垂直距离为d=L/Rm，速度的无量纲表达式为v/(Rm/t)，以此研究初始距离与角度对气泡运动特性的影响规律。

3.2.1 刚性面角度对气泡的影响

本节首先研究气泡与刚性壁面距离与角度变化时对气泡射流角度与射流速度的影响。图5为在气泡与刚性壁面初始距离为Rm且刚性壁面与水平面的夹角分别为0°、45°、90°、135°和180°时气泡膨胀过程。从图5可以看出，在Bjerknes力作用下，气泡在膨胀时会被刚性壁面吸引，且气泡最大半径不会随着刚性壁面角度的变化而变化，表现为不同工况下气泡脉动最大半径及膨胀到最大半径所用的时间基本相同。

图6为与图5各工况相对应的气泡收缩过程。从图中可以看出，当气泡与刚性壁面初始距离较近时(特征距离为Rm)，气泡受刚性壁面的影响较大，表现为气泡坍塌发生于远离结构一侧，且最终的射流方向基本垂直于刚性壁面，从射流最终击穿气泡的时间来看，射流发生于气泡收缩至最小体积并开始第二次膨胀阶段。此外，当刚性壁面与水平面夹角大于90°时(如图6(d)、(e)所示)，由于此时Bjerknes力沿竖直方向投影与浮力方向相反，会对气泡坍塌与射流的形成有影响，导致射流速度降低，表现为射流击穿气泡所需时间较其它工况要长。

无量纲时间: (0)0.00, (1)0.21, (2)0.42, (3)0.62, (4)1.09

图7为特征距离为Rm时，气泡射流方向与壁面法向的夹角的关系曲线。从图中可以看到，由于特征距离较近，Bjerknes力效果明显占优，最终气泡射流基本都是垂直射向结构。当刚性面垂直放置时，气泡射流与壁面法向夹角值最大，约为4°。

图7 特征距离为Rm时射流角度与刚性面角度之间的关系

图8为特征距离为Rm时，气泡射流速度与刚性面角度的变化关系曲线。从图中可以看到，随着刚性面与水平面夹角的增加，气泡射流速度会逐渐降低，其主要原因就是气泡坍塌与射流主要受Bjerknes力与浮力的共同作用，Bjerknes力方向主要是气泡与结构连线方向，而浮力就是沿着竖直方向，当Bjerknes力沿竖直方向分量与浮力作用方向相同时，所形成的合力能使气泡产生高速射流，而当二者方向相反时，合力导致气泡射流速度降低。

图8 特征距离为Rm时射流速度与刚性面角度之间的关系

图9、10分别为特征距离为2Rm时，气泡射流角度、速度和刚性面放置角度之间的变化关系。从图9可以看出，射流方向受刚性面角度影响明显，射流与刚性面夹角会经历先增加后减小的过程。对比图7与图9可以发现，当距离增加时，由于Bjerknes力作用效果减弱，导致气泡射流方向与结构夹角会有不同程度地增加。从图10中可以看出，当特征距离不变时，射流速度会随着刚性面角度的增加而减小，对比图8与图10可以发现，当距离增加时，气泡与结构之间的吸引力减小，导致射流速度降低，且在相同角度情况下，特征距离为Rm所产生的气泡射流速度约为2Rm时的3倍左右。因此，特征距离对气泡射流方向与速度影响明显。

图9 特征距离为2Rm时射流角度与刚性面角度之间的关系

图10 特征距离为2Rm时射流速度与刚性面角度之间的关系

3.2.2 特征距离对气泡运动的影响

前述计算结果表明，特征距离对气泡射流特性影响明显，且根据图7所示，当特征距离较小时，气泡与壁面最大夹角发生在壁面垂直放置时，为此，本节将研究刚性壁面垂直放置时，特征距离对气泡运动特性的影响。图11表示气泡中心与垂直刚性壁面距离分别为2Rm、3Rm和4Rm时气泡的收缩过程(特征距离为Rm工况参见图6(c))。从图中可以看到，当气泡与刚性壁面较近时，由于Bjerknes力作用效果显著，气泡被结构强烈吸引，随即产生指向结构的射流。当气泡与结构距离为2Rm时，Bjerknes力作用效果有所减弱，表现为气泡最终形成的射流方向与刚性面法向夹角增加，而当特征距离大于3Rm时，气泡无明显指向结构的射流产生，此时气泡受结构作用影响较小，基本可以忽略结构对气泡运动的影响。

无量纲时间：(0)1.09 (1) 1.33 (2) 1.58 (3) 1.80 (4) 2.03

图12 气泡与垂直刚性壁面特征距离与射流角度之间的关系

图12为气泡在垂直刚性壁面附近运动，其射流角度与特征距离之间的关系曲线。从图中可以看到，随着特征距离的增加，气泡与刚性面的夹角也随之增加，特别是当特征距离大于3Rm时，射流角度与垂直刚性面法向的夹角基本为90°，射流方向基本沿着竖直方向，与气泡在自由场中运动特性类似，表明此时刚性面对气泡运动影响较小，基本可以忽略。

3 结 论

采用一维有限元模型对爆轰与冲击波传播过程进行了仿真分析，以此为基础，将一维计算结果映射到三维流场中，确定了气泡脉动的初始条件，并对气泡与结构相互作用全过程进行了仿真分析，通过与实验结果的对比，验证了气泡初始条件确定方法、子程序开发与模型建立的正确性，以此模型为基础，通过仿真计算，研究了特征距离与刚性面角度对气泡运动特性的影响，分析总结相关规律，获得了如下一些主要结论：

(1) 在一维有限元模型利用JWL状态方程中模拟炸药爆轰与冲击波传播过程，并将其计算结果作为气泡脉动的初始条件，在三维仿真模型中采用Tait状态方程描述气泡，通过与经验公式和已有实验测试结果的对比，验证了此方法在三维气泡运动数值模拟的可行性与正确性。

(2) 由于Bjerknes力对气泡的作用力主要是沿着气泡与结构的连线方向，而浮力作用是沿着竖直方向，因此，刚性面的角度对气泡运动、射流方向与射流速度均有不同程度的影响。表现为当二者作用方向相同时，气泡射流速度较高，反之，气泡射流速度较低。

(3) 当特征距离在3Rm以内时，由于Bjerknes力作用占优，可以产生指向结构的射流。表现为当初始距离在Rm时，气泡射流角度与结构法向夹角基本重合，且射流速度较高，射流速度随着角度的增加而减小；当初始距离在2Rm时，气泡射流角度与结构法向夹角增加，最大角度出现在刚性面与水平面成120°时，且射流的速度随着角度的增加而减小。

(4) 随着特征距离的增加，射流方向与结构法向之间的夹角随之增加，而射流速度随着距离的增加而减小。表现为当特征距离在3Rm以外时，射流方向基本与竖直方向重合，即气泡受结构的影响较弱，基本可以忽略其对气泡运动的影响。

[1]Klaseboer E, Khoo B C, Hung K C. Dynamics of an oscillating bubble near a floating structure[J]. Journal of Fluids and Structures, 2005, 21:395-412.

[2]李健, 荣吉利. 水下爆炸圆柱壳塑性动态响应实验及数值计算[J]. 北京理工大学学报, 2008, 28(8): 559-562.

LI Jian, RONG Ji-li. Experimental and numerical investigations of dynamic response for cylindrical shell subjected to underwater explosion[J]. Transactions of Beijing Institute of Technology, 2008, 28(8): 559-562.

[3]张阿漫, 肖巍, 王诗平,等. 不同沙粒底面下气泡脉动特性实验研究[J]. 物理学报, 2013, 62(1): 014703, 1-8.

ZHANG A-man, XIAO Wei, WANG Shi-ping, et al. Experimental study of the interactions between a pulsating bubble and sand particles with different diameters[J]. Acta Physica Sinica, 2013, 62(1): 014703, 1-8.

[4]徐玲君, 陈刚, 邵建斌, 等. 不同直径气泡在静水中运动特性的研究[J]. 水动力学研究与进展，2012, 27(5): 582-588.

XU Ling-jun, CHEN Gang, SHAO Jian-bin, et al. Study on motion behavior of bubbles with different diameters in still water[J]. Chinese Journal of Hydrodynamics, 2012, 27(5): 582-588.

[5]Klaseboer E, Manica R, Chan Y C. BEM simulations of potential flow with viscous effects as applied to a rising bubble[J]. Engineering Analysis with Boundary Elements, 2011, 35: 489-494.

[6]张阿漫, 姚熊亮. 水下爆炸气泡与海底相互作用研究，振动与冲击, 2008,37(3), 92-98.

ZHANG A-man, YAO Xiong-liang. Interaction between underwater explosion bubble and seabed[J]. Journal of Vibration and Shock, 2008, 37(3), 92-98.

[7]张海鹏, 岳永威, 苏罗青等. 水下爆炸作用下舰船总纵强度模型实验方案研究[J]. 振动与冲击, 2012, 31(6): 175-180.

ZHANG Hai-peng, YUE Yong-wei, SU Luo-qing, et al. Model experimental scheme for longitudinal strength of a warship subjected to underwater explosion[J]. Journal of Vibration and Shock, 2012, 31(6): 175-180.

[8]Cox P E, Blake J R, et al. Bubble interactions near a free surface[J]. Engineering Analysis with Boundary Elements, 2004, 28:295-313.

[9]袁建红, 朱锡, 张振华. 水下爆炸近距气泡作用下船体箱型梁损伤特性数值模拟方法研究[J]. 振动与冲击, 2010, 29(4): 111-116.

YUAN Jian-hong, ZHU Xi, ZHANG Zhen-hua. Numerical simulation method for damage properties of a box-beam with action of underwater near-distance explosion bubble[J]. Journal of Vibration and Shock, 2010, 29(4): 111-116.

[10]李健, 荣吉利, 杨荣杰, 等.水中爆炸冲击波传播与气泡脉动的实验及数值模拟[J]. 兵工学报, 2008, 29(12): 1437-1443.

LI Jian, RONG Ji-li, YANG Rong-jie, et al. Experiment and numerical simulation of shock wave propagation and bubble impulse of underwater explosion[J]. Acta Armamentarii, 2008, 29(12): 1437-1443.

[11]黄晓明, 朱锡, 牟金磊等. 近距水下爆炸作用下箱形梁模型中垂破坏试验研究[J]. 振动与冲击, 2011, 30(2): 20-23.

HUANG Xiao-ming, ZHU Xi,MU Jin-lei, et al. Sagging damage test of box-beam models subjected to close range underwater explosion[J]. Journal of Vibration and Shock, 2011, 30(2): 20-23.

[12]Rayleigh J W. On t he pressure developed in a liquid during the collapse of a spherical cavity [J]. Philos Mag, 1917 , 34: 94-98.

[13]Price R S. Simulitude equations for explosives fired underwater[R]. NSWC: Technical Report NSWC TR, 1979, 80-299.

[14]Hunter K S. Global-shape-function models of an underwater explosion bubble[D]. Faculty of Graduate School of the University of Colorado, Department of Mechanical Engineering, 2001, 16.