爆炸烟云中粒子运动的数值模拟＊

袁 涛，罗永锋，尚爱国，罗昆升，高洪泉

（1.第二炮兵工程大学，陕西 西安 710025；2.第二炮兵装备研究院第二研究所，北京 100085）

爆炸装置爆炸会产生烟云［1］，烟云的主要成分是被气溶胶化的颗粒物。这些微小的颗粒物可以长时间悬浮在空气中，并随着环境而扩散或迁移。在这个过程中，它们极易被动植物吸附在表面或吸入体内，根据颗粒物的种类不同而造成不同程度的危害。研究这些颗粒物运动和沉积特性是评估爆炸产物对环境影响的理论基础。实地实验从事颗粒沉积的研究非常困难，这主要是因为烟云初始温度高、延迟时间短、变化性强，并且对于其中的物理过程和化学过程知道的非常少［2－3］。爆炸产物中颗粒污染物的排放依据颗粒尺寸可以分为2大类：第1类是粒径小于50μm的颗粒，烟云在大气运动时它们将随气相排放物而分散；第2类是粒径大于50μm的较大颗粒，它们由于重力作用而沉降于地面上。本文中拟通过数值计算，分析爆炸烟云近场扩散的粒子运动规律，选取烟云气溶胶中具有代表性的粒径大小1、10、50、100μm粒子来研究这2类颗粒的运动规律，以期得到的量化结果可用于实验设计，评估爆炸烟云对空气、水体和土壤等的影响。

1 爆炸源项计算

图1 t＝1ms时，源项的压力分布Fig.1 Pressure distribution of source at t＝1ms

本文的研究对象是10kg TNT爆炸装置在开阔无风环境下爆炸的烟云运动问题。根据源项表征模型［4］，将高能炸药的爆炸过程划分为爆炸、浮力气团上升、大气散布3个理化过程。第1阶段爆炸过程是一个极其复杂的理化过程，整个特征时间尺度为几皮秒到几微秒［5］。就本文而言，第1阶段的爆炸过程是第2阶段浮力气团上升过程的源项，所以必须要得到源项的物理特征参数。

为了简化源项计算，本文中做了如下设定：

（1）环境温度为300K，气压为标准大气压，空气为理想气体；

（2）炸药为球状装药，放置于地面中央，装药密度ρe＝1.63g／cm3，半径r＝11.36cm；

（3）第1阶段时间尺度为1ms［6］，即源项的初始时刻为1ms。

Fluent软件不能直接计算爆炸装置爆炸后的物理参数，利用Autodyn软件计算爆炸第1阶段后压力、温度等物理参数，为源项第2阶段在Fluent软件中的计算提供初始参数［7］。接下来在Autodyn软件中建立几何模型，输入初始参数条件，计算得到1ms时源项压力分布如图1所示。

图2 爆炸源项双层模型Fig.2 Adouble－layer model for explosion source

从图1中可以看出，1ms时，爆炸产物的外缘为击波区域，压力p高，形成了明显的分层状态，如果在源项设定中使用简单的球状或椭球状单层结构是不合适的。基于Autodyn的计算结果，本文中建立了如图2所示的爆炸源项双层模型：内层为短轴a1＝1.466m，长轴b1＝1.806m，压力 p1＝18kPa的负压椭球；外层为短轴a2＝1.8 m，长轴b2＝2.5m，压力p2＝530kPa的超压壳状结构。

在源项几何模型中填充粒子要根据爆炸后的物质密度ρ分布。图3表明：密度ρ分布基本服从压力的分布，其规律是内部稀疏、外层紧实，并且外层粒子密度是内层的4～6倍。结合源项双层模型，本文中设置外层粒子密度为内层的5倍。

图4是1ms时粒子速度v与距离爆心距离s关系图。爆炸后产生的能量将给予粒子初始运动速度，速度分布的规律即：击波区域速度最高，击波区域到爆心为线性衰减。

至此，本文源项所需的物理参数已全部得到。

图3 t＝1ms时，源项密度分布Fig.3 Density distribution of source at t＝1ms

图4 t＝1ms时，源项速度曲线Fig.4 Velocity curve of source at t＝1ms

2 几何模型建立及控制方程设置

