核材料炸药化学爆炸条件下烟云源项的实验与仿真

段中山，过惠平，罗永锋，罗昆升，尚爱国，龚朋彬

(1.火箭军工程大学，陕西 西安 710025；2.陆军勤务学院，重庆 401331；3.火箭军研究院，北京 100094)

在核材料与炸药组成的体系中，一旦炸药发生爆炸将导致核材料气溶胶化并产生大量放射性烟云，对周围人员及环境造成严重辐射危害，常见核材料炸药爆炸事件主要包括脏弹恐怖袭击、核武器化学爆炸以及粗糙核装置爆炸。目前的研究大多运用高斯烟云模式、随机游走算法等风场扩散模型开展外场核材料化学爆炸事故评估，以获取事故后放射性污染强度和范围[1-3]。不管采用哪种方法，都需先获取核材料参数、炸药参数、稳定烟云高度H、稳定烟云半径R、空载气溶胶释放率(ARF)、可吸入比率(RF)和气溶胶可吸入释放率(RRF，RRF=ARF×RF)等源项参数。由于外场烟云稳定状态大多在高空位置，H、R、ARF、RF难以实验测量，而大多只能采用经验估算。当前外场烟云高度与半径计算大多采用Church公式(H=76M0.25、R=13.13M0.25，M为当量)[4-5]，而ARF大多取值地面气化率，即假设爆轰生成的气溶胶全部进入到稳定烟云中，RF取值于地面粒径测量结果[6-8]。事实上Church公式主要由53 kg以上TNT爆炸及核爆炸实验结果拟合得到，未讨论其用于计算小尺度TNT爆炸烟云高度误差，而从外场炸药爆炸到形成稳定烟云这段时间已有一定比例气溶胶颗粒就地发生了驱散或沉降，并非所有地面气溶胶全部进入烟云并辐射下游地区，即ARF、RF值与地面气化率、粒径测量值都存在明显差距，简单且粗糙的源项估算方法导致事故评价结果误差较大。

为解决以上问题，本文拟开展烟云扩散抬升阶段外场实验与仿真，建立一种可用于烟云扩散抬升阶段的数值仿真方法，并用烟云实验测量值验证仿真方法的正确性，最后在仿真结果中提取实验测量难以获取的烟云参数及气溶胶参数变化规律。

1 爆炸烟云气溶胶扩散理论与实验

1.1 核化学爆炸烟云扩散理论

图1 SCM示意图Fig.1 Schematic of SCM

源项模型(SCM)将炸药爆轰过程分为3个阶段[9]，具体如图1所示。第1阶段是爆轰热力学阶段。爆炸瞬间，高温、高压的气体产物迅速膨胀，同时温度、压力不断下降，当火球膨胀至地面烟团直径最大时刻，高温、低密度且包含放射性气溶胶的初始烟团形成[6]。第2阶段是爆炸烟云主动抬升阶段。爆炸膨胀完成后，烟团在大气浮力作用下将脱离地面开始扩散上升，同时空气不断进入、密度变大、温度降低，期间大部分气溶胶颗粒将随烟云扩散上升，大粒径颗粒由于重力作用而沉降，小粒径颗粒由于碰撞而凝并，沉降和凝并则会导致ARF和RF发生变化，该阶段是一个复杂的流体力学过程，可采用计算流体动力学(CFD)方法开展数值仿真。第3阶段是大气被动扩散阶段。烟团到达顶高后基本稳定，烟团和烟柱中的气溶胶随风场逐步扩散、沉降并辐射下游区域。

1.2 爆炸烟云时空分布实验

为获取烟云扩散时空分布过程并检验数值仿真结果的正确性，本文实施了TNT爆炸烟云扩散外场实验。实验设计了1 kg密度为1.63 g/cm3的TNT裸装药爆炸装置，装置上围均匀包裹少量惰性金属粉末用于烟云示踪，大气稳定度、风速、地面条件、炸药壳体等都将影响烟云扩散，实验选择在大气稳定(阴天、傍晚)、无风(v<0.5 m/s)、硬质水泥地面条件下开展。实验简易流程如图2所示，首先根据经验公式对烟云的扩散范围进行大致估算[4]，标定爆炸烟云扩散视野并调整照相机视场，使烟云在可视阶段的飘散范围全部包含在视场范围内；然后选取理想实验条件时刻起爆装置，同时用照相与摄像方式记录爆炸烟团初始形态及动态演变过程；最后根据视场标定结果，使用像素计算软件对实验视频及图片进行分析，获取爆炸烟云扩散时空分布参数。

图2 烟云扩散实验简易流程Fig.2 Simple process of cloud diffusion experiment

2 爆炸烟云气溶胶扩散的理论方程与仿真计算

2.1 烟云扩散流体动力学基本方程

根据烟云浮力扩散理论描述，爆炸烟云源项主要受重力、上升阻力和大气密度差带来的浮力等影响，烟云在浮力作用下的扩散上升过程中遵循质量、动量和能量守恒方程。

质量守恒方程为：

(1)

