高空强爆炸X 射线辐照铝靶板动响应的数值模拟

摘要： 高空强爆炸所产生的X 射线辐照至导弹壳体结构时产生的汽化反冲冲量（ blow-off impulse， BOI）及热激波，能够引起目标的动响应破坏。现有的Whitener、BBAY 和MBBAY 理论模型仅能给出一维近似BOI 值，无法处理复杂三维情况并给出对应的热激波峰值压力p，因此，对该问题的研究非常依赖数值计算。利用X 射线热激波数值计算程序TSHOCK3D 对矩形铝靶板在0.1～3.0 keV 范围的普朗克黑体温度和220～400 J/cm2 辐射能通量下的汽化反冲冲量及峰值压力进行计算，并与理论模型作了对比分析。结果表明，TSHOCK3D 程序可以得到可靠的结果，正辐照靶板中心处近似一维工况下的BOI 与Whitener、BBAY 和MBBAY 三个理论模型下的BOI 基本相符。通过单变量分析可得，靶板中BOI 和峰值压力p 均与入射能通量呈近似线性关系；而对于不同的黑体温度，BOI 和峰值压力则在1.5～2.0 keV处存在极大值。

关键词： X 射线；汽化反冲冲量；热激波；TSHOCK3D

中图分类号： O347.3; O434.19 国标学科代码： 13015 文献标志码： A

高空强爆炸实现弹道导弹中段拦截主要依靠其所释放的高能X 射线。不同于大气内空爆及地下实验，高空强爆炸中70% 左右的能量是以X 射线脉冲而非冲击波的形式向外释放。当这一脉冲辐照至导弹壳体时，X 射线光子能量在极短的时间内迅速以指数形式在壳体内发生衰减沉积，转化为壳体材料的内能并引发烧蚀相变，从而实现对导弹壳体表层结构的烧蚀破坏。随着汽化层向近真空的剧烈自由膨胀，所形成的汽化反冲冲量（blow-off impulse， BOI） 会对剩余结构产生压缩反冲的力学加载，形成一个三角形的纯压缩应力波。同时在剩余壳体结构中，由于沉积能量的非均匀分布会产生非均匀的热应力，应力沿梯度下降方向发展传播并最终形成一个近似为正弦曲线的先压缩后拉伸的应力波，其中拉伸波的作用会引起受辐照面的层裂。上述两种机理产生的应力波在传播过程中互相叠加，形成如图1[1] 所示的具有高压缩峰值并伴随有一个较低拉伸波的复杂热激波，即X 射线热激波。

早期的高空强爆炸实验部分验证了X 射线热激波产生与后续动力学破坏效应的物理机理[2]，初步建立起了对X 射线辐照产生破坏效应基本规律的认识。此后则主要利用Z-pinch 类设备模拟产生X 射线源，对缩比尺寸靶板开展辐照实验并精确测量BOI 等相关物理量。这一方面代表性的工作是Sandia 实验室的Remo 等[3-6] 利用Saturn Z-pinch 设备开展的一系列模拟强爆炸的软X 射线辐照实验，靶板涵盖了从多孔岩土到金属材料。近年来，NIF（national ignition facility） 大装置在激光间接驱动惯性约束聚变内爆问题上获得了突破性进展，实现了能量增益大于1 的突破[7]。间接驱动的核心物理问题仍是基于激光辐照加热黑腔后，产生的X 射线烧蚀球形DT 靶丸产生BOI 及球形均匀热激波，从而实现靶丸的等效均匀压缩[8]，相关研究的一个重点在于对热激波的精确控制从而克服其中的流体不稳定性[9]。与实验相对应，在BOI 理论研究方面，最早由Whitner 基于指数递减型能量沉积剖面及材料瞬时汽化假设给出BOI 的大小；Hans Bethe 等提出的BBAY 模型中考虑了汽化区域的复杂相变，并在此基础上通过对液化相变过程及其对冲量的贡献的修正建立起了MBBAY 模型[10]。Lawrence[11] 对上述三种模型在矩形、指数递减形能量沉积剖面下的BOI 进行了理论计算对比，证明了三种模型的结果没有显著差异，液相对BOI 贡献不明显。受限于实验设备，同时期国内大多采用电子束辐照实验进行代替[12-13]，其原理是，利用电子束辐照模拟X 射线辐照产生类似的能量沉积曲线，并假定后续动力学响应规律一致。张朝辉等[14]利用PTS（primary test stand）实验装置达到了文献[3] 的实验加载水平，测量得到了BOI 的值，为后续的数值计算研究提供了非常重要的对比数据样本。

