水介质初始参数设置对水下爆炸载荷的影响

郑永辉; 魏继锋

doi:10.11883/bzycj-2021-0485

水介质初始参数设置对水下爆炸载荷的影响

doi: 10.11883/bzycj-2021-0485

郑永辉,
魏继锋^,

北京理工大学爆炸科学与技术国家重点实验室，北京 100081

基金项目: 国防科研基金（1020904）

详细信息

作者简介:
郑永辉（1996－　），男，硕士研究生，3038169560@qq.com

通讯作者:
魏继锋（1977－　），男，博士，副教授，weijifeng@bit.edu.cn

中图分类号: O389；TJ02
计量
- 文章访问数: 422
- HTML全文浏览量: 152
- PDF下载量: 114
- 被引次数: 6
出版历程
- 收稿日期: 2021-11-18
- 修回日期: 2021-12-25
- 网络出版日期: 2022-04-25
- 刊出日期: 2022-05-27

Effect of initial parameter setting of water on load characteristics of underwater explosion

ZHENG Yonghui,
WEI Jifeng^,

State Key Laboratory of Explosion Science and Technology, Beijing Institute of Technology, Beijing 100081, China

摘要

摘要: 针对数值计算中水介质初始参数设置对水下爆炸载荷特性的影响开展了深入分析。基于参考状态参数确定了水介质状态方程形式；从热力学角度分析了常用的两种初始参数设置方式，提出了一种按等温假设设置初始参数的方式，并对LS-DYNA中INITIAL_EOS_ALE关键字给出的参数设置结果进行了分析；采用LS-DYNA程序进行一维球形装药水下爆炸数值计算，分析了3种设置方式下爆炸载荷特性的差异，并与已有研究成果进行了对比。结果表明：当仅改变水介质内能项时，参数按等容过程变化，流场压力源于外界传热，与实际深水环境严重不符；INITIAL_EOS_ALE关键字给出的参数设置结果与仅改变水介质密度（等内能过程）接近，水温变化规律与真实环境不符；按等内能过程和等温过程设置初始参数时，水下爆炸载荷特性计算结果基本相同，与已有成果吻合；综合分析认为，按等温形式进行初始参数设置方式较优。研究成果可为水下爆炸尤其是深水爆炸数值仿真提供参考。
- 水下爆炸 /
- 初始压缩比 /
- 初始内能 /
- 冲击波 /
- 气泡
Abstract: The effect of the initial parameter setting of water medium on underwater explosion load characteristics in numerical simulation is studied. Firstly, based on the parameters under reference state, a kind of polynomial equation of state (EOS) is chosen as the EOS of water medium. Secondly, from the perspective of thermodynamics, two existing setting modes of initial parameters are analyzed, and a new setting mode following isothermal process is proposed. In addition, the results of initial compression ratio, initial internal energy and acoustic velocity of water under different depths given by INITIAL_EOS_ALE keyword in LS-DYNA program are investigated and compared with other three modes. Finally, by using the LS-DYNA program, numerical simulations of underwater explosion under different depths are conducted with a one-dimensional spherical charge model. The differences of shock wave load and bubble pulsation characteristics among the first three setting modes are discussed, which are also compared with previous studies. The results show that the parameters of water medium change with depths according to isochoric process if only the internal energy term of water medium is changed. Hence, it indicates that the pressures under different water depths are caused by thermal conduction from external environment, which is seriously inconsistent with the actual deep water condition. Initial parameters given by INITIAL_EOS_ALE keyword are close to the results obtained by only changing the density of water (i.e., following isometric energy process), but the changing laws of temperature for these two modes are both inconsistent with the real environment. When the parameters follow equal internal energy process or isothermal process, the calculated load characteristics are close to each other, which are consistent with existed studies. It is concluded, therefore, that initial parameter setting mode based on isothermal process is better than other three modes. This conclusion can provide an important reference to ensure the accuracy of underwater explosion numerical simulation, especially for deep water explosion.
- underwater explosion /
- initial compression ratio /
- initial internal energy /
- shock wave /
- bubble

HTML全文

随着数值计算方法的进步和计算能力的提升，数值仿真技术为解决水下爆炸机理分析、结构响应机制规律揭示等科学技术难题提供了有力的手段。在水下爆炸数值仿真研究中，研究人员已对水介质状态方程、人工黏性系数、网格密度等影响计算精度的因素开展了较为细致的探讨^[1-3]。

在进行水下爆炸问题计算时，需设定水介质的初始参数以设置流场初始压力。已有研究中，大多采用以下两种方式：一是认为水介质密度变化很小，可忽略不计，仅通过改变内能设置流场初始压力^[4-6]；二是认为流体静压力的增大导致水的密度增大，仅通过调整比容来设置初始压力^[7-8]。当采用前一种方式时，冲击波超压峰值随水深的增大略有减小^[6]，后一种方式下的变化趋势则相反^[7]。可以看出，初始参数的设置方式会带来水下爆炸载荷计算结果的差异，而对此的细致分析尚未见报道。

