王晔， 白春华， 李建平， 陈风云

(1.北京理工大学 爆炸科学与技术国家重点实验室， 北京 100081；2.交通运输部水运科学研究所， 北京 100088)



非对称云雾爆炸超压场数值模拟

王晔1， 白春华1， 李建平1， 陈风云2

(1.北京理工大学 爆炸科学与技术国家重点实验室， 北京 100081；2.交通运输部水运科学研究所， 北京 100088)

针对非对称云雾爆炸超压场的分布特性问题，利用LS-DYNA程序对非对称云雾爆炸超压过程进行数值模拟。将数值模拟超压值与空投试验超压结果进行对比，得到了0°、90°和180°共3个方向的峰值超压随距离的变化规律以及不同云雾倾角对超压场分布的影响。研究结果表明：非对称云雾与地面存在一定的倾斜角度，在起爆后地面冲击波轨迹呈现倾斜的椭球形，且各个方向上的超压衰减特性不同，在0°方向冲击波峰值超压较大，且衰减速率最大，在180°方向冲击波峰值超压较低，且衰减速率较小；当云雾倾角由0°增加到8.27°时，传播距离在13～30 m范围内的0°方向上，峰值超压平均提高了7%，在90°方向上，超压值基本无变化，而在180°方向上，峰值超压平均降低了8%，因此增加云雾倾角使得峰值超压强度降低。

兵器科学与技术； 非对称云雾爆轰； 峰值超压； 云雾倾角； 数值模拟

0 引言

无约束云雾的爆炸超压场分布规律，是云雾爆轰[1-3]和石化工业事故预防[4-6]等领域研究的热点问题。峰值超压对人员的伤害程度、周边设备的破坏情况及其随距离变化的衰减特性，是计算安全距离与评估事故灾害的主要依据。

最早关于云雾爆轰的研究是Gey等[7]于1963年提出的，自此学者们对云雾爆炸场相继展开研究。Raju[8]、Liu等[9]和刘庆明等[10]对云雾爆轰超压场进行了研究，分别进行理论与数值计算，拟合不同燃料组分云雾爆炸超压与比距离的关系，得到不同距离处超压的TNT当量值，发现了多相爆轰波的多峰结构压力波，并认为云雾区内峰值超压是恒定的结论。

过去研究的云雾状态多为圆柱对称云雾[11-14]，而针对非对称云雾爆轰超压场的文献相对较少。为模拟更加符合实际情况的云雾爆炸超压场，非对称云雾超压场的研究开始起步[15]。由于云雾在外形上的不规则性，导致了爆炸超压场不同方向上的超压衰减速率和分布差异。本文根据由空投试验确定的云雾形态，对非对称云雾的超压场分布进行了数值模拟，并在试验压力测试结果对比的基础上分析了云雾倾角对爆炸超压场的影响规律。

1 试验与模拟参数

1.1 试验

为获得非对称条件下的云雾模型，利用空投方法实现云雾的形成和爆轰过程，试验现场布置如图1所示。利用热气球将战斗部吊至320 m空中，利用地面上的3个锚点调整战斗部的投放位置。以投放位置的地面投影为中心，在0°、90°和180°共3个方向布置高速摄像系统，记录整个试验过程，每台高速摄像拍摄参数设置为500帧/s，记录时长不少于2 s，型号为V12. 试验中采用了24个压力传感器和24个等效压力罐进行云雾爆炸超压场的测试，传感器和压力罐均为标定过的压力测试装置[16]，每个测点间距为12 m，以投放点地面投影为中心，成方形布置，如图2所示。云雾起爆前的云雾形态如图3所示，基于试验结果确定，非对称云雾形状近似为具有一定角度的“倾斜圆柱体”。这是由于一次起爆时，战斗部与地面成一定角度进行燃料抛撒。图3中红色实线为云雾形状辅助线，并以起爆云雾的时刻为0 ms.

图1 空投试验现场布置Fig.1 Experimental layout

图2 地面超压测试场布置Fig.2 Pressure sensors arrangement for experiment

图3 试验中高速摄像拍摄的云雾形成(0 ms)Fig.3 Cloud formation recorded by high-speed cameras in experiment (0 ms)

1.2 数值模型