考虑到实验开展成本高昂，目前，开展X 射线辐照材料产生动响应问题的主要研究方法是数值模拟。对于该问题的数值模拟主要是由美国Sandia 实验室开发的CTH、PUFF-TFT[3-4] 等程序，结合上述实验进行的等效一维模拟。Huang 等[15] 采用有限元方法对碳酚醛各向异性材料进行了二维数值模拟；张昆等[16] 将其扩展至三维，编写了TSHOCK3D 程序并完成了X 射线辐照铝板的汽化反冲冲量和碳纤维增强树脂（carbon fiber reinforced polymer， CFRP） 平板的热激波三维数值模拟。Wang 等[17] 和Lin 等[18] 采用一维计算程序RAMA，模拟了碳酚醛材料由X 射线引起的二维热激波，通过实验将数值模拟结果与相同强电流电子束产生的热激波应力进行了比较，验证了数值模拟的正确性，讨论和分析了热冲击波在碳/酚醛材料中的传播规律，并就单能X 射线和软/硬X 射线以不同入射角照射LY-12 铝靶的BOI 开展了等效计算。

综上所述，目前数值计算研究聚焦于结合具体辐照实验，对部分特定工况下的X 射线辐照动力学过程进行计算，研究结论对理解物质在极端条件下的响应具有参考价值，但数值计算工作没有涵盖实际应用中可能遇到的各种辐照工况，计算结果缺乏与BOI 理论模型的系统对比分析。本文中使用TSHOCK3D 程序，对黑体温度及入射能通量这两个核心辐照参数进行多工况遍历计算，将计算得到的BOI 与理论模型得到的进行对比，并对黑体温度-能通量与辐照后BOI-峰值压力的变化规律进行定性分析。

1 高空强爆炸X 射线基本理论

高空强爆炸形成的火球可视为理想黑体辐射源，产生的X 射线存在连续波长的谱分布，各个波长光子所占能量组分f （λ，T） 可由普朗克公式描述：

式中：c1 和c2 分别为第一和第二辐射常数，

c1 = 2πhc2 = 3.743 5×10-12 J ·"cm2/s （2）

c2 = hc／k = 1.439 cm ·"K （3）

式中：h=6.626×10−34 J·s 为普朗克常数；c 为光速；λ 为黑体辐射的波长；kT 为黑体温度，通常用能量来表示，单位为keV，其中k=1.38×10−23 J/K 为玻尔兹曼常数。

如果辐射源表面温度并非均匀分布，其发射谱并非黑体谱，则可以由若干个不同温度的黑体谱加权平均来得到等效的黑体温度k-T：

式中：αi 是温度为kTi 的黑体谱辐射能量与总能量的比值。对于不同种类的核武器，其等效黑体温度不尽相同，但普遍在1～5 keV 之间。表1给出了常见核武器高空强爆炸下的黑体温度和占比的虚拟参考值。

由表1 可知，裂变武器的辐射谱可以由黑体温度为1.0 keV 的辐射谱表示，普通热核武器存在高能的硬X 射线谱，但0.8～1.0 keV 和3.0～4.0 keV 黑体温度所占能量高达97%，增强辐射武器有75% 的能量由0.8～1.0 keV 和3.0～5.0 keV 的黑体提供。因此，本文在能谱工况设计上主要考虑0.1～3.0 keV 段。