除直接设置水介质初始参数外，也可通过施加重力、采用程序提供的特定功能等方式实现流场压力初始化^[9-10]。这些方式简化了操作，但程序仍需提供水介质的初始参数使其具有一定的压力。以LS-DYNA提供的关键字INITIAL_EOS_ALE为例，其说明文档指出，程序将采用一种迭代方法根据用户提供的压力计算介质的初始内能和比容^[11]。简而言之，该类方法其本质上也是对水介质的初始参数进行设置。

水介质初始参数设置本质上是对水介质参数随水深变化的规律的反映。本文中将首先从常用的状态方程中选定符合参考状态时参数的水介质状态方程；然后从热力学角度分析现有两种设置方式对应的热力学过程，并给出第3种参数设置方式，同时以关键字INITIAL_EOS_ALE为例，对程序给出的初始化结果进行分析；随后采用LS-DYNA程序进行一维球形装药水下爆炸仿真研究，细致分析前三种初始参数设置方式对水下爆炸载荷特性的影响，并与已有研究成果进行对比；最后确定适当的初始参数设置方式。

1. 水介质状态方程

数值仿真软件中常见的水介质状态方程为Mie-Grüneisen状态方程和Polynomial状态方程。Mie-Grüneisen状态方程以冲击绝热线作为参考线，最终形式为：

