真实空气中TNT装药爆炸近区冲击波传播特性研究

周丰峻，郑磊，孙云厚，张鏖，夏明

（1. 中国人民解放军军事科学院国防工程研究院，北京 100850；2. 中国人民解放军92611部队，广州 510715）

一、TNT装药强爆炸自模拟解

点源爆炸理论是研究球形装药爆炸冲击波最简单的模型之一，其假定爆轰产物质量无限小，而装药爆炸释放的能量有限，忽略运动中静止空气的影响，因此，不存在任何一个可以表征爆炸的初始线尺度参数，任何一个物理量在空间上的分布随时间的变化是自相似的，即离爆心不远处冲击波后的气体运动满足自相似条件 [1]。Taylor [2]、Von Neumann [3]、Bach和Lee [4]以及Sedov [5]各自独立得到了点源强爆炸自模拟的理论解，与实验比对可以很好地描述爆炸初始阶段的流场特征 [6]。本文基于点源爆炸理论，对平面、柱面以及球面的一维变几何对称面（面源、线源和点源）装药强爆炸解的自模拟解流场进行了分析。

（一）控制方程

爆炸冲击波的传播是一个非常复杂的不定常过程，但满足质量守恒、动量守恒和能量守恒关系。由于爆炸力学所考虑的运动变化快，以至流线之间不会发生明显的动量和能量迁移，因此可以忽略粘性和热传导 [7]。

1. 欧拉坐标下的控制方程

欧拉坐标下流场的状态量主要包括压力P、密度ρ、质点流动速度u、内能E以及温度T，这些量值均为自变量x, y , z , t的函数。

质量守恒：

动量守恒：

式（2）中，F为质量力。

能量守恒：

物态方程：

对于一维对称TNT装药爆炸可将介质的运动看作是一维绝热运动，即介质中各量值只依赖时间和某单一几何坐标而变化，去掉质量力，即F = 0，可将其在平面直角坐标系、柱坐标系和球形坐标系下建立的欧拉坐标下的控制方程组统一写为：，u , v , w分别为

式（5）中，r为欧拉位置坐标；t为时间；N为常数，等于0，1和2时分别对应一维平面对称、一维柱面对称和一维球面对称情形。

2. 拉格朗日坐标下的控制方程

为了提高冲击波阵面局部计算精度，可将方程组写为拉格朗日形式，控制方程为：

再补充一个运动关系式：

式（6）中，R为拉格朗日位置坐标；u为质点速度；P为压力；E为内能；ρ为密度；ρ0为爆炸波前初始空气密度；t为时间；N为相应平面、柱面和球面坐标指数。

为了计算方便，引入质量坐标方程式：

则可将拉格朗日坐标方程组转化为无量纲形式：

式（9）中，N为常数，等于0，1和2时分别为包含一维平面、柱面和球面的一体无量纲方程组。

（二）真实空气条件下TNT装药爆轰产物状态方程

1. 真实空气无量纲形式状态方程

试验表明，当温度大于等于2 000 K或冲击波超压大于等于4 MPa时，实际空气不再满足理想气体状态方程 [7]。随着温度与压力的不断升高，空气中各种分子之间会发生化学反应，如氮氧分子开始离解，经反应生成NO和NO2。空气组分在不同温度下是不同的，当温度高于8 000 K后，空气中各种原子会发生电离反应，空气组分变得更加复杂；当温度达到300 000 K后，任何密度下气体的电离过程都可以完成。在温度升高过程中发生的空气中各原子电离、分子离解及激发的各种分子振动均需要消耗能量，因此，在高温条件下，空气的状态方程与理想气体状态方程有显著不同 [8,9]。俄罗斯学者Kuznetzov在试验基础上，采用统计物理学和冲击波气体动力学的理论方法得到了空气热力学函数表和冲击绝热线的计算表格，Brode [10]对表格数据进行拟合得到了真实空气状态方程的拟合公式 ：

式（10）～式（16）中，e，θ，π，η分别为无量纲内能、无量纲温度、无量纲压力和无量纲密度。

