强爆炸与化学爆炸是当前研究较多的2种爆炸现象^[1]，也是国防安全中重要的研究课题。它们在能量转化方式和爆炸毁伤方面存在显著的差异。由于能量密度较低，化学爆炸产生的温升效应较弱，介质的热辐射影响基本可以忽略，毁伤效应以冲击波破坏为主。强爆炸的爆炸过程中，大量能量被瞬间释放^[2]，产生强烈的温升效应，并发射X射线加热周围冷空气形成高温高压火球，其毁伤作用除冲击破坏外还伴随较强的热辐射毁伤。热辐射毁伤是强爆炸火球现象的典型问题，开展热辐射规律研究能提升对强爆炸过程的认识，并支撑毁伤评估。

目前，关于强爆炸火球问题，理论、数值和实验研究均取得了显著进展。强爆炸火球的热辐射特征强烈依赖于火球尺度和辐射自由程特征尺度。辐射自由程是在一定温度和密度下爆炸场介质的固有属性，而火球尺度受辐射自由程和爆炸能量的影响。对于不同爆炸高度的强爆炸问题，随着爆炸高度的增加，大气介质密度逐渐稀薄，辐射自由程占主导作用。对于不同初始火球半径的强爆炸问题，随着爆炸能量的减小，火球尺度的影响逐渐增强，以至于可显著改变强爆炸火球的热辐射特征。强爆炸火球问题中的尺度效应受到了学者们的广泛关注。

基于几何相似原理，乔登江^[3]详细介绍了强爆炸理论中的立方根相似律。Zhang等^[4]理论分析了辐射热传导条件下辐射流体方程的相似变换，推导了辐射能和穿透深度的尺度律。Fournier等^[5]采用NIF（National Ignition Facility）实验设备研究了空中爆炸和地面爆炸时的冲击波演化规律，指出由于实验设备尺度的影响，需加入惰性气体来增加实验气体的吸收截面，改变辐射自由程后，实验结果更符合空中爆炸和地面爆炸的结果。Bouquet等^[6]系统分析了实验室中辐射流体的尺度律，并介绍了它在天体物理和辐射激波中的应用。Koenig等^[7]和Vinci等^[8]开展了高能量密度氙气介质实验，研究了辐射激波的温度和膨胀规律；赵多等^[9]采用数值模拟方法研究了该实验条件下氙气的发光特征。这些研究对于理解不同尺度下的火球现象具有重要意义。

对于热辐射自由程占主导的情形，学者们也开展了较多的研究。辐射自由程是火球辐射流体力学中的一个重要参量，它是绝对线度，使得火球的某些参量不满足立方根相似律；随着高度的增加，当辐射自由程增长到始终大于火球半径时，火球现象将消失^[3]。乔登江^[3]总结了强爆炸火球极值特征问题中时间、半径、有效温度随爆炸当量变化的经验公式，爆炸高度通过大气介质的稀薄程度（密度）来描述。孙景文^[10]指出：爆炸高度小于30 km时，空中爆炸火球将形成2个热脉冲；爆炸高度高于30 km时，由于空气稀薄，热激波不足以形成阻挡层，火球的热辐射呈现单脉冲特征。Brode^[11-12]采用数值模拟方法研究了不同当量和不同高度条件下的立方根相似律，模拟结果不满足立方根相似律。Svettsov^[13]基于球对称方法对低空、中空等高能量密度爆炸问题开展了数值模拟研究，结果显示：随着爆炸高度的增加，极小亮度和第2极大亮度出现的时间提前，且后者提前的时间更长；爆炸高度等于30 km时，随着当量的增加，火球热辐射的极值特征更加明显。田宙等^[14]研究了不同爆炸高度下火球阵面密度、火球中心温度、计算区域辐射能量随时间的变化规律，火球阵面参量和内部参量不满足立方根相似律。随后，田宙等^[15]研究了爆炸高度为40～60 km时火球的演化规律，不同爆炸高度下，火球中心温度的时间分布存在明显差异；中、低空爆炸时，火球热辐射呈现显著的双脉冲极值特征；高空爆炸时，火球热辐射呈现单脉冲特征。由此可见，不同爆炸高度下，火球现象具有显著的尺度效应，通过介质密度描述热辐射极值特征时，其适用范围仍难以界定，不同的文献给出了不同的结论。

综上所述，强爆炸火球问题中的2种尺度效应难以统一，缺乏系统研究。本文中，基于辐射流体热传导近似模型，理论推导包含尺度效应的相似参数来界定热辐射极值特征的适用域，选取火球特征尺度和辐射自由程特征尺度差异较大的2类典型问题来验证尺度效应相似参数的有效性，并采用高精度Euler辐射流体计算程序来模拟火球热辐射对相似参数的依赖性。