其中：ux、uy、uz为x、y、z方向的速度分量；t为时间；ρ为密度。

烟云扩散过程中满足x、y、z3个方向的动量方程，以x方向为例：

(2)

能量守恒方程为：

(3)

其中：E=h-p/ρ+u2/2为流体微团总能量，包含内能、动能和势能；h为焓，J/kg；hj为组分j的焓，J/kg；keff=k+kt为有效热传导系数，W/(m·K)，k为流体本征热导率，kt为湍流热传导系数，根据所用的湍流模型来确定；Jj为组分j的扩散通量；Sh为化学反应热及定义的其他体积热源项[10]。

2.2 烟云扩散仿真建模方法与计算

根据SCM对烟云扩散过程与机理的描述，采用CFD方法对烟云扩散开展数值仿真[11]。通过在Fluent求解器中编写UDF源项，结合火球现场测量与公式R=1.93M0.32/(T/3 600)1/3(T=2 861 K)确定烟团直径和高度，1 kg TNT爆炸将形成半椭球形态烟团，水平直径约4.17 m、高约3.13 m。火球膨胀到最大尺寸后压力略低于或等于空气压力，由ρ=ρ0(T0/T)可得烟团整体密度范围为0.08～0.12 kg/m3。烟团组分按TNT化学爆炸式理想气体产物体积分数设定，由于气溶胶颗粒整体体积分数小于1/10，采用离散相模型反复创建出空气动力学直径为1、10、30、50、70、100 μm的球形颗粒来描述气溶胶颗粒。由于气溶胶比例较小，在外场开放空间碰撞并凝并概率低，计算中暂不考虑其对ARF和RF的影响[12]。考虑到温度下降迅速，将烟团混合物平均温度设为约1 000 K，结合实验现象设置烟团的初始翻滚抬升速度约为5 m/s。

计算采用60 m×60 m二维空间。网格从中间到两边、从下到上以1.02倍步进递增方法划分，该方法划分的网格在烟云扩散区域密集度高、计算更准确，下边界为墙面，其他边界为自由界面。选取非定常、2阶方程计算方式及标准k-ε湍流方程。对x方向和y方向速度、能量、各物相分数及k值等进行收敛设定和残差监视，计算中残差线下行收敛并达到计算要求精度[10]。

3 实验与仿真结果讨论

3.1 放射性烟云演变形态分析

图3a～h较完整地记录了不同时刻实验烟云形态，图中烟云抬升速度由快到慢，浓度不断降低直至稳定状态。图3a～d显示，烟云扩散前期由于示踪剂驱散导致形态不够规整。图3e～h显示，后期烟云在浮力主导下其扩散的整体呈现椭球形态，烟团上端分布均匀、下端烟柱明显。

图3 实验烟云扩散过程Fig.3 Experimental cloud diffusion process

图4 气溶胶烟云扩散仿真结果Fig.4 Simulated aerosol cloud diffusion

图4为烟云扩散到约25 m高度的气溶胶分布计算结果(图3f时刻)，图片清晰显示出此时气溶胶随烟云扩散的整体分布与形态，由于烟云密度梯度差产生的浮力将气团由传统的高斯形态浓度分布转变成由两个反旋涡环流组成的涡流分布，横截面上展现出由两个反旋涡环流所组成的“肾形状”模式，即“蘑菇云”现象。涡环围绕其涡核由内向外翻卷，直径逐渐变大且上升过程非常稳定，仿真计算结果与实验观察吻合，表明气溶胶上升驱动力主要来自于低密度高温烟团浮力并非完全依赖初始爆轰驱散。

3.2 放射性烟云时空分布参数分析

烟云扩散高度、直径(宽度)和宽高比随时间变化的实验和仿真结果示于图5。

从图5a可见，实验烟云前期较仿真结果略高，考虑为实验过程存在示踪剂爆轰抛撒及地面反弹，而仿真模型则是理想烟团浮力模型，中后期烟云变化规律基本一致，实验和仿真结果相互验证、数据可信度高。烟云高度随时间变化整体呈现前快后慢的幂函数扩散规律，与文献[13]的结果一致。采用Church模型计算的1 kg TNT爆炸高度高达76 m，而仿真和实验结果均显示烟云约在30 s时刻几乎不再上升，此时烟云稳定高度约38 m，证明Church模型在计算小当量爆高时误差过大，小当量爆炸烟云顶高可采用CFD仿真方法直接计算得到或采用小当量实验确定。

烟云宽度是计算烟云体积的重要依据，Church模型给出了烟云宽度约为高度0.4倍的结论。图5b显示，实验和仿真烟团宽度变化曲线规律较为一致，烟云宽度随时间几乎呈线性增长，但由于实验存在抛撒导致测量宽度略大，最终烟云宽度约为22 m。图5c为烟云宽度与高度数据对比，斜率显示其宽高比先小后大，最终宽高比可达0.58，这一结果与德国爆炸实验结果[14]吻合。

3.3 烟云中气溶胶扩散参数变化规律

1) 不同粒径气溶胶在烟云中的分布规律