2. TNT爆轰产物状态方程

Brode [10]从能量角度出发，对TNT装药爆轰产物的状态方程进行了详细研究，并给出了与真实空气状态方程形式相同的描述方程：

式（19）中，α1= 1.76，b1= 52.4892，c1= 7.3。

（三）面源、线源和点源TNT装药爆轰阶段流场参数计算

利用量纲分析方法中的Π定理，可得到由主定参量初始密度ρ0、大气声速C0、爆炸所释放的能量E、球形坐标r和时间t构成的面源、线源和点源TNT装药强爆炸的独立自模拟解变量分别为ξ =若冲击波阵面的位置为rs，则ξ取值为ξ0，则可得到面源、线源和点源TNT对称装药强爆炸自模拟解的控制方程组分别如下：

面源装药：

线源装药：

点源装药：

边界条件：在波阵面上，f(1) = g(1) = h(1) = 1；在中心处，f(0) = 0。其中，f，g和h分别为无量纲速度、无量纲密度和无量纲压力。

式（21）～式（23）为封闭方程组，可利用其进行求解。

通过计算可得到面源、线源和点源装药强爆炸自模拟数值解，结果如图1～图3所示。图中横、纵坐标的数值区间（0～1）均为无量纲形式，其中横坐标的相对位置，即r/rs，为无量纲形式的距爆炸中心的距离。

从图1～图3可以看出，面源、线源和点源装药强爆炸解参数随波阵面相对位置分布的特征基本相同，均呈现类似的变化趋势。质点速度随着r/rs下降基本上呈线性下降趋势；压力在r/rs从1降低到0的过程中先下降比较快，当r/rs降到0.5之后直到爆炸中心，压力基本稳定在0.36～0.37；密度随r/rs降低很快接近至0。

试验表明，自模拟解 [7]可以很好地描述爆炸的初始阶段。因此，利用真实空气条件下面源、线源和点源装药强爆炸的自模拟解，作为TNT装药爆炸近区流场特性数值计算的初始条件是精确合理的。

二、球形TNT装药爆炸近区流场特性

（一）爆炸流场物理模型

选取爆轰波传至装药表面的时刻为初始计算时刻，爆轰波的波阵面将整个流场分为两个部分（见图4）。

根据守恒定律，在激波阵面上必须满足如下关系：

质量守恒：

动量守恒：

图1 面源装药强爆炸解参数分布图

图2 线源装药强爆炸解参数分布图

图3 点源装药强爆炸解参数分布图

能量守恒：

式（24）～式（26）中，P0、u0、ρ0、E0分别表示激波前未扰动空气的压力、质点速度、密度和内能；用P、u、ρ、E表示激波阵面后扰动空气的压力、质点速度、密度和内能；D为冲击波速度。

图4 初始计算时刻流场分布示意图

对式（24）～式（26）无量纲化后，可以解出单值对应的激波关系。对偏微分方程组（9）进行差分近似，为了保持接触间断面处和激波跳跃间断面处计算的贯通，采用人工粘性法，在出现压力项处附加人工粘性项ζ，并将其无量纲化：ρ = ζ/η，于是将方程组（9）（N = 2时为球形装药爆炸）变为无量纲差分格式：

式（27）中，l为差分网格整格点，n为时间整步长，e为空间半格l–1/2和时间整格点n + 1处。

其中粘性项为：

关于能量方程的差分，由于参数比较复杂，因此，采用新的方程式，从无量纲状态方程出发，将e进行代换可得到新的形式：

式（28）中，ξ = lnη。

从而得到能量的差分格式：

因上式是一个隐式的差分格式，因此采用二次迭代进行求解。建立激波后装配方法：首先进行激波探查，用三次样条插值求ξ的最大值点，将此位置定为激波阵面位置λs，求得激波阵面上的压力，然后调用波阵面参数的求解子程序求得波阵面上的其他参数。数值计算的时间步长由Courant条件和激波区域的差分方程稳定性条件共同确定。

（二）算例计算与结果