基于空投试验结果，确定模拟计算中的云雾形状为一定倾角(4.27°)的圆柱形云雾，云雾半径为11.5 m，中心高度为3.3 m. 为方便分析，定义由于倾角导致最高一侧的云雾高度为0°方向，中心高度为90°方向，较低一侧高度为180°方向，如图4所示。经试验结果记录并结合前期试算，同时考虑云雾的长径比以及试验验证中的需要，确定空气域为直径104 m，高度为20 m的圆柱形空气域。

图4 计算模型Fig.4 Computation model

1.3 计算方法

为简化模型便于计算，假设云雾区域内的化学反应过程极短，反应时间可忽略不计，最终表现为爆轰压力。因此，将云雾爆轰产物简化为复合γ定律的理想气体，并利用压力- 能量多项式状态方程描述爆炸后冲击波压力与气体内能之间的关系。采用任意拉格朗日- 欧拉(ALE)多物质组分模型进行计算[17-18]，为了节省运算时间，采用1/2模型进行计算分析。

1.4 材料参数

计算区域包括空气和云雾两个部分。空气采用MAT_NULL材料模型和线性多项式LINEAR_POLYNOMIAL状态方程，主要计算参数如表1所示。状态方程的表达形式为

p=c0+c1μ+c2μ2+c3μ3+(c4+c5μ+c6μ2)E，

(1)

式中：p是压力；c0、c1、c2、c3、c4、c5、c6是常数；E是比内能；μ为动力黏性系数，取值1.79×10-5Pa·s.

由于云雾区内的燃料颗粒远大于空气分子尺寸，可将云雾等同于抛撒燃料与周边空气的均匀混合物。利用空气密度、燃料质量(130 kg)以及云雾体积，可计算出云雾质量密度。利用压力- 能量多项式态方程描述爆炸后冲击波压力与气体内能之间的关系，表达形式如(1)式所示。

假设云雾为理想气体，常数c0=c1=c2=c3=c6=0，c4=c5=γ-1，对于理想气体，多方气体指数γ取1.4[19]，可知c4=c5=0.4. 根据云雾尺寸可确定云雾内爆炸压力平均值为3.0 MPa. 云雾爆轰超压的主要计算参数如表1所示。

表1 空气和云雾的基本计算参数Tab.1 Calculation parameters of air and cloud

1.5 初始条件和边界条件

根据简化模型，以爆轰时刻为0 ms的云雾区域内初始超压3.0 MPa为初始条件。由于爆炸超压场传播区域为无限边界的空气，故边界条件采用无反射边界条件，地面采用反射边界条件。

1.6 网格精度

为了尽可能减少计算精度带来的误差，对云雾区域的网格大小对计算结果的影响进行讨论。针对区域单元边长取0.005 m、0.01 m、0.02 m和0.03 m共4种条件进行对比计算，取距离为50 m处的压力值p50进行对比分析。计算结果如表2所示。

表2 不同网格尺寸的超压计算结果Tab.2 Overpressure computation of different mesh size

由表2可知，不同网格尺寸同一监测点的计算结果不同。随着网格尺寸的细化，计算结果越精确，但是计算时间却成倍增加，因此综合计算时间和计算结果，选取网格尺寸为0.01 m进行模拟。

1.7 网格的选取和计算时长

终止计算条件由计算时间决定，计算时间为100 ms. 每隔0.01 m对1/2模型进行映射网格划分，体单元形状为六面体，网格总数约为300万。模型总计算时长约为10 h.

2 模拟结果与试验验证

2.1 模拟超压分布

根据超压场等压变化情况，选取2 ms、10 ms、20 ms和40 ms时刻对应的压力变化情况。为了直观观测出超压场变化的情况，将左视图进行俯角20°的旋转，模拟压力云图结果如图5所示。

图5 云雾爆轰的超压分布Fig.5 Overpressure distribution for cloud detonation

由图5可知，当云雾完全爆轰后，2 ms时刻，冲击波波阵面前沿外形与倾斜云雾初始形态近似，呈现一端高一端低的轮廓。到10 ms时，冲击波波阵面的形态逐渐演化为半球形，同时不同方向上的压力分布不均匀，出现局部较强超压，即图5(b)中红色区域。到20 ms时，冲击波轨迹保持半球形态继续向远处传播，压力衰减，存在方向差异性。随着时间的推移，到40 ms以后，冲击波压力不断衰减，扩大作用区域，最终消除压力分布的方向差异性，形成均匀分布的压力场。

2.2 试验验证