1.   计算方法

1.1   控制方程

采用基于Euler框架的一维球几何辐射流体力学方程组求解强爆炸火球问题，方程组的形式为：

式中：t为时间；r为半径；$\rho$为介质密度；$u$为介质速度；$E=\rho(e_{\mathrm{m}}+e\mathrm{_r})+\rho u^2/2$为流场单位体积总能，${e_{\mathrm{m}}}$为介质内能，${e_{\mathrm{r}}}$为辐射内能；$p = {p_{\text{m}}} + {p_{\text{r}}}$为流场压力，下标m表示流场介质，下标r表示辐射；${S_\rho }$、${S_{\rho u}}$和${S_E}$分别为质量、动量和能量守恒方程源项；$\delta$为几何因子，对于柱几何，$\delta = 1$，对于球几何，$\delta = 2$。

在强爆炸火球的数值模拟中，忽略重力、黏性等源项的影响，即${S_\rho } = 0$，${S_{\rho u}}{\text{ = }}0$，能量方程源项${S_E}$主要为辐射源。采用Zinn^[16]发展的辐射输运近似方法求解式(1)～(2)。辐射能群离散支持单群灰体近似和多群近似，在流体力学方程组的数值离散中，对流项的求解采用五阶WENO（weighted essentially non-oscillatory）格式^[17]，对流通量求解采用LF（Lax-Friedriches）方法，时间项求解采用Euler方法。采用自适应网格计算方法^[18]提高计算效率，高精度计算方法可有效提升辐射波阵面的分辨率，从而提高热辐射的模拟精度。

1.2   数值方法验证

参考文献[3]中表征热辐射特征的经验公式，验证数值方法和物性参数的可靠性。爆炸当量为1 kt时，爆炸高度为海平面高度。Brode等^[19]和Symbalisty等^[20]的研究表明，当关注的火球特征不同时，设置的爆炸初始条件也应当不同。对于热辐射问题，计算中初始温度为6.0×10⁶ K，高温区内外介质的密度相同。在数值方法方面，Svettsov^[13]和Zinn等^[21]的研究表明，采用拉式激波捕捉方法，数值振荡较为严重，第1极大亮度难以分辨。基于Euler方法求解辐射流体方程，可高精度地获得计算场的温度和密度分布，从而精细地描述热辐射输出规律。

图1给出了爆炸当量为1 kt时，爆炸火球的特征有效半径（effective radius）和有效温度（effective temperature）随时间的变化规律，其中辐射输运采用42群罗西兰平均自由程。模拟计算的火球第1极大亮度（first maximum brightness，FMB）、极小亮度（minimum brightness，MB）和第2极大亮度（second maximum brightness，SMB）时间分别为0.09、3.84和37.90 ms，与采用经验公式^[3]计算的理论值（0.09、4.56和38.00 ms）符合良好。多群近似能够产生高精度的可信结果，但多群辐射参数的每一群都存在特征尺度，理论分析十分困难，本研究拟采用单群灰体近似模型进行理论和数值模拟研究。模拟的初始温度为2.0×10⁶ K，高温区内外介质的密度相同，单群灰体近似的计算结果如图2所示。可以看出，相较于多群近似，单群灰体近似难以表征第1极大亮度，极小亮度和第2极大亮度的时间略微后移，但大体趋势基本一致。单群灰体模型可用于描述热辐射特征。

图  1  爆炸当量为1 kt、爆炸高度为海平面高度时爆炸火球有效半径和有效温度随时间的变化规律(42群)

Figure  1.  Effective temperature and radius of yield 1 kt explosion at sea level (42 group)

下载: 全尺寸图片 幻灯片

图  2  爆炸当量为1 kt、爆炸高度为海平面高度时爆炸火球有效半径和有效温度随时间的变化规律（单群）

Figure  2.  Effective temperature and radius of yield 1 kt explosion at sea level (single group)

下载: 全尺寸图片 幻灯片

2.   结果与分析

2.1   尺度效应理论分析

热辐射极值（如极小亮度、第2极大亮度等）特征是强爆炸火球的典型现象。已有研究通过引入密度来表征不同的爆炸高度对热辐射极值特征的影响，但未将辐射自由程的特征尺度与物理现象的特征尺度联系起来^{[4-5, 22]}。密度变化难以同时反映火球尺度和辐射自由程的尺度效应，所以在不同高度、不同当量火球问题中，热辐射极值特征经验公式难以统一。为了解决该问题，本文中，从基本控制方程出发，引入火球尺度和热辐射自由程尺度以表征强爆炸火球中的尺度效应。