由图4可见，大部分气溶胶颗粒随向上运动的双涡环集聚于爆炸烟云的顶部和茎部形成蘑菇云及烟柱，少量颗粒物将随四周横向运动的气流发生驱散或由于重力作用沉降于地面。为更好地观察烟云中颗粒物的分布规律，截取了烟云中1、10、50、100 μm颗粒物在烟云稳定时刻的分布(图6)。从图6可见，此时涡环高度和宽度已达最大，组成烟云顶端的主要成分是1 μm和10 μm的颗粒物，50 μm颗粒物主要分布在烟云茎部或少量到达顶端，由于大粒径颗粒物稳定上升至顶端需要流场向心力最大位置才能完成，故卷入顶端的少量颗粒物主要位于外环，实验后期也显示烟云外环位置能见度低于中心位置，也证明外环处颗粒物浓度高或较多大粒径颗粒物汇聚在此。100 μm颗粒物几乎不能到达涡环顶高，大多集聚在烟柱中或沉积于地面。

2) 不同粒径气溶胶的扩散沉降规律

烟云在扩散过程中颗粒物将发生水平驱散和垂直沉降[15]，根据气溶胶在100 s时间内变化的统计结果(图7)，可将其运动过程分为3个阶段：第1阶段发生在前10 s，烟云抬升时少量颗粒物随周围小涡环发生水平驱散或未卷入涡环，这些颗粒物在重力作用下迅速沉降于地面并造成烟云中颗粒物数目快速下降，颗粒物粒径越大，沉降越快，比例越高，如100 μm颗粒物4 s时减少约40%，而可吸入粒径颗粒则90%以上随烟云上升；第2阶段为颗粒物停留阶段，颗粒物随烟云上升过程中将经历浮力加速、重力与阻力减速、逐步从烟顶分开经烟柱在空气中悬浮、沉降过程，该阶段颗粒物数明显进入一个停留期，粒径越大停留期越短，如100 μm颗粒物在烟云中停留约40 s之后开始沉积于地面，而70 μm颗粒物在80 s时刻开始沉积于地面，其他粒径颗粒物停留期均大于100 s；第3阶段为悬浮、沉积阶段，经过不同停留期后将呈现大颗粒沉积、小颗粒悬浮现象，该规律与文献[9]经验计算结果基本吻合。

图5 烟云高度、直径和宽高比对比Fig.5 Comparison of cloud height, diameter and diameter and height ratio

图6 稳定时刻不同粒径的颗粒物分布Fig.6 Particle distribution of different sizes at steady time

图7 颗粒物在烟云中停留时间Fig.7 Residence time of particulate matter in cloud

3) 烟云中气溶胶源项参数变化规律

图8 烟云抬升过程中颗粒物的留存比Fig.8 Retention ratio of particles during rising of cloud

图8为烟云抬升过程中颗粒物在空气中的留存比(RR)。烟云主动扩散到稳定状态，由于气溶胶运动造成其ARF和RF与地面值相比发生了变化，1 μm和10 μm颗粒物在空气中的留存比约90%，大部分分布在烟云上部涡环中；30 μm和50 μm颗粒物分别有63%和52%在空气中，大多位于顶部烟云的烟柱或外围；70 μm和100 μm颗粒物驱散或沉降较多，分别只有46%和36%留在空气中。以Operation Roller Coaster(ORC)实验质量归一化分布为例，地面实验测量ARF为100%、RF为20%、RRF为20%，粒径区间为<10、10～30、30～50、50～70、70～100、>100 μm时的质量比(MR)分别为20%、25%、15%、10%、10%、20%[5]。烟云中ARF计算公式为ARF=∑RR×MR，RR值取粒径留存比区间期望值，计算得到ARF约58%，由于ARF变化导致烟云中可吸入比例升高，RF由20%增至31%，RRF由20%变为18%。这一结果说明烟云抬升阶段颗粒物驱散与沉降会导致爆炸点附近辐射强度增加，进入大气扩散模式的气溶胶数量整体减少，但可吸入部分气溶胶数量几乎保持不变。

4 结论

本文通过对核材料化学爆炸烟云的实验与仿真分析，得到以下结论。

1) CFD方法仿真烟云扩散的过程、形态和时空分布参数与实验吻合。该方法可用于核材料化学爆炸烟云扩散实验设计、源项计算、近场污染评估等。

2) 烟云流场由两个反旋涡环流组成，横截面展现出一个“肾形状”模型，烟云中可吸入粒径气溶胶大部分分布在烟云涡环中且颗粒粒径越小离涡环中心越近，粒径在50 μm以上颗粒在重力作用下沉降于地面或集聚在烟柱中。1 kg TNT爆炸烟云高度约38 m、宽度约22 m，最终宽高比可达0.58。

3) 由于驱散与沉降效应烟云中气溶胶大多经过初期减少、颗粒物停留、逐步沉降3个阶段，沉降导致稳定烟云ARF与RF均发生一定变化，以ORC实验质量归一化分布为例，其烟云稳定时刻ARF约为地面气化率的58%，RF由20%增至31%，RRF约为18%。