2 数值模拟计算流程

采用FORTRAN90 语言编写的基于显式拉格朗日动力学有限元程序的TSHOCK3D[16] 可以对X 射线辐照动力学过程进行三维数值模拟。TSHOCK3D 采用模块化设计，将计算控制文件与靶板模型文件分离。本文中的计算主要考虑辐照源参数的影响，因此，设定靶板参数不变。为了提高多次迭代计算效率，实现样本自动循环生成，在编译好的TSHOCK3D 程序的基础上，采用Python 语言编写了辐照工况遍历计算自动化脚本程序，分两层循环更改控制文件中黑体温度和辐射能通量的数值，并对计算结果自动抓取汽化反冲冲量和峰值压力数值，供后续分析。计算流程也会留存当前计算工况下的等效质量吸收系数μ，用于后续与理论模型的对比中。计算流程如图2 所示。

3 结果对比与分析

3.1 工况设置及理论模型对比

设定X 射线为方波脉冲，脉宽为0.1 μs。黑体温度范围为0.1～3.0 keV，间隔0.1 keV，共30 组；辐射能通量范围为220～400 J/cm2，间隔20 J/cm2，共10 组，双因素遍历组合合计300 组工况，其中单次计算耗时约4 h。

如图3[ 1 6 ] 所示，铝靶板双边长4 mm、厚1 mm，X 射线正辐照加载。考虑靶板及正辐照的对称性，采用1/4 模型进行计算。单元划分采用六面体实体单元，单元在靶板厚度方向的长度为0.01 mm，宽和高方向的单元长度取为 0.16 mm，单元合计62 500 个。铝的相关物态参数见表2[19]，其中E 为杨氏模量，G 为剪切模量，v 为泊松比，σy为屈服应力，ρ0为密度，c0为铝中的声速，s 为Hugoniot 系数，E0 为升华能，Γ0 为Grüneisen 系数。将程序设置为在10、20、50 和100 ns 时记录并输出靶板中的压力三维分布并观察汽化层的运动情况。本文中取黑体温度为1 keV、辐射能通量400 J/cm2的计算结果进行展示，在照射100 ns 后单元内的能量沉积分布云图如图4 所示。

可以看到，该工况下的X 射线对于金属铝靶的穿透能力较差，能量沉积集中在距靶板表层0.1 mm 的厚度范围内，最外层单元获得的比内能远高于升华能，发生相变并向外自由膨胀，产生BOI。表层汽化层的演化及剩余固体靶板中热激波向内传播的过程见图5。

从图5 中可以看到，表层汽化后的物质沿x 轴方向快速膨胀，在靶板边缘自由面处，汽化物质沿y、z 轴也产生了一定的自由飞散运动。靶板中心区域受到了临近区域物质的约束，可视为一维加载，提取该处的BOI 用于后续分析。辐照结束时刻（t=100 ns） 的热激波压力峰值可达26 GPa 左右，随着激波的向内传播，压力峰值逐渐降低，并且当热激波到达靶板后自由面即发生自由面等值反射，所产生的等值反向拉伸波会导致材料产生层裂破坏。事实上，这一现象在靶板自由边界处已发生典型的边侧稀疏破坏，与实验定性相符。上述结果定性证明了本文采用的TSHOCK3D 程序的可靠性，而文献[1，20] 采用的TSHOCK3D 程序，对文献[4] 的1310#实验进行了定量的详细验证，证明了本文使用程序的可靠性，因此，本文研究重点在BOI 及对激波压力同辐照参数的关系。

对BOI 的理论计算，通常需要对输入能通量F0 及汽化反冲冲量I 进行无量纲化处理：