基于模拟结果，进行试验验证。由于战斗部在下落过程中出现锥形摆动现象，当一次起爆时，可能出现战斗部与地面不垂直，最终导致云雾形态出现一定角度的倾斜。由于倾斜方向不可控，为便于分析压力不同方向变化情况，定义云雾高度较高的一侧为0°方向，中心为90°方向，较低方向为180°方向。为保持视角与模拟结果的一致性，采用90°方向的高速摄像的图片进行分析，如图6所示。0 ms时刻云雾起爆(见图2)，经过4 ms，云雾完全爆轰，火球形状为半球形状，之后，火球体积不断增大，亮度保持不变，到64 ms时，爆轰产物增多，火球的亮度开始衰减，其体积达到最大值，云雾区内的反应物已完全消耗。虽然试验中的冲击波波阵面轨迹为无色的，但是可利用冲击波轨迹与其周边空气折射率的差异来判断波阵面。图6(a)较明显地显现出波阵面为半球型，随着时间的推移，以半球形状向四周扩散。图6中黑色虚线为波阵面辅助线，在4 ms时，冲击波波阵面为明显的半球形，轮廓明显，在24 ms时，冲击波已传播到画面边缘，且轮廓不明显，到44 ms时冲击波已完全传播出画面，无法直接观察。

图6 高速摄像系统记录的试验结果(90°方向)Fig.6 Experimental results recorded by high-speed camera(90°)

为了对比试验与数值模拟结果的差异，将0°、90°和180°方向上的试验超压值和模拟压力值进行对比，如图7所示。图7中模拟压力值，是利用指定监测单元来获取的。

图7 试验与数值模拟超压结果对比Fig.7 Comparison of simulated and experimental overpressures

由图7可知，随着传播距离的增加，在不同的方向上，模拟和试验的冲击波超压值均不断衰减，虽然存在一定的误差，但在整体趋势上，冲击波传播区域峰值超压的模拟结果与试验值具有较好的一致性。导致误差的原因，是由于模拟中的假设条件造成。为计算简便，在模拟中假设倾斜圆柱云雾侧面状态完全相同。但在试验中，由高速摄像拍摄的结果(见图2)可以明显发现，云雾一端侧面存在上浮、鼓包状的不规则湍流云团，另一端云团模糊，较为稀疏，从而导致初始爆轰超压值的不均匀分布。由于试验中初始爆轰超压分布不均匀，而模拟中假设成均匀爆轰，因此导致在传播距离为15～30 m之间(近场区域)的误差较为明显。当超压进行衰减后，渐渐趋于平衡，因此在30 m以外的远场传播区域符合良好。

3 不同倾角的云雾爆炸超压场模拟

由于战斗部与地面存在的角度不确定，为分析云雾倾角对爆炸超压场等距离处峰值超压的影响，分别采用云雾倾角θ为0°、4.27°和8.27°模拟爆炸超压场的峰值超压变化过程。获得0°、90°和180°方向上不同距离处对应的峰值超压，结果如图8所示。

图8 不同云雾倾角压力分布Fig.8 Overpressure distribution at different angles of inclined cloud

由图8可知，随着云雾倾角的增加，3个方向上的峰值超压均表现为随着冲击波传播距离的增加而逐渐衰减的趋势。当云雾倾角由0°增加到8.27°时，传播距离在13～30 m范围内的0°方向上，峰值超压平均提高了7%，在90°方向上，超压值基本无变化，而在180°方向上，峰值超压平均降低了8%.在30～50 m的区域内，3个方向上的超压变化受云雾倾角的影响不显著，基本无变化。

根据数值模型可知，3个方向对应的峰值超压变化差异与云雾截面存在直接关系[20]。0°方向对应的是倾斜云雾的较高一端；90°方向对应的是倾斜云雾基本无变化的一侧；180°方向对应的是倾斜云雾较低的一端。所以，云雾区域的不同形态在很大程度上影响了冲击波传播区域的峰值超压。

通过试验结果可知，虽然增加云雾倾角在0°方向上能够提高峰值超压强度，但是在180°方向上却表现为明显降低，因此从总体上看，增加云雾倾角使得峰值超压强度降低。为了确保获得相对较强的云雾爆炸地面超压场，应保持战斗部与地面垂直，使云雾不发生倾斜。

4 结论