对于中低空爆炸，火球的光学基本结构为：火球阵面前和核心部分为光学薄区域，阵面附近为光学厚区域。这种基本结构从第1极大亮度前后开始，一直维持到第2极大亮度^[23]。辐射热传导模型可用于表征光学厚区域的相似性和存在性。

一维球对称辐射流体热传导近似方程为：

式中：c为光速，${a_{\text{R}}}$为辐射常数，${\lambda _{\text{R}}}$为辐射平均自由程。

为导出相似参数并建立它与特征尺度之间的关系，选用海平面温度${T_0}$、海平面高度的声速${a_{\text{s}}}$、高度为h_b处的大气密度$\rho_{h_{\mathrm{b}}}$、大气密度$\rho_{h_{\mathrm{b}}}$对应的最小辐射平均自由程${L_{{\text{MFP}}}}$、最小亮度时火球半径${L_{{\text{FB}}}}$作为特征参量，对式(1)～(3)中的部分物理量进行无量纲化处理：$\overline T = T/{T_0}$，$\overline u = u/{a_{\text{s}}}$，$\overline{\rho}=\rho/\rho_{h_{\mathrm{b}}}$，$\overline{t}=t/(L_{\text{FB}}/a_{\mathrm{s}})$，$\overline r = r/{L_{{\text{FB}}}}$，${\overline \lambda _{\text{R}}} = {\lambda _{\text{R}}}/{L_{{\text{MFP}}}}$，$\overline{p}=p/(\rho_{h_{\mathrm{b}}}a_{\text{s}}^2)$，$\overline{E}=E/(\rho_{h_{\mathrm{b}}}a_{\text{s}}^2)$。选用海平面温度作为温度的特征参量，是因为不同海平面高度下密度和压力的变化比温度剧烈，温度的影响基本可以忽略。

无量纲化的辐射流体热传导近似控制方程可以表示为：

式中：${\rho _0}$为海平面高度的气体介质密度。

由式(4)可知，等式右端存在2个相似参量，即$ca_{\text{R}}T_0^4/(3\rho_0a_{\mathrm{s}}^3)$和$\rho_0L_{\text{MFP}}/(\rho_{h_{\mathrm{b}}}L_{\text{FB}})$，定义后者为尺度效应相似参数f_scale，前者约为常数1.27×10⁻⁵。${L_{{\text{MFP}}}}$可从辐射不透明度物性参数中直接提取，而火球特征尺寸${L_{{\text{FB}}}}$是唯一需要确定的相似参量，这是本文理论分析的难点。

在冲击波扩张阶段，冲击波半径^[3]为：

式中：Q为TNT当量，kt；ρ₀为介质初始密度，kg/m³；时间t的单位为s；R_s的单位为m。

对于海平面高度爆炸的火球现象，在极小亮度时刻的火球半径^[24]为：

由于极小亮度时间与密度基本无关^[3]，不同海平面高度下极小亮度时刻的火球半径与密度之比呈0.21的指数关系，即

辐射流体热传导近似控制方程的相似性参数由此确定，相似参数中成功引入了火球特征尺度与辐射自由程特征尺度。对于实验室特征尺度的相似参数，温度变化不大时，将不同高度下的介质密度等效为实验室介质密度。表1给出了爆炸当量为1 kt时不同爆炸高度下的尺度效应参数。当f_scale<1时，火球特征尺度占优，辐射输运满足光学厚假设；而f_scale>1时，辐射自由程特征尺度占优，光学厚假设不再满足，辐射热传导近似也将不再成立，此时的尺度效应参数仅标识辐射热传导模型偏离光学厚假设的程度。

表  1  爆炸当量为1 kt时不同爆炸高度下的尺度效应参数

Table  1.  Scale effect similarity parameters of yield 1 kt at different altitudes

hb/km	$\rho_{h_{\mathrm{b}}}$/(kg·m−3)	LMFP/m	LFB/m	fscale
0	1.22×100	5.98×10−5	32.00	1.87×10−6
10	4.13×10−1	3.87×10−4	40.19	2.85×10−5
20	8.84×10−2	5.70×10−3	55.55	1.42×10−3
30	1.82×10−2	9.32×10−2	77.45	8.10×10−2
35	8.30×10−3	3.70×10−1	91.30	5.97×10−1
40	3.90×10−3	1.44×100	107.01	4.22×100
50	9.92×10−4	1.67×101	146.62	1.44×102
60	2.93×10−4	1.49×102	184.19	3.37×103
70	7.58×10−5	1.71×103	244.72	1.12×105