式中：E0 = 10.89 kJ/g 为表2 中铝的升华能，μ 为铝对该黑体温度下X 射线的等效质量吸收系数，该数值与具体光子波长分布及受辐照材料的原子序数相关。以黑体温度为1.0 keV 工况为例，对于本文所研究的纯铝靶板，μ = 814.99 cm2/g。相应的Whitener、BBAY 和MBBAY 三个理论模型的无量纲解分别为：

TSHOCK3D 程序计算结果与上述理论模型的对比如图6 所示。从图6 中可以看到，三种理论模型中的I*与F0*均呈现出比例系数近似的线性增长关系，而Whitener 模型的最高，在低能通量下达到MBBAY 模型的近两倍。TSHOCK3D 程序的计算结果非常接近BBAY 模型，而略高于MBBAY 模型。值得注意的是，上述理论模型均为一维情况理论解，模型均假设X 射线辐照与物质相互作用的能量沉积瞬时完成，即完全不考虑X 射线时间波形和时间脉宽。在BOI 形成机制上忽略材料汽化相变物理过程。因此，理论模型在处理具有很大不确定性的实际工况时非常依赖能量沉积曲线的计算，并且材料的升华能及材料对不同黑体谱的等效质量吸收系数的准确度对理论模型的计算也存在很大影响。三种模型在高能量通量下可通过适当的“等效归一化”[11] 使结果近似相等，但在中低能量通量下仍存在差异。为此，将BBAY 和MBBAY 两模型带有的修正系数α 作为线性比例经验性修正参数（取值为1≤α≤根号下2），上下限取值的不同对结果的影响接近50%，本文中α 设定为1.2。

3.2 TSHOCK3D 结果分析

利用300 组工况的计算结果，将BOI/峰值压力与黑体温度/能通量的关系，在双因素分析下使用Tecplot 软件绘制成三维映射曲面图，如图7 所示。

由图7 可知，BOI 和峰值压力随黑体温度和辐射能通量的变化都有较大的变化，存在明显的全局和局部极值。对于BOI 分布曲面，在硬谱低通量条件（对应黑体温度3.0 keV，辐射能通量 220 J/cm2）下BOI 为零，这是由于在该工况下，X 射线穿透能力较强，沉积能量分布平缓，且入射总能量较低，没有在局部形成足够的能量分布峰值，因此，不足以使得表层物质发生汽化，继而导致没有BOI 产生。但由于仍存在能量沉积形成了热应力，靶板中的压力并不为零。

图8 为单变量分析下BOI 和峰值压力p 随辐照参数的变化。可以看出，BOI 与峰值压力均随着能通量的增加近似呈线性变化，这与2.1 节中的结论一致。对于2.5 keV 等较高黑体温度的计算结果，可以看到，BOI 曲线在部分位置并不平滑，发生了一定的转折和阶跃现象，分析发现，发生这种现象的主要原因，是由于数值计算中，汽化层质量是单位单元质量的整数倍而非连续值，因此，在该工况下汽化层质量发生了阶跃增加，导致BOI 和峰值压力计算结果发生突变，但BOI 与能通量的基本线性关系仍是成立的。对于峰值压力，可以看到，黑体温度为0.5、1.0、1.5 keV 的三种工况，在低能通量下计算得到的峰值压力近乎相同，这一点在三维云图中也可以看出。而更高能通量下，2.0 keV 以上的黑体温度对应的峰值压力曲线有更大的斜率，更低的黑体温度对应的峰值压力曲线斜率则近似相同。BOI 和p 随黑体温度的变化如图9 所示。