基于有限元分析法，对非对称条件下的爆炸超压场在3个方向的分布规律进行了数值模拟，并与试验超压测试结果进行了对比验证。在此基础上，分析了云雾倾角对爆炸超压场的影响。主要结论有：

1)非对称云雾与地面存在一定的倾斜角度，在起爆后地面冲击波轨迹呈现倾斜的椭球形，且各个方向上的超压衰减特性不同。在0°方向冲击波峰值超压较大，且衰减速率最大，在180°方向冲击波峰值超压较低，且衰减速率较小。可通过改变云雾倾向实现特定方向的超压强弱分布，实现不同的能量输出需要。

2)在不同云雾倾角的模拟条件下，峰值超压强度随着云雾倾角的增加而相对减弱。当云雾倾角由0°增加到8.27°时，传播距离在13～30 m范围内的0°方向上，峰值超压平均提高了7%，在90°方向上，超压值基本无变化，而在180°方向上，峰值超压平均降低了8%. 在30～50 m的区域内，3个方向上的超压变化受云雾倾角的影响不显著，基本无变化。为了确保获得相对较强的云雾爆炸地面超压场，应保持战斗部与地面垂直，使云雾不发生倾斜。

References)

[1] 许会林,汪家骅.燃料空气炸药[M].北京:国防工业出版社,1980:10-18. XU Hui-lin, WANG Jia-hua.Fuel air explosion [M].Beijing: National Defense Industry Press, 1980:10-18.(in Chinese)

[2] 白春华,梁慧敏,李建平,等.云雾爆轰[M].北京:科学出版社,2012:4-12. BAI Chun-hua, LIANG Hui-min, LI Jian-ping, et al. Cloud detonation[M].Beijing: Science Press, 2012:4-12.(in Chinese)

[3] 李秀丽,惠君明,王伯良. 云爆剂爆炸/冲击波参数研究[J].含能材料,2008,16(4):410-414. LI Xiu-li,HUI Jun-ming,WANG Bo-liang. Blast shock wave parameters of single-event FAE[J].Chinese Journal of Energetic Materials,2008,16(4):410-414.(in Chinese)

[4] 宫广东,刘庆明,白春华.管道中瓦斯爆炸特性的数值模拟[J].兵工学报,2010,31(增刊1):17-21. GONG Guang-dong, LIU Qing-ming, BAI Chun-hua. Numerical simulation for gas explosion in tubes[J]. Acta Armamentarii, 2010,31(S1):17-21.(in Chinese)

[5] 陈明生,李建平,白春华.非圆截面云雾爆炸超压场数值模拟[J]. 含能材料,2015,23(5):484-489. CHEN Ming-sheng, LI Jian-ping, BAI Chun-hua. Simulation of explosion overpressure distribution for non-circular cross-section cloud[J]. Chinese Journal of Energetic Materials, 2015,23(5): 484-489.(in Chinese)

[6] Bai C, Fan B, Li B. Improvement of a gas intake device and a gas mixing study in a large-scale vessel[J]. Chemical, 2014, 37(10):1751-1756.

[7] Gey W, Mygaard M. Feasibility study of FAX explosive, NAVWEPS-8065[R]. China Lake, CA, US: US Naval Ordnance Test Station, 1963.

[8] Raju M S N. The blast waves from unconfined axisymmetric vapour cloud explosion[D]. Urbana, IL, US: University of Illinois, 1975.

[9] Liu G, Hou F, Cao B, et al. Experimental study of fuel-air explosive[J]. Combustion, Explosion, and Shock Waves, 2008, 44(2): 213-217.

[10] 刘庆明, 白春华, 李建平. 多相燃料空气炸药爆炸压力场研究[J]. 实验力学, 2008, 23(4): 361-370. LIU Qing-ming,BAI Chun-hua,LI Jian-ping. Study on blast field characteristics of multiphase fuel air explosive[J]. Journal of Experimental Mechanics, 2008, 23(4): 361-370.(in Chinese)

[11] Leyer J C. An experimental study of pressure fields by exploding cylindrical clouds[J]. Combust Flame, 1982, 48(3):251-263.

[12] Borisov A A, Gelfand B E, Gubin S A, et al. Detonation of fuel-air mixtures above the surface of the earth[J]. Combustion, Explosion, and Shock Waves,1988,24(2):238-240.

[13] Fishburn B, Slagg N, Lu P. Blast effect from a pancake shaped fuel drop-air cloud detonation(theory and experiment)[J]. Journal of Hazardous Materials, 1981, 5(1):65-75.