1. 空气和TNT装药参数

静止空气的参数见表1。

表1 静止空气参数

TNT装药参数见表2。

表2 TNT装药参数

2. 结果

采用本文建立的计算方法可得到下述结果：

（1）超压

图5为无量纲峰值超压ΔPs、峰值动压qs随无量纲激波阵面距离位置的变化曲线，从图中可看出，本文的计算结果与Baker [11]的计算结果吻合度较好。

（2）接触间断面和主激波阵面的时–空曲线

图6为主激波（MS）、二次激波（2S）以及接触间断面（CS）无量纲位置λ随无量纲时间τ的变化规律。从图中可看到Friedman [12]发现的二次激波现象，而出现的多个后续激波是由于波的反射产生的。从图6可得到冲击波的影响范围远远大于爆轰产物接触间断面的扩散范围。

（3）峰值超压、峰值动压随激波阵面距离变化曲线

图7和图8分别为不同装药量球形TNT装药在真实空气中爆炸产生冲击波的主激波峰值超压及峰值动压随激波阵面距离的变化曲线图（采用的是指数坐标）。由图7、图8可知，对于同一装药量，球形TNT装药随着距离的不断增大，超压峰值和动压峰值不断下降。不同装药量的炸药爆炸超压、动压随距离变化符合相似律关系。随着装药量的增加，相同位置处的超压峰值和动压峰值相应增加。

图5 峰值超压、峰值动压随激波半径的变化曲线（实心点为本文结果，空心点为Baker的计算结果）

图6 冲击波与接触间断面的时–空曲线

（4）冲击波到达时间随激波阵面距离变化的曲线

图9为不同装药量球形TNT装药在真实空气中爆炸产生的冲击波到达时间随激波阵面距离的变化规律。可见随着装药量的增加，球形TNT装药爆炸冲击波到达相同位置处的时间不断降低。不同装药量的冲击波到达时间变化规律大致相同，可见同样满足相似律关系。

（5）冲击波超压冲量随激波阵面距离变化的曲线

图8 不同装药量球形TNT装药爆炸峰值动压随激波阵面距离变化的曲线

图9 不同装药量球形TNT装药爆炸冲击波到达时间随激波阵面距离变化的曲线

图10为不同装药量球形TNT装药爆炸冲击波超压冲量随激波阵面距离的变化曲线。从图中可知，对于同一装药量，球形TNT装药随着激波阵面距离的不断增大，冲量先逐渐降低，当到达一定距离后，在此区域会出现一个平缓甚至稍微增加的阶段，之后再开始衰减，这是由于二次激波追赶上了主激波使得正压作用冲量增加造成的；另一方面，随着装药量的增加，相同位置处的冲量也逐渐增加。

（6）冲击波的正压作用时间随激波阵面距离变化的曲线

图11是不同装药量的球形TNT装药空中爆炸冲击波的正压作用时间随激波阵面距离变化的曲线。从图中可知，对于同一装药量，爆炸冲击波随激波阵面距离的变化规律与图10类似，正压作用时间也存在一个转折点，它先随距离的增加而降低，但在某个位置处又开始随距离的增大而增大。这主要是因为存在二次激波的原因，随着冲击波不断向外传播，主激波的峰值超压不断降低，且主激波越晚到达的点的正压作用时间就会越少，但是由于二次激波的不断追赶，当到达某个位置后，主激波的正压还未衰减到负压，而二次激波就追赶上来，从而使得正压作用持续时间增加。对于1～50 kg的TNT装药，转折点的位置彼此相差不大，基本在4 m左右；而对于100～1 000 kg的大装药量药球，转折点离爆心的距离随装药量的增大而增大。另外，对于装药量为1～50 kg的药球，在固定的位置点，转折点前后均有正压作用时间随装药量的增大而增大的趋势；但对于装药量为100～1 000 kg的药球，在转折点前，正压作用时间同样随装药量的增大而增大，但是在转折点后则相反，正压作用时间随装药量的增加反而降低（见图11中的箭头显示）。