下载: 导出CSV 

| 显示表格

爆炸高度不同时，尺度效应相似参数之比与密度之比呈指数关系，如图3所示。通过指数拟合方法，可以确定指数关系为：

图  3  爆炸高度不同时尺度效应相似参数与密度的关系

Figure  3.  Scale effect similarity parameter vs. density at different altitudes

下载: 全尺寸图片 幻灯片

式中：f_scale,0为海平面高度的尺度效应相似参数，$f_{\text{scale},h_{\mathrm{b}}}$为高度h_b处的尺度效应相似参数。

将式(8)代入第2极大亮度时刻经验公式中，可得第2极大亮度时刻：

式中：${t_{\max ,2}}$的单位为ms。

孙景文^[10]的研究表明，单脉冲和双脉冲热辐射特征是爆炸高度为30 km以上和以下爆炸的显著区别。然而，采用极小亮度和第2极大亮度时间来表征单脉冲和双脉冲热辐射特征^[13]时，区分高度约为45 km，与孙景文^[10]的研究成果（30 km）有较大差异，需要引入其他参数来表征热辐射特征。已有的研究通过密度之比来表征爆炸高度，与式(8)～(9)的表述相同，因此，可以引入f_scale来表征热辐射特征。

对于实验室尺度的火球现象，尺度效应的分析较为复杂。实验室条件下，火球的尺度较小，且常常通过混合气体来改变介质的辐射自由程尺度以满足相似律。典型案例中，火球尺度会缩小，假定实验气体的温度和密度均不变。表2给出了海平面高度（h_b=0）时不同火球尺度（R₀）所对应的尺度效应参数。可以看出，随着火球尺度的减小，f_scale逐渐变大，但即使对于0.1 mm量级的火球，火球尺度依旧占优，火球的热辐射呈现双脉冲特征， f_scale可为实验设计提供理论支持。表3给出了爆炸高度（h_b）为20 km时不同火球尺度所对应的尺度效应参数，当R₀≥1.0 mm时， f_scale>1，辐射自由程尺度占优，热辐射的双脉冲特征消失。

表  2  海平面高度时不同火球尺度所对应的尺度效应参数（$h_{\mathrm{b}}=0$）

Table  2.  Scale effect similarity parameter in different laboratory scales ($h_{\mathrm{b}}=0$)

R0/m	$\rho_{h_{\mathrm{b}}}$/(kg·m−3)	LMFP/m	LFB/m	fscale
2.81	1.225	5.98×10−5	32.00	1.87×10−6
1.00	1.225	5.98×10−5	9.86	6.07×10−6
1.00×10−1	1.225	5.98×10−5	7.14×10−1	8.37×10−5
1.00×10−2	1.225	5.98×10−5	5.17×10−2	1.16×10−3
1.00×10−3	1.225	5.98×10−5	3.75×10−3	1.60×10−2
5.00×10−4	1.225	5.98×10−5	1.70×10−3	3.52×10−2
2.00×10−4	1.225	5.98×10−5	5.98×10−4	1.00×10−1

下载: 导出CSV 

| 显示表格

表  3  $h_{\mathrm{b}}=20\;{\mathrm{km}}$时不同火球尺度所对应的尺度效应参数

Table  3.  Scale effect similarity parameter in different laboratory scales ($h_{\mathrm{b}}=20\;{\mathrm{km}}$)

R0/m	$\rho_{h_{\mathrm{b}}}$/(kg·m−3)	LMFP/m	LFB/m	fscale
2.81	8.84×10−2	5.70×10−3	55.55	1.42×10−3
1.00	8.84×10−2	5.70×10−3	17.11	4.61×10−3
1.00×10−1	8.84×10−2	5.70×10−3	1.24	6.37×10−2
1.00×10−2	8.84×10−2	5.70×10−3	8.98×10−2	8.79×10−1
1.00×10−3	8.84×10−2	5.70×10−3	6.51×10−3	1.21×101
5.00×10−4	8.84×10−2	5.70×10−3	2.95×10−3	2.67×101
2.00×10−4	8.84×10−2	5.70×10−3	1.04×10−3	7.60×101

下载: 导出CSV 

| 显示表格

2.2   尺度效应对辐射自由程占优问题中热辐射特征的影响