从图9 可以看到，相同能通量下，BOI 和p 随黑体温度并没有单调发展，而是在0.8～2.6 keV 区间中存在局部极值点。对于这一现象，从能量分配的角度定性分析认为：X 射线辐照输入的能量，一部分抵消材料相变潜能将固体靶板汽化，提供汽化层的质量m，而另一部分能量转变为汽化层的动能，赋予汽化物质非均匀的初始速度v，两者共同构成了汽化层的反冲冲量即BOI。而相同能通量下，不同的黑体谱及辐照靶板材质，其等效质量吸收系数不同，能量用于抵消潜能及转变动能的比例不同，造成了最终合成mv 即BOI 的大小产生了明显差异。峰值压力也存在类似现象，但各曲线极大值对应的黑体温度相较于BOI 峰值更高，并随着能通量下降，极大值逐渐左移。这是由于在较低黑体温度下，软X 射线能量更多沉积在表面，少量汽化部分获得更高动能有利于提高局部峰值压力。对于这一现象，本文仅进行了定性分析，有待后续深入研究。

对上述确定能通量下存在局部极值现象的推论是：对于强爆炸X 射线的毁伤问题，相同当量的核武器（相同的能通量），不同的爆炸类型（黑体温度）对不同的目标材质的毁伤明显不同，因此，在核武器当量和相对爆距确定的情况下，选取恰当爆炸类型的核武器可产生最优的毁伤效果。

4 结 论

采用TSHOCK3D 程序对不同黑体温度和辐射能通量下辐照铝靶板产生的BOI 和峰值压力进行了数值模拟，将其与Whitener、BBAY 和MBBAY 三个理论模型作对比，并对数值结果做了初步分析，得到以下主要结论：

（1） 对于汽化反冲冲量，TSHOCK3D 程序的计算结果介于三个模型的理论解之间，与BBAY 模型吻合程度最高，TSHOCK3D 程序可以有效可靠地给出模拟结果；

（2）BOI 和峰值压力相对于辐射能通量的变化近似呈线性关系，辐射能通量越高，BOI 和峰值压力越高；

（3）BOI 和峰值压力相对于黑体温度的变化存在局部极大值，能通量越高，极大值越明显。在研究X 射线对导弹壳体材料的毁伤效应时，通过选用合适的核武器战斗部，产生最高BOI 和峰值压力以得到最佳的毁伤效果，是后续值得深入研究的问题。

参考文献：

[1]ZHANG K， TANG W H， FU K K. Modeling of dynamic behavior of carbon fiber-reinforced polymer （CFRP） composite"under X-ray radiation [J]. Materials， 2018， 11（1）： 143. DOI： 10.3390/ma11010143.

[2] 王建国. 高空核爆炸效应参数手册 [M]. 北京： 原子能出版社， 2010： 6–9.

[3]REMO J L， FURNISH M D， LAWRENCE R J. Soft X-ray shock loading and momentum coupling in meteorite and planetary"materials [J]. AIP Conference Proceedings， 2012， 1426（1）： 879–882. DOI： 10.1063/1.3686418.

[4]REMO J L， FURNISH M D， LAWRENCE R J. Plasma-driven Z-pinch X-ray loading and momentum coupling in meteorite"and planetary materials [J]. Journal of Plasma Physics， 2013， 79（2）： 121–141. DOI： 10.1017/s0022377812000712.

[5]REMO J L， LAWRENCE R J， JACOBSEN S B， et al. High energy density soft X-ray momentum coupling to comet analogs"for NEO mitigation [J]. Acta Astronautica， 2016， 129： 384–388. DOI： 10.1016/j.actaastro.2016.09.026.

[6]REMO J L， FURNISH M D. Analysis of Z-pinch shock wave experiments on meteorite and planetary materials [J].International Journal of Impact Engineering， 2008， 35（12）： 1516–1521. DOI： 10.1016/j.ijimpeng.2008.07.075.

[7]LIBERATORE S， GAUTHIER P， WILLIEN J L， et al. First indirect drive inertial confinement fusion campaign at laser"megajoule [J]. Physics of Plasmas， 2023， 30（12）： 122707. DOI： 10.1063/5.0176446.

[8]HURRICANE O A， PATEL P K， BETTI R， et al. Physics principles of inertial confinement fusion and U. S. program"overview [J]. Reviews of Modern Physics， 2023， 95（2）： 025005. DOI： 10.1103/RevModPhys.95.025005.

