球形装药动态爆炸冲击波超压场计算模型*

聂 源,蒋建伟,李 梅

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

球形装药动态爆炸冲击波超压场计算模型*

聂 源,蒋建伟,李 梅

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

为获得球形装药动态爆炸冲击波超压场计算模型，对静态爆炸冲击波超压Baker计算公式加入修正因子进行修正，并建立了构造包含装药运动速度、对比距离和方位角的修正因子函数的方法。为获得修正因子的函数表达式，采用高精度显式欧拉流体动力学软件SPEED针对具有典型运动速度的球形装药空中爆炸过程进行了数值模拟，得到了沿装药不同对比距离和方位角处的动态爆炸冲击波超压峰值。在对数值模拟结果处理的基础上，经过数据拟合获得了动态爆炸冲击波超压场计算模型。校验结果表明，该模型能较准确描述动态爆炸冲击波超压分布，具有普适性。

冲击波超压；动态爆炸；装药运动速度；超压计算模型

爆炸冲击波超压是衡量弹药战斗部爆炸威力的重要指标之一，由于装药随载体运动爆炸产生的冲击波超压测试较复杂，目前主要以地面静爆实验为主[1]。动态与静态爆炸冲击波场存在明显差异[1]，静爆实验无法准确测量战斗部的实际威力，因此研究装药在动态爆炸条件下冲击波超压场及其计算模型非常必要。

目前，关于静态爆炸的冲击波超压计算模型较多，如根据爆炸相似律建立的Henrych公式[2]、Baker公式[3]和Sadovskyi公式[4]等，而动态爆炸的超压计算模型却不多。在文献[5]中，根据能量相似原理，将炸药随载体运动的动能等效为附加装药量，并与实际装药量累加后代入静爆超压计算模型，该模型获得的冲击波超压值与方位角θ(运动方向为θ=0)无关；张光莹等[1]对具有典型运动速度的装药在正向(θ=0)和反向(θ=π)处的冲击波超压峰值进行了数值计算，得到动态爆炸冲击波超压峰值具有方向性，但并未得出其计算模型。

本文中，在对数值模拟结果数据分析的基础上，通过建立装药运动速度、方位角和对比距离的修正因子，得到基于静态条件冲击波超压Baker计算模型的球形装药动态爆炸冲击波超压场计算模型，并采用数值模拟方法对不同尺寸的球形装药动态爆炸冲击波超压进行验证，试图表明该模型具有普适性。

1 模型建立

记Δp0为静态爆炸条件下(装药运动速度v=0时)，A点的冲击波超压，MPa；Δpv,Z,θ为动态爆炸条件下(装药运动速度v≠0时)，A点的冲击波超压，MPa。

定义修正因子δ为动态爆炸冲击波超压相对静态爆炸冲击波超压的增量，即：

(1)

用修正因子δ计算动态爆炸过程A点处的冲击波超压峰值Δpv,Z,θ，变换式(1)得：

Δpv,Z,θ=(1+δ)Δp0

(2)

静态爆炸冲击波超压Δp0采用Baker公式进行计算[3]：

(3)

修正因子δ是与装药运动速度v、方位角θ和对比距离Z等3个因子有关的量，即δ=f(v,Z,θ)，假设3个因子v、Z和θ相互独立[6]，则δ可表示为：

δ=f(v)f(Z)f(θ)

(4)

f(v)为在特定的方位角θ*和对比距离Z*处(θ*和Z*均为任取)、不同装药运动速度v时的δ值，也就是：

f(v)=δv,Z*,θ*

(5)

f(Z) 为在特定方位角θ*，取一系列装药运动速度v时、不同对比距离Z处的δ值与Z=Z*时δ值之比，也就是：

(6)

式中：v=[v1,v2,…,vn]，v1、v2、 …、vn为不同的装药运动速度。

f(θ)为取一系列装药运动速度v和对比距离Z时，不同方位角θ处的δ值与θ=θ*时的δ值的比值，也就是：

(7)

式中：Z=[Z1,Z2,…,Zn]，Z1、Z2、 …、Zn为不同的对比距离。

以上为基于Baker公式的动态爆炸冲击波超压场计算模型的建立方法。该计算模型得到的动态爆炸冲击波超压，可以作为入射超压用于计算反射超压，并且同样适用于不同炸药类型、装药形状的修正模型。

2 修正因子的获取

式(5)～(7)为构造修正因子函数的

方法。为获得具体表达式，采用数值模拟方法对具有不同运动速度的标准球形装药在空气中动态爆炸进行数值模拟分析，根据式(5)～(7)采用数据拟合方法得到修正因子的函数表达式。

数值模拟采用高精度显式欧拉流体动力学软件SPEED(shock physical explicit eulerian dynamic)，图2为1 kg的TNT球形装药在空气中动态爆炸的二维轴对称数值模型，其中装药运动速度v=[400,800,1 200,1 600] m/s。

利用SPEED的自适应网格技术[7]，并设置空气域边界条件为透射边界。当冲击波到达边界网格时，边界单元内材料的速度达到一定范围，触发网格按一定放大系数向四周扩张[7]，以确保爆轰产物在欧拉区域内。

沿装药典型方位角θ=[0,π/6,π/4,π/3,π/2,2π/3,3π/4,5π/6,π]和Z=[0.25,0.5,0.75,1,2,3,4,5,6] m/kg1/3处设置了观测点。

TNT炸药采用JWL状态方程，其一般压力形式[8]为：

(8)

TNT炸药的JWL状态方程参数分别为：ρ0=1.63 g/cm3,A=373.8 GPa,B=3.747 GPa,R1=4.15,R2=0.90,ω=0.35,E0=6.00 GPa。