$p = \left\{ \begin{array}{*{20}{l}} \dfrac{{{\rho _0}c_0^2\mu \left[ {1 + \left( {1 - \displaystyle\frac{{{\gamma _0}}}{2}} \right)\mu - \dfrac{a}{2}{\mu ^2}} \right]}}{{{{\left[ {1 + \mu -\displaystyle \sum\limits_{i = 1}^3 {{S_i}\dfrac{{{\mu ^i}}}{{{{\left( {1 + \mu } \right)}^{i - 1}}}}} } \right]}^2}}} + \left( {{\gamma _0} + a\mu } \right){e_V}&\mu \text{≥} {\text{0}} \\ {\rho _0}c_0^2\mu + \left( {{\gamma _0} + a\mu } \right){e_V}&\mu {{\text{＜} 0}} \end{array} \right.$

(1)

式中：c₀为参考状态下的声速；ρ₀为参考状态下的密度；μ为压缩比，μ = ρ/ρ₀−1；γ₀为Mie-Grüneisen系数；a为体积修正系数；S₁、S₂和S₃为试验拟合系数；e_V为体积内能增量。

从Mie-Grüneisen状态方程的推导过程可以看出，该方程忽略了参考压力p₀项，内能项e_V实际上是相对于参考状态下的内能增量^[12]。在参考状态下，μ = 0，e_V = 0，进而求得p = 0，这与真实参考压力p₀ = 101325 Pa不符。现有仿真计算中，大多将e_V设置为一个非零值，以使p = p₀^[8]。

Polynomial状态方程一般是对Mie-Grüneisen状态方程或者试验数据拟合得到^{[4, 8]}，常见形式为：

$p = \left\{ \begin{array}{*{20}{l}} {A_1}\mu + {A_2}{\mu ^2} + {A_3}{\mu ^3} + \left( {{B_0} + {B_1}\mu } \right){\rho _0}{e_M}&\mu \text{≥} {\text{0}} \\ {T_1}\mu + {T_2}{\mu ^2} + {B_0}{\rho _0}{e_M}&\mu \text{＜} 0 \end{array} \right.$

(2)

$p = \left\{ \begin{array}{*{20}{l}} {a_1}\mu + {a_2}{\mu ^2} + {a_3}{\mu ^3} + \left( {\displaystyle\sum\limits_{i = 0}^3 {{b_i}{u^i}} } \right){\rho _0}{e_M}&\mu \text{≥} 0 \\ {a_1}\mu + \left( {{b_0} + {b_1}\mu } \right){\rho _0}{e_M}&\mu \text{＜} 0 \end{array} \right.$

(3)

$p = \left\{ \begin{array}{*{20}{l}} {C_0} + {C_1}\mu + {C_2}{\mu ^2} + {C_3}{\mu ^3} + \left( {\displaystyle\sum\limits_{i = 4}^6 {{C_i}{\mu ^{i - 4}}} } \right){e_V}&\mu \text{≥} 0 \\ {C_0} + {C_1}\mu + {C_3}{\mu ^3} + \left( {{C_4} + {C_5}\mu } \right){e_V}&\mu \text{＜} 0 \hfill \\ \end{array} \right.$

(4)

式中：A₁、A₂、A₃、B₀、B₁、T₁和T₂均为状态方程系数，e_M为质量内能增量，a₁、a₂、a₃、b₀、b₁、b₂和b₃均为状态方程系数，C₀、C₁、C₂、C₃、C₄、C₅和C₆均为状态方程系数。

与Mie-Grüneisen状态方程类似，式（2）～（3）也存在与实际压力不符的问题。式（4）在C₀ = p₀时，可满足μ = 0，e_V = 0以及p = p₀。因此，宜选择Polynomial状态方程式（4）作为水介质状态方程，具体参数见表1^[13]。

表 1 水介质状态方程参数

Table 1. EOS parameters of water

C₀/Pa	C₁/GPa	C₂/GPa	C₃/GPa	C₄	C₅	C₆
101325	2.2	9.54	14.57	0.28	0.28	0

下载: 导出CSV

| 显示表格

2. 初始参数设置

2.1 设置方式及状态参数变化

根据水介质状态方程，设置初始压力时可改变参量μ和e_V。现今大多数仿真研究是通过改变e_V（方式Ⅰ）来进行流场初始压力设置^[4-6]。从热学角度看，方式Ⅰ中不同水深下的水介质参数按等容过程（dV=0）变化，此时外界做功为零，流场的内能增量全部由外界传导的热量提供（dE_i = dQ）。这表明采用方式Ⅰ时，流体静压力的增大单纯源于外界传热，与重力等外力做功等因素无关，这与实际深水环境以及加压模拟深水环境严重不符。

少数研究人员通过调整μ（方式Ⅱ）来进行初始压力设置^[7-8]。此时，内能变化为零（dE_i = 0），即等内能形式，外界对流体做功与外界向流体传导的热量之和为零（dE_i = dQ − pdV = 0），水介质被压缩后向外放出热量。为分析温度变化趋势，引入含温度T的内能表达式^[14]：

$\dfrac{{{e_V}}}{{{p_0}}} = 10\;160\dfrac{{\Delta T}}{{{T_0}}} + 0.63 \times {10^5} \times \dfrac{{\left( {1 - {{{\rho _0}}}/{\rho }} \right)\left( {0.71 - {{{\rho _0}}}/{\rho }} \right)}}{{{{\left( {{{{\rho _0}}}/{\rho }} \right)}^{{4 /3}}}}}\left( {1 - 2\dfrac{\rho }{{{\rho _0}}}{{\rm{e}}^{ - {{\left( {{\rho }/{{{\rho _0}}}} \right)}^2}}}} \right)$

(5)

式中：T₀ = 273 K。当0 ≤ μ ≤ 0.025时，T随μ的增大而逐渐增大（见图1）。在实际深水环境中，随着水深的增大，水温逐渐减小^[15]；在通过加压构造的模拟深水环境中，水温与环境温度基本相同，即温度可认为保持不变。因此，从温度变化趋势的角度来看，方式Ⅱ与实际情况有一定的差距。

图 1 不同设置方式下的μ、e_V和ΔT值

Figure 1. Values of μ, e_V and ΔT in different setting modes

下载: 全尺寸图片幻灯片

热力学研究结果表明，温度到处相同是重力场的热平衡条件^[16]。据此，本文中提出第3种初始参数设置方式，即假设水介质参数随水深按照等温过程（dT = 0，方式Ⅲ）变化。对于该方式，ΔT = 0，将式（5）代入式（4），可得等温过程对应的状态方程。显然，对于实际深水环境，方式Ⅲ只是一种粗糙的近似，不过理论上优于方式Ⅰ和方式Ⅱ；对于模拟深水爆炸试验，近似程度较高。

除上述设置方式外，数值仿真软件往往也具备根据重力和水深信息或者设定的流体静压力计算水介质初始参数的功能，只是用户往往无法直接获知其遵循的规律。现以关键字INITIAL_EOS_ALE为例，通过读取不同水深下的流场压力初始化结果来初步了解LS-DYNA程序中μ和e_V随水深的变化规律（方式Ⅳ），结果如图1所示，可以看出，方式Ⅳ对应的热力学过程接近方式Ⅱ。

图1中给出了4种设置方式在不同水深条件时的压缩比μ、内能增量e_V和温度变化ΔT的计算结果。方式Ⅱ和方式Ⅳ的计算结果几乎重合，二者均与方式Ⅲ较接近；方式Ⅰ的计算结果与其他3种相去甚远。采用方式Ⅱ、Ⅲ、Ⅳ时，μ随水深增加而增大；采用方式Ⅰ时μ则始终为零；在同一水深处，μ的值从小到大依次为：方式Ⅰ、方式Ⅳ、方式Ⅱ、方式Ⅲ。采用方式Ⅰ和方式Ⅳ时，e_V随水深增加而增大，采用方式Ⅱ时e_V始终为零，采用方式Ⅲ时e_V随水深的增加而减小；在同一水深处，e_V的值从小到大依次为：方式Ⅲ、方式Ⅱ、方式Ⅳ、方式Ⅰ。采用方式Ⅰ、方式Ⅱ和方式Ⅳ时，ΔT随水深增加而增大；采用方式Ⅲ时ΔT始终为零；在同一水深处，ΔT的值从小到大依次为：方式Ⅲ、方式Ⅱ、方式Ⅳ、方式Ⅰ。

2.2 流体可压缩性

流体的可压缩性决定流体内微弱扰动波的传播速度，即流体内的声音传播速度，因此常用声速c来表示流体的可压缩性。声速c的表达式为：

${c^2} = {\left( {\frac{{{\text{d}}p}}{{{\text{d}}\rho }}} \right)_S}$

(6)

式中：下标S表示等熵过程。

当初始状态为( μ_∞ , e_V∞ )时，等熵压缩状态下的压力可表示为：

$p = {C_0} + {C_1}\mu + {C_2}{\mu ^2} + {C_3}{\mu ^3} + \left( {{C_4} + {C_5}\mu + {C_6}{\mu ^2}} \right)\left( {{e_{V\infty }} + \int\nolimits_{{\mu _\infty }}^\mu {\frac{p}{{{{\left( {1 + \mu } \right)}^2}}}{\text{d}}\mu } } \right)$

(7)

于是，声速可表示为：

${c^2} = \frac{{{\text{d}}p}}{{{\text{d}}\rho }} = \frac{{{\text{d}}p}}{{{\text{d}}\mu }}\frac{{{\text{d}}\mu }}{{{\text{d}}\rho }} = \frac{{{C_1} + 2{C_2}\mu + 3{C_3}{\mu ^2}}}{{{\rho _0}}} + \frac{{{C_5} + 2{C_6}\mu }}{{{\rho _0}}}\left( {{e_{V\infty }} + \int\nolimits_{{\mu _\infty }}^\mu {\frac{p}{{{{\left( {1 + \mu } \right)}^2}}}{\text{d}}\mu } } \right) + \frac{{{C_4} + {C_5}\mu + {C_6}{\mu ^2}}}{{{\rho _{0}}}}\frac{p}{{{{\left( {1 + \mu } \right)}^2}}}$

(8)

当μ = μ_∞时，声速为：

${c^2} = \frac{{{C_1} + 2{C_2}{\mu _\infty } + 3{C_3}\mu _\infty ^2}}{{{\rho _0}}} + \frac{{{C_5} + 2{C_6}{\mu _\infty }}}{{{\rho _0}}}{e_{V\infty }} + \frac{{{C_4} + {C_5}{\mu _\infty } + {C_6}\mu _\infty ^2}}{{{\rho _0}}}\frac{{{p_\infty }}}{{{{\left( {1 + {\mu _\infty }} \right)}^2}}}$

(9)

图2中给出了4种设置方式在不同水深条件时的声速。从图2中可以看出，随着水深的增大，4种方式下的水介质声速逐渐增大，即可压缩性均出现减弱，但减弱幅度不同，由小到大的排序为：方式Ⅰ、方式Ⅳ、方式Ⅱ、方式Ⅲ。可见，声速的变化与压缩比的变化情形一致。

图 2 不同设置方式下的声速变化

Figure 2. Acoustic velocity in different modes

下载: 全尺寸图片幻灯片

3. 初始参数设置对水下爆炸载荷的影响

3.1 数值仿真模型

建立了一维球对称水下爆炸数值仿真模型，并应用LS-DYNA数值仿真软件进行计算。其中，炸药采用MAT_HIGH_EXPLOSIVE_BURN材料模型以及JWL状态方程，参数取自文献[17]。水介质的状态方程采用式（4），参数在表1中列出。

炸药选用1 kg TNT，水深范围0 ～ 5 km，由于方式Ⅱ和方式Ⅳ对应的初始参数极为接近，因此只对比前3种方式，各自对应的μ和e_V值取自图1。为保证计算精度，网格尺寸为装药半径的1/100；为减小边界反射冲击波对气泡脉动的影响，水域半径取237.2 m，大于水深为0 m时声波在一个脉动周期内传播距离的1/2。

3.2 冲击波载荷特性

选取水深为5 km，爆距R = 55R₀（R₀为装药半径），计算得到了3种方式下的冲击波压力-时间曲线，如图3所示，对应的冲击波到达时间t_a、超压峰值ΔP_m和正压持续时间t_c列于表2中。可以看出，3条曲线并不完全重合，除方式Ⅰ外，其他2种方式的仿真结果相差较小；方式Ⅲ、Ⅱ、Ⅰ的t_a及t_c值依次减小，而ΔP_m依次增大。该规律与流体可压缩性结果基本一致，即流体可压缩性越弱，冲击波到达时间越短，冲击波超压峰值越大，正压持续时间越短。

图 3 H=5 km及R=55R₀时3种方式下的冲击波压力-时间曲线

Figure 3. Shock wave pressure-time curves in three modes when H=5 km and R=55R₀

下载: 全尺寸图片幻灯片

表 2

$H=5\;{\rm {km}}$ 及

$R=55R_0$ 时3种方式下的

$t_{\rm a}$ 、

$\Delta P_{\rm m}$ 及

$t_{\rm c}$

Table 2. Values of

$t_{\rm a}$ ,

$\Delta P_{\rm m}$ and

$t_{\rm c}$ in three modes when

$H=5\;{\rm {km}}$ and

$R=55R_0$

方式	t_a/ms	ΔP_m/MPa	t_c/ms
Ⅰ	1.7745	12.080	0.2676
Ⅱ	1.6578	13.603	0.2587
Ⅲ	1.6518	13.692	0.2581

下载: 导出CSV

| 显示表格

图4中给出了R = 55R₀时3种设置方式在不同水深下获得的冲击波超压峰值、冲量I和能流密度e_s，其中，I和e_s分别由下式给出；

图 4 R=55R₀时不同水深下3种方式所对应的冲击波载荷参数

Figure 4. Values of shock wave load in three modes when R=55R₀

下载: 全尺寸图片幻灯片

$I = \int\nolimits_0^{{t_{\text{c}}}} {\left( {p - {p_\infty }} \right){\text{d}}t}$

(10)

${e_{\text{s}}} = \frac{1}{{{\rho _\infty }{c_\infty }}}\int\nolimits_0^{{t_{\text{c}}}} {{{\left( {p - {p_\infty }} \right)}^2}{\text{d}}t}$

(11)

式中：e_s为冲击波能流密度，ρ_∞、c_∞分别为未扰动流体的密度和声速，为便于比较，统一取ρ_∞ = ρ₀=1000 kg/m³，c_∞ =c₀ = 1483.2 m/s。

表3中给出了冲击波载荷参数在5 km水深处相对于0 m的变化幅度。可以看出，随着水深的增大，采用方式Ⅰ时ΔP_m逐渐减小，其他2种方式则逐渐增大，且变化幅度略大于方式Ⅰ。对于I和e_s，3种方式给出的仿真结果均随水深的增大而减小，方式Ⅱ和Ⅲ的变化幅度略小于方式Ⅰ。

表 3 冲击波载荷在5 km水深处相对于0 m的变化幅度

Table 3. Changing amplitudes of shock wave load when depth changes from 0 m to 5 km

方式	变化幅度/%
方式	ΔP_m	I	e_s
Ⅰ	−1.202	−74.304	−32.896
Ⅱ	11.255	−73.227	−21.877
Ⅲ	11.981	−73.171	−21.247

下载: 导出CSV

| 显示表格

3.3 气泡脉动特性

分析H=5 km时3种设置方式下气泡的脉动曲线，如图5所示，相应的气泡最大半径R_max和第1次脉动周期T列于表4。可以看出，除方式Ⅰ外，其余两条曲线基本重合；方式Ⅲ、方式Ⅱ、方式Ⅰ对应的R_max和T依次增大。显然，该规律与水介质压缩比的变化一致。由于气泡脉动与水介质的惯性密切相关，当体积不变时，水介质密度的增大导致其惯性增大，在密度变化较小时惯性增幅有限，因此，3种方式下气泡脉动曲线的差别并不明显。

图 5 H=5 km时3种方式下气泡的半径-时间曲线

Figure 5. Bubble radius-time curves in three modes when H=5 km

下载: 全尺寸图片幻灯片

表 4

$H=5\;{\rm{km}}$ 时3种方式下对应的

$R_{\max}$ 和

$T$

Table 4. Values of

$R_{\max}$ and

$T$ in three modes when

$H=5\;{\rm{km}}$

方式	R_max/mm	T/ms
Ⅰ	178.722	1.7129
Ⅱ	180.618	1.7445
Ⅲ	180.731	1.7460

下载: 导出CSV

| 显示表格

4. 对比分析

4.1 冲击波载荷对比分析

自20世纪60年代起，研究人员开展了较多的水深对冲击波载荷特性影响的研究。Baum等^[18]在模拟深水爆炸容器中进行了深水爆炸试验，当1 < R/R₀ < 7时，ΔP_m在0 ～ 4 km的模拟水深范围内基本保持不变。Vanzant等^[19]给出了14.75 g Pentolite炸药在更大爆距范围内的试验结果，本文对其试验数据进行了线性拟合，结果表明，ΔP_m随水深的增大而略有增大，且在较大爆距处的变化更明显，如图6所示。钟帅等^[20]开展了0～150 m模拟水深水下爆炸试验，指出ΔP_m也随水深增大而略有增加。Slifko^[21]、Xiao等^[22]也进行了类似研究，不过试验在大洋中进行，爆源和测点并不位于同一深度，而且爆距较大，与本文分析工况不符。

图 6 Vanzant等^[19]的试验结果及线性拟合曲线

Figure 6. Test results from Vanzant et al^[19] and linear fitting results

下载: 全尺寸图片幻灯片

理论上，Baum等^[18]从相似性和量纲理论出发，在浅水经验公式的基础上，给出了考虑流体静压力时ΔP_m的计算式：

$\Delta {P_{\text{m}}} = A{\left( {\frac{{{W^{{1 \mathord{\left/ {\vphantom {1 3}} \right. } 3}}}}}{R}} \right)^\alpha }{\left( {\frac{{{p_\infty }}}{B} + 1} \right)^{{{\left( {N - \alpha } \right)} \mathord{\left/ {\vphantom {{\left({N - \alpha }\right)} N}} \right. } N}}}$

(12)

式中：A和α为常数，当6 ≤ R/R₀ < 12时，A = 44.1 MPa，α = 1.5；当12 ≤ R/R₀ < 240时，A = 52.4 MPa，α = 1.13；B = 304.9 MPa为Tait状态方程中的常数；N为几何指数，对于球面冲击波，N = 3。

式（12）表明，ΔP_m随水深的增大而缓慢增大；同时随着爆距的增大，α逐渐减小，流体静压力对ΔP_m的影响逐渐增大。

从前述研究中可以看到，随着水深的增大，采用方式Ⅱ和方式Ⅲ时ΔP_m会逐渐增大。从图7中可以看到，5 km水深条件下当7R₀ ≤ R≤ 65R₀时，随着爆距的增大，水深对ΔP_m的影响也在逐渐增大。图8为H = 5 km时3种方式下ΔP_m计算结果相对于式（12）的误差。综合来看，采用方式Ⅱ和方式Ⅲ获得的计算结果与试验值及理论值均较为吻合，方式Ⅲ的相对误差更小。

图 7 ΔP_m在5 km水深处相对于在0 m处的变化幅度

Figure 7. Changing amplitudes of ΔP_m when depth changes from 0 m to 5 km

下载: 全尺寸图片幻灯片

图 8 H=5 km时ΔP_m相对于计算式(8)的误差

Figure 8. Relative errors of ΔP_m between simulation results and calculation formulas when H=5 km

下载: 全尺寸图片幻灯片

4.2 气泡脉动参数对比分析

Cole^[23]基于不可压缩理论模型推导了R_max和T与水深和装药量的关系，Swift等^[24]基于试验数据给出了计算气泡脉动参数的经验公式：

${R_{{\text{max}}}} = 3.36\frac{{{W^{{1 \mathord{\left/ {\vphantom {1 3}} \right. } 3}}}}}{{{{\left( {H + 10.3} \right)}^{{1 \mathord{\left/ {\vphantom {1 3}} \right. } 3}}}}}$

(13)

$T = 2.11\frac{{{W^{{1 \mathord{\left/ {\vphantom {1 3}} \right. } 3}}}}}{{{{\left( {H + 10.3} \right)}^{{5 \mathord{\left/ {\vphantom {5 6}} \right. } 6}}}}}$

(14)

Baum等^[18]、Slifko^[21]和梁浩哲等^[25]借助试验进一步研究了气泡脉动参数随水深的变化规律，其中，Baum等^[18]验证了R_max ∝ p_∞^−1/3及T ∝ p_∞^−5/6在大水深（H > 1 km）下的适用性。

图9中给出了3种方式下气泡脉动参数仿真结果相对于经验值的误差。可以看出，在同一水深下，方式Ⅱ、Ⅲ误差较小，与经验值更接近，且方式Ⅲ的相对误差最小。

图 9 气泡脉动参数相对于经验值的误差

Figure 9. Relative errors of R_max and T between simulation and empirical formula

下载: 全尺寸图片幻灯片

4.3 设置方式合理性探讨

从热力学角度看，2.1节中已指出，采用方式Ⅰ时，流体静压力的升高单纯由外界传热引起，这与真实情形不符；采用方式Ⅱ和方式Ⅳ时，温度随水深的增大而略有增大，与实际情况有一定差距。采用方式Ⅲ时，水介质温度不随水深变化，与模拟深水试验条件吻合，与实际深水环境的吻合性较方式Ⅰ、Ⅱ和Ⅳ要好。从与已有成果的对比分析看，方式Ⅱ、Ⅲ的结果与已有成果吻合较好；相较而言，采用方式Ⅲ时载荷值的相对误差最小。综上分析认为，初始参数设置宜选用方式Ⅲ，方式Ⅱ和方式Ⅳ次之，在讨论水深影响时不建议选用方式Ⅰ。

5. 结　论

分析了4种水介质初始参数设置方式对应的热力学过程，并采用有限元数值仿真方法，深入探讨了初始参数设置方式对水下爆炸载荷特性的影响，并与相关试验与理论研究成果进行了对比分析，主要结论如下。

（1）根据参考状态下水介质参数，选择式（4）形式的状态方程作为水介质状态方程；其他常用状态方程均忽略了参考压力项，无法满足参考状态下的水介质参数。

（2）仅改变水介质内能增量e_V，流体压力源于外界传热，与实际环境不符；仅改变压缩比μ时，水介质温度随水深的增大而略有增加；LS-DYNA中INITIAL_EOS_ALE关键字给出的初始参数计算结果与等内能假定非常接近；与上述3种方式相比，本文提出的等温假定更符合真实深水环境以及模拟深水试验条件。

（3）按等内能假定和等温假定设置初始压力时，水下爆炸载荷特性结果接近，与相关试验与理论研究成果吻合良好；综合分析实际深水环境以及模拟深水试验条件，在数值仿真计算的初始压力设置时，宜选用本文提出的水介质参数等温变化假定。

图 1 不同设置方式下的μ、e_V和ΔT值

Figure 1. Values of μ, e_V and ΔT in different setting modes

下载: 全尺寸图片幻灯片

图 2 不同设置方式下的声速变化

Figure 2. Acoustic velocity in different modes

下载: 全尺寸图片幻灯片

图 3 H=5 km及R=55R₀时3种方式下的冲击波压力-时间曲线

Figure 3. Shock wave pressure-time curves in three modes when H=5 km and R=55R₀

下载: 全尺寸图片幻灯片

图 4 R=55R₀时不同水深下3种方式所对应的冲击波载荷参数

Figure 4. Values of shock wave load in three modes when R=55R₀

下载: 全尺寸图片幻灯片

图 5 H=5 km时3种方式下气泡的半径-时间曲线

Figure 5. Bubble radius-time curves in three modes when H=5 km

下载: 全尺寸图片幻灯片

图 6 Vanzant等^[19]的试验结果及线性拟合曲线

Figure 6. Test results from Vanzant et al^[19] and linear fitting results

下载: 全尺寸图片幻灯片

图 7 ΔP_m在5 km水深处相对于在0 m处的变化幅度

Figure 7. Changing amplitudes of ΔP_m when depth changes from 0 m to 5 km

下载: 全尺寸图片幻灯片

图 8 H=5 km时ΔP_m相对于计算式(8)的误差

Figure 8. Relative errors of ΔP_m between simulation results and calculation formulas when H=5 km

下载: 全尺寸图片幻灯片

图 9 气泡脉动参数相对于经验值的误差

Figure 9. Relative errors of R_max and T between simulation and empirical formula

下载: 全尺寸图片幻灯片

表 1 水介质状态方程参数

Table 1. EOS parameters of water

C₀/Pa C₁/GPa C₂/GPa C₃/GPa C₄ C₅ C₆

101325 2.2 9.54 14.57 0.28 0.28 0

下载: 导出CSV

表 2 $H=5\;{\rm {km}}$ 及 $R=55R_0$ 时3种方式下的 $t_{\rm a}$ 、 $\Delta P_{\rm m}$ 及 $t_{\rm c}$

Table 2. Values of $t_{\rm a}$ , $\Delta P_{\rm m}$ and $t_{\rm c}$ in three modes when $H=5\;{\rm {km}}$ and $R=55R_0$

方式 t_a/ms ΔP_m/MPa t_c/ms

Ⅰ 1.7745 12.080 0.2676
Ⅱ 1.6578 13.603 0.2587
Ⅲ 1.6518 13.692 0.2581

下载: 导出CSV

表 3 冲击波载荷在5 km水深处相对于0 m的变化幅度

Table 3. Changing amplitudes of shock wave load when depth changes from 0 m to 5 km

方式变化幅度/%
ΔP_m I e_s

Ⅰ −1.202 −74.304 −32.896
Ⅱ 11.255 −73.227 −21.877
Ⅲ 11.981 −73.171 −21.247

下载: 导出CSV

表 4 $H=5\;{\rm{km}}$ 时3种方式下对应的 $R_{\max}$ 和 $T$

Table 4. Values of $R_{\max}$ and $T$ in three modes when $H=5\;{\rm{km}}$

方式 R_max/mm T/ms

Ⅰ 178.722 1.7129
Ⅱ 180.618 1.7445
Ⅲ 180.731 1.7460

下载: 导出CSV

参考文献(25)

[1]	杨坤, 陈朗, 伍俊英, 等. 计算网格与人工粘性系数对炸药水中爆炸数值模拟计算的影响分析 [J]. 兵工学报, 2014, 35(S2): 237–243. YANG K, CHEN L, WU J Y, et al. The effects of computing grid and artificial viscosity coefficient on underwater explosion numerical simulation [J]. Acta Armamentarii, 2014, 35(S2): 237–243.
[2]	WANG G H, WANG Y X, LU W B, et al. On the determination of the mesh size for numerical simulations of shock wave propagation in near field underwater explosion [J]. Applied Ocean Research, 2016, 59: 1–9. DOI: 10.1016/j.apor.2016.05.011.
[3]	胡亮亮, 黄瑞源, 李世超, 等. 水下爆炸冲击波数值仿真研究 [J]. 高压物理学报, 2020, 34(1): 015102. DOI: 10.11858/gywlxb.20190773. HU L L, HUANG R Y, LI S C, et al. Shock wave simulation of underwater explosion [J]. Chinese Journal of High Pressure Physics, 2020, 34(1): 015102. DOI: 10.11858/gywlxb.20190773.
[4]	CHISUM J E, SHIN Y S. Explosion gas bubbles near simple boundaries [J]. Shock and Vibration, 1997, 4(1): 11–25. DOI: 10.3233/SAV-1997-4102.
[5]	詹发民, 姜涛, 任佳宁, 等. 深水爆炸载荷数值仿真研究 [J]. 工程爆破, 2013, 19(6): 9–12. DOI: 10.3969/j.issn.1006-7051.2013.06.003. ZHAN F M, JIANG T, REN J N, et al. Numerical simulation of explosion load in deep-water [J]. Engineering Blasting, 2013, 19(6): 9–12. DOI: 10.3969/j.issn.1006-7051.2013.06.003.
[6]	张亮. 典型装药深水爆炸特性研究 [D]. 北京: 北京理工大学, 2020: 56–57.
[7]	BRETT J M. Numerical modelling of shock wave and pressure pulse generation by underwater explosions: AD010558 [R]. Australia: Aeronautical and Maritime Research Laboratory, 1998.
[8]	SHIN Y S, LEE M, LAM K Y, et al. Modeling mitigation effects of watershield on shock waves [J]. Shock and Vibration, 1998, 5(4): 225–234. DOI: 10.1155/1998/782032.
[9]	鲁忠宝, 南长江, 步相东, 等. 不同水深爆炸气泡运动特性仿真 [J]. 鱼雷技术, 2009, 17(5): 15–18. DOI: 10.3969/j.issn.1673-1948.2009.05.004. LU Z B, NAN C J, BU X D, et al. Simulation of moving features of underwater explosion bubble in different water depths [J]. Torpedo Technology, 2009, 17(5): 15–18. DOI: 10.3969/j.issn.1673-1948.2009.05.004.
[10]	董琪, 韦灼彬, 唐廷, 等. 爆炸深度对浅水爆炸气泡脉动的影响 [J]. 高压物理学报, 2018, 32(2): 201420. DOI: 10.11858/gywlxb.20170580. DONG Q, WEI Z B, TANG T, et al. Influence of explosion depth on bubble pulsation in shallow water explosion [J]. Chinese Journal of High Pressure Physics, 2018, 32(2): 201420. DOI: 10.11858/gywlxb.20170580.
[11]	LSTC. LS-DYNA keyword user’s manual Volume 1 [M]. LSTC, 2019: 2259-2260.
[12]	李晓杰, 张程娇, 王小红, 等. 水的状态方程对水下爆炸影响的研究 [J]. 工程力学, 2014, 31(8): 46–52. DOI: 10.6052/j.issn.1000-4750.2013.03.0180. LI X J, ZHANG C J, WANG X H, et al. Numerical study on the effect of equations of state of water on underwater explosions [J]. Engineering Mechanics, 2014, 31(8): 46–52. DOI: 10.6052/j.issn.1000-4750.2013.03.0180.
[13]	PONSONNAILLE J P, EDJTEMAI N. Numerical simulation and analysis of fluid surge due to the rupture of a storage tank [C]//First European LS-DYNA Conference. England, 1997:1–5.
[14]	ОРЛЕНКО Л П. 爆炸物理学(下册) [M]. 孙承纬, 译. 北京: 科学出版社, 2011: 1311–1312. ОРЛЕНКО Л П. Explosion physics[M]. SUN C W, translate. Beijing: Science Press, 2011: 1311–1312.
[15]	霍文娟, 韩震. 印度洋北部海域垂直剖面温度结构的参数模型 [J]. 海洋环境科学, 2013, 32(3): 368–372,380. HUO W J, HAN Z. Parametric model of seawater temperature profile in the Northern Indian Ocean [J]. Marine Environmental Science, 2013, 32(3): 368–372,380.
[16]	王竹溪. 热力学 [M]. 北京: 北京大学出版社, 2011: 274–276.
[17]	LEE E L, HORNIG H C, KURY J W, et al. Adiabatic expansion of high explosive detonation products: UCRL-50422 [R]. Livermore, California, USA: Lawrence Radiation Laboratory, 1968.
[18]	BAUM F A, SANASARYAN N S. Effect of hydrostatic pressure on the parameters of an underwater explosion [J]. Combustion, Explosion and Shock Waves, 1965, 1(4): 33–38. DOI: 10.1007/BF00748810.
[19]	VANZANT B W, DEHART R C. Effect of hydrostatic pressure on shock waves from underwater explosions [J]. Journal of Applied Physics, 1965, 36(10): 3116–3117. DOI: 10.1063/1.1702934.
[20]	钟帅. 模拟深水爆炸装药输出能量的研究 [D]. 安徽淮南: 安徽理工大学, 2007: 44–45. ZHONG S. Study on energy output of explosive charge under simulated deep-water explosion[D]. Huainan, Anhui: Anhui University of Science and Technology, 2007: 44–45.
[21]	SLIFKO J F. Pressure-pulse characteristics of deep explosions as functions of depth and range: AD661804 [R]. Maryland: U. S. Naval Ordnance Lab, 1967.
[22]	XIAO P, YANG K D. Experimental results for peak pressure and sound exposure level in deep-sea explosions [J]. Acoustics Australia, 2015, 43(2): 175–178. DOI: 10.1007/s40857-015-0020-9.
[23]	COLE R H. Underwater explosion [M]. Princeton, USA: Princeton University Press, 1948: 270–276.
[24]	SWIFT E, DECIUS J C. Measurement of bubble pulse phenomena Ⅲ: Radius and period studies in underwater explosion research [R]. Washington D C: Office of Naval Research, 1950.
[25]	LIANG H Z, ZHANG Q M, LONG R R, et al. Pulsation behavior of a bubble generated by a deep underwater explosion [J]. AIP Advances, 2019, 9(2): 025108. DOI: 10.1063/1.5086361.

施引文献

期刊类型引用(5)

1.	周刚，孔阳，崔洋洋，钱新明，傅砺烨，张琦. 城市地下排水管道中燃气爆炸及气-液两相耦合作用规律. 爆炸与冲击. 2024(03): 90-104 . 本站查看
2.	张勇，肖正明，段浩，伍星，卢敏，王浩. 水下中远场爆炸冲击波作用下航行体表面动态响应分析. 兵工学报. 2024(07): 2341-2350 . 百度学术
3.	严侃. 典型水下航行器爆炸毁伤动力学特性研究. 水下无人系统学报. 2024(06): 1108-1116 . 百度学术
4.	严侃. 典型水下航行器爆炸毁伤动力学特性研究. 水下无人系统学报. 2024(06): 1108-1116 . 百度学术
5.	闫秋实，常松. 水下爆炸三维数值模拟特征参量敏感性分析. 北京工业大学学报. 2023(10): 1099-1108 . 百度学术