[14] 赵传荣,孔德仁,王胜强,等.飞片厚度对冲击波压力峰值衰减特性的影响分析[J]. 振动与冲击,2016, 35(3):135-138. ZHAO Chuan-rong, KONG de-ren, WANG Sheng-qiang,et al. Influence of flyer’s thickness on the attenuation characteristics of shock wave’s pressure peak[J]. Journal of Vibration and Shock, 2016, 35(3):135-138.(in Chinese)

[15] 陈明生,白春华,李建平.多点云雾爆炸波相互作用的数值模拟[J].爆炸与冲击,2016, 36(1):81-86. CHEN Ming-sheng, BAI Chun-hua, LI Jian-ping. Simulation of blast waves interaction for multiple cloud explosion[J].Explosion and Shock Waves, 2016, 36(1):81-86. (in Chinese)

[16] 陈昌明,李建平,白春华.一种新的冲击波超压测试方法:压力响应膜片[J].中国水运,2013,13(1):104-106. CHEN Chang-ming, LI Jian-ping, BAI Chun-hua. A new shock wave test method: pressure response diaphragm[J].China Water Transport, 2013,13(1):104-106.(in Chinese)

[17] 时党勇,李裕春,张胜民. 基于ANSYS/LS-DYNA 8.1进行显式动力分析[M]. 北京:清华大学出版社,2005:184-249． SHI Dang-yong, LI Yu-chun, ZHANG Sheng-min. Explicit dynamic analysis based on ANSYS/LS-DYNA 8.1[M].Beijing:Tsinghua University Press,2005:184-249.(in Chinese)

[18] 白金泽. LS-DYNA 3D基础理论与实例分析[M]. 北京:科学出版社,2005:74-103． BAI Jin-ze. Basic theory and case analysis for LS-DYNA 3D[M]. Beijing: Science Press, 2005: 74-103.(in Chinese)

[19] 张宝钅平,张庆明,黄风雷. 爆轰物理学[M].北京:兵器工业出版社,2001. ZHANG Bao-ping, ZHANG Qing-ming, HUANG Feng-lei. Detonation physics[M]. Beijing:Publishing House of Ordnance Industry, 2001.(in Chinese)

[20] Bai C H, Wang Y, Li J P, et al. Influences of the cloud shape of fuel-air mixtures on the overpressure field[J].Shock and Vibration, 2016(1):1-7.

Simulation of Explosion Overpressure Distribution forNon-symmetry Cloud Detonation

WANG Ye1, BAI Chun-hua1, LI Jian-ping1, CHEN Feng-yun2

(1.State Key Laboratory of Explosion Science and Technology, Beijing Institute of Technology, Beijing 100081, China；2.China Waterborne Transport Research Institute, Beijing 100088, China)

In order to study the distribution characteristics of overpressure field for non-symmetry detonation, the process of non-symmetry cloud detonation is simulated by LS-DYNA code. The simulated and experimental results are compared. The changing rules of overpressure with distance in three directions, including 0°, 90°and 180°, are obtained. The influences of different angles of inclined cloud on overpressure are studied. Results show that there is an angle between cloud and ground. The track of shock wave presents inclined spheroidicity. The attenuation properties of overpressure in different directions are different. The peak overpressure in the direction of 0° is larger than those in other directions, and the peak overpressure in the direction of 180° is the smallest. The attenuation velocity in the direction of 0° is higher than those in other directions, and the attenuation velocity in the direction of 180° is lower than those in other directions. When the angle of inclined cloud increases from 0° to 8.27°, the peak overpressure is decreased by 7% averagely in the range from 13 m to 30 m in the direction of 0°. The peak overpressure remains nearly unchanged in the direction of 90°, but is increased by 8% averagely in the range from 13 m to 15 m in the direction of 180°. The peak overpressure is decreased by increasing the angle of inclined cloud.

ordnance science and technology; non-symmetry detonation; overpressure; angle of inclined cloud; numerical simulation

2016-06-15

国家部委预先研究重点项目(9140A05080507)

王晔(1986—)， 女， 博士研究生。 E-mail: wangye_0422@126.com

白春华(1959—)， 男， 教授， 博士生导师。 E-mail: chbai@bit.edu.cn

O383+.1

A

1000-1093(2017)05-0910-07

10.3969/j.issn.1000-1093.2017.05.010