[9]DO A， CASEY D T， CLARK D S， et al. Measurements of improved stability to achieve higher fuel compression in ICF [J].Physics of Plasmas， 2023， 30（11）： 112703. DOI： 10.1063/5.0167424.

[10]LONGLEY R W. Analytical relationships for estimating the effects of X-rays on materials： AFRPL-TR-74-52 [R]. 1974. DOI：10.21236/ad0786926.

[11]LAWRENCE R J. The equivalence of simple models for radiation-induced impulse [C]// SCHMIDT S C， DICK R D，FORBES J W， et al. Shock Compression of Condensed Matter-1991. Amsterdam： North Holland， 1992： 785–788. DOI： 10.1016/B978-0-444-89732-9.50179-5.

[12]李清源， 王国庆， 吴军， 等. 脉冲电子束产生喷射冲量的实验研究 [J]. 爆炸与冲击， 1991， 11（4）： 339–345. DOI： 10.11883/1001-1455（1991）04-0339-7.

LI Q Y， WANG G Q， WU J， et al. Experimental studies of blow-off impulse generated by a pulse electron beam [J]. Explosion"and Shock Waves， 1991， 11（4）： 339–345. DOI： 10.11883/1001-1455（1991）04-0339-7.

[13]彭常贤， 胥永亮， 徐建波. 电子束辐照平板靶产生喷射冲量的实验研究 [J]. 高压物理学报， 1994， 8（1）： 23–29. DOI：10.11858/gywlxb.1994.01.004.

PENG C X， XU Y L， XU J B. Experimental studies of the blowoff impulses produced in the flat plate targets bombarded by"electron beam [J]. Chinese Journal of High Pressure Physics， 1994， 8（1）： 23–29. DOI： 10.11858/gywlxb.1994.01.004.

[14]张朝辉， 张思群， 任晓东， 等. 基于Z 箍缩X 射线源的热-力学效应实验 [J]. 爆炸与冲击， 2021， 41（9）： 094101. DOI：10.11883/bzycj-2021-0124.

ZHANG Z H， ZHANG S Q， REN X D， et al. Experiments for thermomechanical effects based on Z-pinch X-ray sources [J].Explosion and Shock Waves， 2021， 41（9）： 094101. DOI： 10.11883/bzycj-2021-0124.

[15]HUANG X， TANG W H， JIANG B H. A modified anisotropic PUFF equation of state for composite materials [J]. Journal of"Composite Materials， 2012， 46（5）： 499–506. DOI： 10.1177/0021998311415724.

[16] 张昆， 汤文辉， 冉宪文. X 射线三维热力学效应模拟软件： CN2016SR110024 [P]. 2016.

[17]WANG D W， GAO Y， WANG S， et al. Study on X-ray induced two-dimensional thermal shock waves in carbon/phenolic [J].Materials， 2021， 14（13）： 3553. DOI： 10.3390/ma14133553.

[18]LIN P， CHEN R H， WANG D W. Studies on the thermodynamic properties of C/ph irradiated by intense electron beams [J].Coatings， 2022， 12（8）： 1128. DOI： 10.3390/coatings12081128.

[19]汤文辉， 张若棋. 物态方程理论及计算概论 [M]. 2 版. 北京： 高等教育出版社， 2008.

TANG W H， ZHANG R Q. Introduction to theory and computation of equations of state [M]. 2nd ed. Beijing： Higher"Education Press， 2008.

[20]张昆. 各向异性复合材料的本构关系及其在X 射线辐照下动力学响应的三维有限元模拟 [D]. 长沙： 国防科技大学，2018.

ZHANG K. Constitutive relationship of anisotropic composites and its application in a FEM simulation of the dynamic"response within the X-ray radiation in 3D condition [D]. Changsha： National University of Defense Technology， 2018.

（责任编辑 曾月蓉）