若要建立源项双层模型，必须在Fluent前处理软件GAMBIT中实现。首先在单位上设置以m为标准单位，按照计算好的源项尺寸绘制。这里有几点需要注意，其一是在公共边的设置上要注意属性上设置为Connect部件，否则将会导致网格错误；其二是外环境的总尺寸问题，原则上来说开阔无风环境应该是无边界最合适。但是在有限元计算时，这个是不允许的，必须要设置一个边界。所以尺寸的选定也很关键：太大会极大地增加网格数，增加计算负荷，降低计算精度；太小又会影响源项压力的释放，导致烟云计算不准确。本文中根据压力衰减的规律和烟云的尺寸预先估算出合适的模型大小，再经过几次粗网格快速计算几个尺寸，比较压力最远传输范围，最后取其中的最优值。

整个区域计算的空间设置在一个60m×60m的区域，源项在底部的正中央。按照1∶1的尺寸建立双层模型，根据源项的物理参数控制源项的网格设置，输出网格，导入到Fluent软件中去。

以2ddp模式启动Fluent软件，在File－Read－Case中选择已经创建好的网格文件（.msh）。使用Grid－Check检查网格是否完整，在Grid－Scale中检查是否以m为标准单位，检查无误后进入基本控制方程的设置。

在Define－Models－Solver中选择Time－Unsteady（非定常）和 Unsteady Formulation－2nd Order Implicit（2阶方程），其他保持默认。勾选 Define－Models－Energy中的 Energy Equation（能量守恒），在Define－Models－Viscous Model中选择k－ε湍流方程。最后勾选Define－Models－Discrete Phase Model中的Interaction with Continuous Phase（粒子与空气相互作用）和Brownian Motion（布朗运动）等选项。接下来进行环境的参数设置。

在Define－Operating Conditions中定义压力为标准大气压，重力加速度为9.8m／s2，空气密度为1.225kg／m3。在Define－Boundary Conditions中分别定义外层边界为标准大气压恒压边界，且为粒子的逃逸界面；底边为地面，且为粒子反射界面。在Define－Materials中分别定义爆炸烟云中物质的密度、传热系数等。

源项粒子的初始设定是通过Define－Injects来实现的，根据本文的要求，分别创建不同粒径的源。在Injection Type中选择Surface，在Release from Surface中选择源项的2个面，在Diameter中输入粒径。如此反复，创建出1、10、50、100μm的粒子源。

最后对求解器进行设置，在Solve－Controls－Solution Controls中选择PISO求解方程，收敛因子保持默认。使用Solve－Initialize中的Solution Initialize和Patch填充双层压力，在Solve－Monitor中设置误差限值，并选定监测误差，使用Solve－Animate和Solve－Particle工具记录粒子运动轨迹和图像。在Solve－Iterate中使用0.01s为步长进行模拟。

3 粒子运动模拟结果

烟云中粒子的扩散与沉积特性如表1所示，不同粒径的粒子在运动时呈现出不同的规律。粒径越小的粒子更容易受到浮力涡环的作用，在双涡环的作用下翻滚着向上运动。构成烟云上端的云团主要是粒径小于50μm的粒子，并且在向上翻滚抛洒的过程中，云团的双涡环由内向外是粒径不断增大的。

从压力演变的图像中可以发现，爆炸后初期涡环上升的速度较快。随着时间的推移，涡环上升的速度变慢，逐渐成为一个动态平衡的状态。涡环中的压力释放速度较快，且涡环会不断膨胀，直至压力与大气平衡，失去束缚云团中粒子的能力。涡环消失时，烟云即会进入SCM模型的第3个阶段，随着大气而进行沉积和迁移［8］。

比较粒径为50和100μm的粒子的运动规律可以发现，粒径大的粒子构成了烟云的茎部，并且粒径越大，粒子的分布就越靠下。粒径为50μm的粒子几乎是均匀地构成了烟云的茎部杆状图像，上部的烟云则是分布在外围。而粒径为100μm的粒子在重力的作用下，沉积得更快，绝大部分粒子已经沉降到了底部的区域。与涡环的浮力作用相比，质量大的粒子更易在重力的作用下沉积。