对辐射自由程占优的火球热辐射问题进行数值模拟，分析火球热辐射特征随f_scale的变化规律。随着爆炸高度的增加，辐射自由程显著增长，数值模拟的初始条件也随之变化。热空气条件下，辐射自由程与温度的三次方成正比，假定火球的初始半径随辐射自由程的增长而线性增大，可以获得不同爆炸高度下的初始温度。图4给出了不同爆炸高度（h_b）下火球的热辐射功率（蓝实线）和有效半径（红虚线）。结果表明：h_b≤30 km时，火球的热辐射功率曲线呈现明显的极小亮度和第2极大亮度特征； h_b=35 km时，这种极值特征十分微弱，难以分辨；h_b=40 km时，热辐射功率曲线的极值特征消失。孙景文^[10]的研究表明，对于h_b>30 km的强爆炸，由于空气稀薄，热激波不足以形成阻挡层，热辐射呈现单脉冲特征，与本文的结论一致。

图  4  不同爆炸高度下火球的热辐射功率和有效半径

Figure  4.  Thermal radiant power and effective radius at different altitudes

下载: 全尺寸图片 幻灯片

由表1可知，当爆炸当量较大时，随着火球尺度的增加， f_scale减小。图5给出了h_b=35 km时不同爆炸当量下火球的热辐射功率和有效半径的数值模拟结果。可以看出，随着爆炸当量的增加，小当量爆炸中未出现的极值特征又重新呈现。因此，以30 km作为低空和高空爆炸的区分高度并非完全准确，爆炸当量越大，区分高度将越高。爆炸当量仅反映火球尺度，难以反映辐射自由程特征尺度，仅用它也不能准确表征火球热辐射的极值特征。

图  5  爆炸高度为35 km时不同爆炸当量下火球的热辐射功率和有效半径

Figure  5.  Thermal radiant power and effective radius when the altitude is 35 km for different yields

下载: 全尺寸图片 幻灯片

2.3   尺度效应对火球尺度占优问题热辐射特征的影响

对火球尺度占优的火球热辐射问题进行数值模拟，不同的火球初始半径（R₀）对应着不同的爆炸当量。图6给出了海平面高度（h_b=0）条件下R₀取0.2和0.5 mm时火球的热辐射功率和有效半径。结果表明，火球的热辐射功率呈现出明显的极值特征，相较于R₀取2.81 m的火球极值特征，极小亮度和第2极大亮度时间均大幅提前。由表2可知，h_b=0、R₀取0.2和0.5 mm时，f_scale（0.100和0.035）均小于1，与图6的结论一致。图7给出了h_b=20 km、R₀取10.0和1.0 mm时火球的热辐射功率和有效半径。可以看出，R₀=10.0 mm时热辐射功率曲线的双脉冲特征已较微弱，R₀=1.0 mm时双脉冲特征完全消失。由表3可知，h_b=20 km、R₀取10.0和1.0 mm时，f_scale分别为0.879（小于1）和12.100（大于1），与图7的结论一致。这表明f_scale可以较准确地表征受火球尺度和辐射自由程影响的热辐射特征。

图  6  R₀取0.2和0.5 mm时火球的热辐射功率和有效半径（h_b=0）

Figure  6.  Thermal radiant power and effective radius under different initial radius for laboratory scale fireball (h_b=0)

下载: 全尺寸图片 幻灯片

图  7  R₀取10.0和1.0 mm时火球的热辐射功率和有效半径（h_b=20 km）

Figure  7.  Thermal radiation power and effective radius under different initial radius for laboratory scale fireball (h_b=20 km)

下载: 全尺寸图片 幻灯片

3.   结　论

针对火球热辐射双脉冲极值特征问题，基于辐射热传导近似控制方程，推导了包含尺度效应的相似参数，选取火球尺度占优和辐射自由程占优的2类典型问题来验证尺度效应相似参数的有效性，采用高精度Euler辐射流体计算程序，数值模拟了火球热辐射对于尺度效应相似参数的依赖性，得到的主要结论如下：

(1) 基于尺度效应相似参数，建立了热辐射极值经验公式，该公式可以描述火球热辐射极值特征；

(2) 在辐射自由程占优和火球尺度占优的强爆炸火球问题中，尺度效应相似参数可以较准确地描述火球热辐射演化规律；

(3) 以30 km作为低空和高空爆炸的区分高度并非完全准确，爆炸当量越大，区分高度将越高。

需要指出的是，为了降低理论分析的难度，辐射热传导模型选取了单群灰体模型，并通过数值模拟验证了其有效性，但它在多群模型中的适用性还有待考证，这是下一步工作的重点。