空气采用理想气体状态方程[6]：

p=(γ-1)ρe

(9)

式中：p为气体压力；γ为多方气体指数，γ=1.4；ρ为密度；e为内能，e=206.8 J/g。

对上述工况进行计算，图3为典型工况在典型时刻的压力云图。从图3可见，动态爆炸冲击波超压场与静态爆炸的不同，动态爆炸冲击波超压随方位角的增大而减小。

提取各工况中各观测点的超压曲线，根据式(1)、(5)、(6)和(7)进行数据处理。为便于分析，取θ*=0和Z*=1.5 m/kg1/3，并对参量做归一化处理。

图4为根据式(5)计算得到的f(v)及其拟合曲线。从图4可见，f(v)与装药运动速度v呈现线性关系，经最小二乘法拟合得到：

f(v)=0.53v/c0

(10)

式中：c0为声速，取340 m/s。

图5为根据式(6)计算得到的不同v时的f(Z)曲线和平均值及拟合曲线。从图5(a)可知，不同v时的f(Z)曲线形状相似，随Z的增大均呈现先增长后递减最终趋于零的规律，且在Z=1 m/kg1/3处，f(Z)达到极大值。取不同v时f(Z)曲线的平均值，得到图5(b)，经分段拟合得到：

(11)

图6为根据式(7)计算得到的不同v和Z时的f(θ)曲线和平均值及拟合曲线。从图6(a)可见，不同v和Z时的f(θ)曲线形状相似，f(θ)均随θ的增大而非线性减小。取不同v和Z时f(θ)曲线的平均值，得到图6(b)，经拟合得到f(θ)与θ符合三角函数的变化趋势：

(12)

综上，得到了考虑动态爆炸条件下冲击波超压修正因子δ的计算公式(4)、(10)、(11)和(12)，将该系数应用于式(2)，则得到基于Baker公式的动态爆炸冲击波超压场计算模型：

(13)

3 模型校验

式(13)是以典型药量为参考得到的动态爆炸冲击波超压场计算模型，为检验该模型能否计算不同药量的动态爆炸冲击波超压，同样采用SPEED软件建立两种工况(工况1,M=15 kg,v=900 m/s，工况2,M=230 kg,v=300 m/s)的空中动态爆炸模型，经数值模拟得到典型方位角和对比距离处的动态爆炸冲击波超压峰值。

图7为在工况1、2条件下，数值模拟和计算模型式(13)得到的在不同方位角和对比距离处的动态爆炸冲击波超压峰值对比图。从图中可见，该修正因子修正后的动态爆炸冲击波超压与数值模拟结果基本吻合，所以该计算模型能正确反映动态爆炸冲击波超压特性，具有普适性。

4 结 论

(1)引入包含装药运动速度v、对比距离Z和方位角θ的修正因子，对Baker公式进行修正，并建立了构造修正因子函数的方法，得到球形装药动态爆炸冲击波超压场计算模型，即式(2)、(4)～(7)。

(2)基于数值模拟结果的数据处理，得到了球形装药动态爆炸冲击波超压场计算模型表达式(13)。对该模型的校验结果表明，计算模型能准确计算不同药量、不同装药运动速度的球形装药动态爆炸冲击波超压，具有普适性。

[1] 张光莹，周旭，黄咏政，等．动爆冲击波特性分析方法研究[C]∥第四届全国计算爆炸力学会议论文集．2008：282-287.

[2] Henrych J, Abrahamson G R. The dynamics of explosion and its use[M]. New York: Elsevier Scientific Pub Co, 1979:218.

[3] Baker W E. Explosions in air[M]. Austin: University of Texas Press, 1974:6-10.

[4] Sadovskyi M A. Mechanical action of air shock waves of explosion, based on experimental data[M]. Moscow: Izd Akad Nauk SSSR, 1952:1-2.

[5] 北京工业学院八系．爆炸及其作用[M]．北京：工业出版社，1978.

[6] Izadifard R A, Foroutan M. Blastwave parameters assessment at different altitude using numerical simulation[J]. Turkish Journal of Engineering & Enviromental Sciences, 2010,34(1):25-41.

[7] NUMERICS GmbH. SPPED user’s Manual[Z]. 2012.

[8] Lee E, Finger M, Collins W. JWL equation of state coeffients for high explosives: UCID-16189[R]. 1973.

Abstract: In order to establish a calculation model for the dynamic blast overpressure field, a correction factor was introduced into the Baker formula, a model for calculating the peak overpressure only in a static blast. The method to obtain the correction factor containing the moving velocity, the azimuth and the scaled distance was established. For getting the function of the correction factor, spherical charge models with the typical moving velocity were established using the shock physical explicit Eulerian dynamic (SPEED) to simulate the dynamic blast process in the air. The peak overpressure in the typical azimuth and scaled distance was obtained. Based on the numerical results, a new calculation model was constructed using the data fitted. The result of the numerical simulation indicates that the corrected formula is a universal calculation model, capable of predicting the peak overpressure in dynamic blast.

Keywords: shockwave overpressure; dynamic blast; moving velocity; overpressure calculation model

(责任编辑 丁 峰)

Overpressurecalculationmodelofspherechargeblastingwithmovingvelocity

Nie Yuan, Jiang Jianwei, Li Mei

(StateKeyLaboratoryofExplosionScienceandTechnology,BeijingInstituteofTechnology,Beijing100081,China)

O382.1国标学科代码1303520

A

10.11883/1001-1455(2017)05-0951-06

2016-01-29；

2016-06-12

聂 源(1992— )，男，博士研究生;

蒋建伟,bitjjw@bit.edu.cn。