美国劳伦斯·利弗莫尔实验室利用实验数据和经验公式开发了经典的计算常规爆炸烟云上升的Hotspot程序［9］。这个程序系统地描述了烟云上升后的表征尺度，并得出烟云上升后直径与高度的比值约为0.4。

提取本文计算的结果，烟云在5和10s的表征尺度分别为：Ht＝5＝41.4m，Dt＝5＝13.2m；Ht＝10＝49.8m，Dt＝10＝16.2m，即比值分别为0.32和0.33。这与 Hotspot软件得到的结果相比偏小，主要原因是本文模拟的对象是被气溶胶化的粒径范围，而实际中会产生粒径大于100μm的粒子。粒子质量越大，在上升的过程中就容易被抛洒得更远，使烟云的直径增大。所以，本文的比值与Hotspot软件基本吻合。

表1 不同时间的粒子分布与压力分布Table1 Distribution of particles and pressure at different times

4 结 论

（1）爆炸装置爆炸后会迅速产生双涡环结构，在负压涡环的作用下粒径小的颗粒（50μm以下）极易受涡环的影响而翻滚上升运动。（2）在烟云运动的抛洒过程中，由于惯性的作用，质量大的颗粒分布在烟云的边缘，质量小的颗粒则分布在顶部烟团中。（3）烟云上升的过程中，构成顶层云团的粒子是粒径在50μm以下的粒子，粒径在50μm以上的大颗粒在重力的作用下沉积到地面或形成杆状茎部。

本文研究结果能够揭示爆炸后烟云中粒子运动的规律，可为爆炸后产物对环境的污染和人类的危害的研究提供了理论基础，同时，该结果也能为实验的设计提供参考。

［1］郑波，陈力，丁雁生，等.温压炸药爆炸抛撒的运动规律［J］.爆炸与冲击，2008，28（5）：433－437.Zheng Bo，Chen Li，Ding Yan－sheng，et al.Dispersal process of explosion production of thermobaric explosive［J］.Explosion and Shock Waves，2008，28（5）：433－437.

［2］Morton B R，Taylor G，Turner J S.Turbulent gravitational convection from maintained and instantaneous source［J］.Proceedings of the Royal Society A：Mathematical，Physical and Engingeering Sciences，1956，234（1196）：1－23.

［3］Makhviladze G M，Roberts J P，Yakush S E.Modelling of atmospheric pollution by explosions［J］.Environmental Software，1995，10（2）：117－127.

［4］Brown R C，Kolb C E，Conant J A，et al.A predictive capability for the source terms of residual energetic materials from burning and／or detonation activities［R］.ARI－RR－1384，2004.

［5］Zhang X J，Ghoniem A F.A computational model for the rise and dispersion of wind－blown，buoyancy－driven plumes：Ⅰ：Neutrally stratified atmosphere［J］.Atmospheric Environment：Part A：General Topics，1993，27（15）：2295－2311.

［6］刘庚冉，曹保榆，周凯元，等.FAE野外空中爆炸效应实验研究及其一次起爆机理探讨［J］.防化研究，2007（1）：7－11.Liu Geng－ran，Cao Bao－yun，Zhou Kai－yuan，et al.Experimental investigation of FAE blast field effect in the free air and the discussion of the single ignition mechanism［J］.Anti－chemical Warfare Research，2007（1）：7－11.

［7］周保顺，张立恒，王少龙，等.TNT炸药爆炸冲击波的数值模拟与实验研究［J］.弹箭与制导学报，2010，30（3）：88－90.Zhou Bao－shun，Zhang Li－heng，Wang Shao－long，et al.Numerical simulation and experimental research on TNT explosion shock wave［J］.Journal of Projectiles，Rockets，Missiles and Guidance，2010，30（3）：88－90.

［8］郑毅，应纯同，李磊.稳定大气条件下爆炸烟云的数值模拟研究［J］.力学与实践，2009，31（3）：31－34.Zheng Yi，Ying Chun－tong，Li Lei.Numberical study on explosion clouds in stably stratified atmosphere［J］.Mechanics in Engineering，2009，31（3）：31－34.

［9］HOTSPOT 2.06help manual［M］.Lawence Livermore National Laboratory，2005.