图10 不同装药量球形TNT装药爆炸冲击波超压冲量随激波阵面距离变化的曲线

三、结论

由于空气中爆炸的近区空气与理想气体之间存在较大差异，因此，本文采用真实空气状态方程来描述爆炸近区流场的特性。基于本文所建立的理论，结合人工粘性技术和激波装配技术，对TNT装药真实空气中爆炸流场参数进行了计算，结果与已有文献研究结果相吻合，并得到以下结论：

图11 不同装药量球形TNT装药爆炸冲击波正压作用时间随激波阵面距离变化的曲线

（1）TNT装药爆炸之后会产生一个初始超压为60 MPa的强冲击波（图5中无量纲波阵面位置λs= 1对应的峰值超压），强度随冲击波的不断传播而逐渐衰减，并且在主激波之后存在一个不断向外传播的接触间断面。

（2）爆炸近区最典型的物理特征是存在一个爆轰产物区，在半径小于3倍装药半径的区域内，压力强度远远大于主激波的压力值。

（3）当爆轰波到达装药与空气的接触面时，在形成向外传播的空气冲击波的同时，还会形成一个向内传播的激波，向内传播的激波在爆心处发生碰撞后会得到加强，然后继续向外传播，在到达接触面后又形成一个新的向内传播和向外传播的激波，从而形成定点超压波形负压区的二次、三次激波。

由此可见，本文提出的计算分析方法可以很好地模拟和描述TNT装药爆炸近区在真实空气中的传播规律。本文的研究成果可以对爆炸力学和防护工程研究起到重要的参考作用。

[1] 卢芳云. 一维不定常流体动力学教程 [M]. 北京∶ 科学出版社,2006.

Lu F Y. Tutorial of one-dimensional unsteady fl ow dynamics [M].Beijing∶ China Science Publishing &Media Ltd., 2006.

[2] Taylor G. The formation of a blast wave by a very intense explosion. I. Theoretical discussion [J]. Proceeding of the Royal Society of London. 1950 (201)∶ 159–174.

[3] Von Neumann J. The point source solution [R]. Washington, D.C.∶National Defense Research Committee, 1947.

[4] Bach G G, Lee J H S. An analytical solution for blast waves [J].AIAA Journal, 1970 (8)∶ 271–275.

[5] Sedov L I. Similarity and dimensional methods in mechanics [M].New York∶ CRC Press, 1965.

[6] 汤文辉. 冲击波物理 [M]. 北京∶ 科学出版社, 2011.

Tang W H. Shock wave physics [M]. Beijing∶ China Science Publishing &Media Ltd., 2011.

[7] Л.П.奥尔连科. 爆炸物理学(原书第三版) [M]. 孙承纬译. 北京∶ 科学出版社, 2011.

Л.П.ОРЛЕНКО. Explosion physics [M]. Translated by Sun C W. Beijing∶ China Science Publishing &Media Ltd., 2011.

[8] 乔登江. 核爆炸物理概论 [M]. 北京∶ 原子能出版社, 1988.Qiao D J. Introduction to nuclear explosion physics [M]. Beijing∶Atomic energy Press, 1988.

[9] 周丰峻, 陈叶青, 任辉启. 爆炸近区空气冲击波规则反射和非规则反射 [C].应用力学进展论文集. 北京∶ 中国力学学会, 1999∶1–12.

Zhou F J, Chen Y Q, Ren H Q. The regular reflection and nonregular re fl ection of shockwave in the close-in fi eld of TNT charge explosion [C]. Proceedings of applied mechanics progress. Beijing∶The Chinese Society of Theoretical and Applied Mechanics, 1999∶1–12.

[10] Brode H L. Blast wave from a spherical charge [J]. Physics of Fluids, 1959, 2 (2)∶ 217–229.

[11] Baker W E. Explosives in air [M]. Austin∶ University of Texas Press, 1973.

[12] Friedman M P. A simplified analysis of spherical or cylindrical blast waves [J]. Journal of Fluid Mechanics, 1961, 11 (1)∶ 1–15.