成立, 夏银, 李璐璐, 李百文, 胡晓棉

(北京应用物理与计算数学研究所, 北京100094)

1 引言

小行星撞击地球引发过10 次以上不同程度的生物灭绝事件, 是人类长期共同面临的重大潜在威胁[1]。 截至2022 年10 月17 日, 已发现的近地小行星(轨道近日点小于1.3 个天文单位的小行星) 有30088 颗, 其中直径大于140m 的有10245颗, 直径大于1km 的有857 颗[2]。

为了应对近地小行星撞击地球事件的发生,人们开展了多种防御技术的探索研究[3-5]: 核爆防御技术、 动能撞击技术、 引力牵引技术、 太阳光压技术、 拖船技术、 激光剥蚀技术、 低推力推进技术、 质量驱动技术、 离子束引导技术等。 2022 年9月26 日, 美国国家航空航天局(NASA) 成功实施了双小行星重定向测试任务(DART), 通过动能撞击技术成功使双小行星系统子星的轨道周期缩短了23min[6], 开启了人类主动防御小行星撞击的新篇章。

面对短预警时间(一般小于5 年)、 大尺寸小行星的撞击, 核爆防御是目前唯一可应对的技术手段[3]。 核爆防御技术可分为两种[7]: 一种是利用核爆产生的能量直接摧毁具有潜在威胁的近地小行星, 另一种是利用核爆效应改变小行星的轨道, 避免与地球相撞事件的发生。 由于摧毁小行星后产生的碎片仍然具有撞击地球的风险, 因此大部分核爆防御技术的研究集中在偏转效应的理论分析[8-11]。

核装置在外太空爆炸时, 释放的力学能量将大幅减少, 大部分能量以X 射线的形式释放出来, 还有部分能量由发射出来的中子、 γ 光子所携带。 人们针对核爆X 光、 中子、 γ 光子对小行星的偏转作用已经开展了相关研究。 有研究表明, 核装置产生的γ 射线对小行星的偏转几乎可以忽略, 而中子和X 射线对小行星的偏转起主要作用[12]。 文献[13]指出, 中子对小行星产生的有效动量变化可能更高, 最高可达X 射线的2 倍。 文献[14]针对1MeV 和14.1MeV 两种单能中子对小行星的偏转效应进行了细致的理论分析, 结果表明两种单能中子源的能量沉积分布及能量耦合效率有较大差别, 由此导致偏转速度的差异可达70 %。 文献[12]分析了核爆X 射线和中子对不同物质组成的小行星偏转效应, 结果表明,物质组成对偏转速度的影响可以达到50%。

尽管已经有大量文献开展了核爆防御小行星的理论研究, 但是大部分都是针对单一核爆效应、 假定能量份额及能谱分布条件下开展的物理规律研究。 事实上, 核爆对小行星的偏转效应是X 光、 中子及γ 光子等多物理场对小行星耦合作用的结果。 因此需要针对核爆多物理场耦合作用下, 对小行星的偏转机理开展细致的理论研究。由于核爆是应对大尺寸小行星撞击的唯一技术手段, 本文将针对直径1km 的矿物质小行星, 开展核爆多物理场耦合作用对小行星偏转效应的理论研究。 文章第2 节介绍理论计算模型, 第3 节分析核爆多物理场造成小行星偏转的物理机理, 第4 节为总结。

2 理论模型

本文主要讨论对峙爆炸对大尺寸小行星的偏转效应, 这里取小行星直径1km, 材质为二氧化硅, 密度2.0g/cm3, 计算模型如图1 所示。 核装置在距离小行星表面NG=23m 处爆炸, 所产生的X 光、 γ 光子、 中子在小行星表面进行能量沉积,造成小行星表面升温、 液化甚至汽化, 并向外抛射物质, 通过火箭效应改变小行星动量, 实现对其轨道的偏转。

图1 核爆多物理场对小行星偏转效应计算模型示意图Fig.1 The modeling of simulating the deflection of asteroid due to coupled interactions of nuclear explosion

本文采用非平衡辐射流体力学程序研究X 光及多物理场耦合作用对小行星的偏转效应, 辐射流体力学方程组见文献[15]。 为了研究核爆多物理场对小行星不同区域的影响, 如图1 所示, 以小行星中心点O和爆点N的连线为极轴, 将核爆多物理场对小行星的作用区域按与ON的夹角划分为9 个区域。

利用辐射流体力学程序开展数值模拟, 需要将核爆多物理场对小行星的作用转化为辐射流体力学方程的能量源项。 核爆产生的中子、 γ 光子在小行星表面沉积能量的具体分布可以通过蒙卡粒子输运程序数值模拟给出。 由于中子、 γ 光子的能谱分布、 能量份额对模拟结果具有较大影响, 为了接近真实外太空核爆场景, 计算时中子、 γ 光子的能谱采用如图2 所示的氢弹高空核爆的参考出壳能谱[16], 并且中子、 γ 光子的能量比值取为9。 利用蒙卡模拟给出小行星表面沉积能量的分布情况, 可以开展中子、 γ 光子对小行星偏转效应的理论评估。

针对核爆产生的X 光效应, 可以通过加辐射温度源的方式, 利用多群辐射输运程序考察X 光对小行星表面的烧蚀过程以及造成的动量增量,其中X 射线能谱近似为黑体辐射谱。

针对核爆多物理场的耦合作用, 可以利用多群辐射输运流体力学程序, 同时考察辐射温度源及中子、 γ 光子能量沉积体源的影响, 给出核爆X 光、 中子、 γ 光子多物理场耦合作用对小行星的偏转作用。 本文假定三者作用于小行星表面的时间行为相同, 模拟过程中采用如图3 所示的参考时间谱[16], 其表达式为

式中:β≈0.17ns-1,tp=50ns。

在本文给定核爆距离、 出壳能谱(图2)、 时间谱(图3) 等条件下, 对于100 万吨核装置, X射线出壳能量份额90%, 中子和γ 光子出壳能量份额10%时, 在图1 中G点位置处, 辐射场峰值温度约3 MK, X 射线的表面能注量约57 MJ/cm2,中子的表面能注量约5.7 MJ/cm2, γ 光子的表面能注量约0.63 MJ/cm2。

图2 氢弹高空核爆的中子(左) 和γ 光子(右) 的参考出壳能谱Fig.2 Reference shell energy spectrum of neutrons (left) and γ photons (right)of hydrogen bomb detonated in high altitude

图3 核爆X 光、 中子与γ 光子的参考时间谱Fig.3 Reference time spectra of X-rays, neutrons and γ photons of nuclear explosion

利用多群辐射输运流体力学程序开展多物理场耦合计算时, 需要根据X 光、 中子、 γ 光子在小行星表面的能量沉积特性, 选择合适的径向网格划分。 由于X 光在小行星表面的穿透深度为厘米量级, 而辐射自由程为微米量级, 因此小行星表面数厘米的网格划分需要很细。 本文在计算过程中采用非均匀网格划分方式, 并对网格的收敛性做了检验。 具体来说, 在深度小于10cm 的范围内, 将网格宽度从0.3μm 通过等比方式逐渐放大到1cm, 与X 光的自由程相匹配; 当深度超过10cm 时, 网格宽度从1cm 通过等比方式逐渐放大到20cm, 与中子和γ 光子的平均自由程相匹配。

在计算小行星偏转速度时, 本文沿用文献[14]给出的方法判断某网格内的物质是否被抛离小行星表面: 即当某网格的温度超过熔点且速度大于小行星表面的逃逸速度时, 就认为它被抛离小行星表面。 图1 中区域i对小行星引起的动量增量可以表示为:

式中:j是抛射物质的径向网格下标,mi,j是i区域j网格的质量,vi,j是相应网格的速度。 将各区域的引起的动量增量求和, 可以得到核爆多物理场对小行星总动量的增量。 由此可以求出小行星的速度增量, 其中M为小行星的质量。 由于小行星表面向外抛射的物质质量占比较小, 可以忽略, 因此计算过程中可以不考虑该部分的影响。

3 结果讨论

3.1 核爆中子与γ 光子对小行星的偏转效应

利用蒙卡粒子输运程序模拟了核爆中子与γ光子照射在小行星上产生的沉积能量分布。 图4给出了单位数量源粒子在小行星上单位质量的沉积能(比沉积能) 随角度和距离小行星表面深度的变化规律。

从图4 可以看出, 中子和γ 光子的比沉积能都随深度的变化迅速衰减。 γ 光子的比沉积能在各个角度基本满足随深度变化指数衰减的规律,中子的比沉积能只在小角度(接近垂直入射) 方向满足随深度变化指数衰减的规律。 对于大角度区域, 如[9.513°, 12.03°] 分区, 在深度10cm范围内, 中子比沉积能迅速衰减; 深度大于10cm后, 中子比沉积能衰减速度明显变缓。 这是由于中子的散射截面比光子的散射截面大, 中子在小行星内部输运过程中, 其散射效应比γ 光子强造成的。 对于角度较大的分区, 小行星表面附近的中子比沉积能主要由直穿中子贡献, 其它区域散射过来的中子贡献较小, 因此比沉积能随深度变化呈指数衰减; 到了一定深度之后, 由于直穿中子的贡献衰减迅速, 其它区域散射过来的中子贡献不可忽略, 因此中子比沉积能的衰减速度变缓。

图4 核爆中子和γ 光子在小行星上的比沉积能随角度及深度的分布图中子的比沉积能分布(左图)、γ 光子的比沉积能分布(右图)Fig.4 Distribution of specific energy deposition of neutrons (left) and γ photons (right)on an asteroid with different angles and depths

利用沉积能分布, 可以给出中子和γ 光子在小行星上的能量沉积效率, 即每个源粒子的平均沉积能与初始平均能量的比值。 对于硅岩材质小行星, 文献[8、 17、 18] 中经常采用的中子能量沉积效率为0.70; 文献[14] 的研究表明, 中子能量沉积效率与中子能量相关, 14.1MeV 的中子能量沉积效率约0.69, 1MeV 的中子能量沉积效率约1.09。 本文计算采用的不是单能中子, 而是接近真实太空核爆情形的高空核爆中子出壳能谱[16], 平均中子能量约2.28MeV, 模拟给出的中子能量沉积效率为1.01。 中子的能量沉积效率大于1 是因为部分中子会与原子核发生放能反应产生次级粒子, 次级粒子与原子核碰撞并沉积能量会提高中子能量沉积效率。

采用高空核爆γ 光子出壳能谱[16]计算出的γ 光子能量沉积效率约0.84。 因此, 相同能量条件下, 中子在小行星上的沉积能量比γ 光子多20%, 考虑到近真实情况下, 出壳中子能量约为γ 光子能量的9 倍, 因此中子在小行星上的沉积能量约为γ 光子的11 倍, 中子对小行星的偏转作用显著大于γ 光子。

图5 给出了核爆出壳中子和γ 光子总能量为500kt TNT 炸药当量情况下, 利用辐射输运流体力学程序数值模拟到600μs 时, 小行星各区域速度的径向分布, 其中横坐标为距离小行星表面的深度。 可以看出, 在中子和γ 光子联合作用下,小行星表面向外抛射物质的深度约10 ～100cm。小行星表面向外抛射物质的速度随角度变大而减小, 这与图4 所示的中子和γ 光子沉积能随角度变大而减小的规律是一致的。 对于中子和γ 光子直射的分区, 表面最大物质抛射速度约9cm/μs。

图5 核爆出壳中子和γ 光子总能量500kt TNT 炸药当量, 600μs 时, 小行星各区域的径向速度分布(以小行星表面外法线方向为正)Fig.5 The radial velocity distribution in each region of the asteroid with a total energy of neutrons and photons of 500kt TNT at 600μs (the direction of the outer normal of the asteroid surface is positive)

表1 中给出了不同能量条件下, 核爆中子与γ 光子共同作用下对小行星产生的偏转速度。 其中,Yn+γ为核装置出壳中子和γ 光子携带的总能量,Yin为照射到小行星表面上的能量,Ydep为沉积到小行星上的能量,为能量沉积效率,KErecoil为沉积能量转化的动能,为动能转化效率, Δv为中子和γ 光子共同作用下对小行星产生的速度增量。

表1 核爆中子与γ 光子共同作用下对小行星产生的偏转速度Table 1 The deflection velocity of the asteroid caused by the combination of neutrons and γ photons

受到小行星对核装置立体张角的限制, 核装置释放的能量只有部分可以到达小行星表面。 在本文给定的模型条件下, 约35%的出壳能量可以到达小行星表面。 对于中子和γ 光子而言, 不同能量条件下的能量沉积效率都很高, 接近100%。沉积能量大部分转化为物质内能, 只有约1%的沉积能量可以转化为动能, 由此造成小行星的动量变化。 表1 中最后一列为不同能量中子和γ 光子对小行星造成的速度增量。 可以看到能量越高, 对小行星的偏转越明显, 500kt 出壳中子和γ光子能量, 可以造成直径1km 的硅岩质小行星约0.34cm/s 的速度偏转。

3.2 核爆X 光对小行星的偏转效应

核爆产生的X 光接近黑体辐射源, 根据核爆威力及出壳X 光能量份额可以确定黑体辐射温度。 通过加辐射温度源的方式, 利用多群辐射输运程序可以考察X 光对小行星表面的辐射烧蚀过程以及对小行星造成的速度增量。

图6 给出了核爆出壳X 光能量500kt TNT 炸药当量时, X 光在小行星表面作用200ns 后(辐射源结束时刻) 的径向温度分布。 从图6 中可以看到, 对于图1 中划分的9 片区域, 角度越小的分区, 小行星表面温升越高, 辐射烧蚀深度越深。 这是由于角度越小的分区, 小行星表面距离爆点越近, X 光的源强越强, 因此小行星表面的烧蚀深度越深、 温升越高。 对于X 光直射的分区, 最大烧蚀深度约0.1cm; 而在深度0.1～1.0cm 区域的物质依然有一定温升, 这主要是通过物质热传导造成的, 但是该区域温度相对较低。

图6 核爆出壳X 光能量500kt TNT 炸药当量,200 ns 时, 小行星表面的温度径向分布Fig.6 Radial distribution of the temperature on the surface of the asteroid with the X-ray energy of 500kt TNT at 200ns

图7 给出了核爆出壳X 光能量500kt TNT 炸药当量情况下, X 光在小行星表面作用600μs 后,小行星径向的速度分布。 与温度分布的规律相同, 角度越小的分区, 由于X 光的强度越大, 表面烧蚀出来的物质速度越大。 从图中可以看到,X 光直射分区的表面物质向外最大抛射速度可以达到约80cm/μs, 远高于相同出壳能量条件下,中子和γ 光子造成的小行星表面最大物质抛射速度(9cm/μs)。 但是X 光使小行星向外抛射物质的深度不到1cm, 远小于中子和γ 光子的作用。

图7 核爆出壳X 光能量500kt TNT 炸药当量, 600 μs 时,小行星各分区物质速度的径向分布(以小行星表面外法线方向为正)Fig.7 Radial distribution of material velocities in different regions of the asteroid with the X-ray energy of 500kt TNT at 600μs (the direction of the outer normal of the asteroid surface is positive)

表2 给出了不同能量核爆X 射线对小行星产生的偏转速度,YX-ray为核装置出壳X 光能量,Yin为照射到小行星表面的X 光能量。 与中子和γ光子类似, 受小行星对核装置立体张角的限制,只有部分X 光可以到达小行星表面, 在本文给定的模型条件下, 约35%的X 光能量可以到达小行星表面。 表2 其它各列的物理含义与表1 相同。与中子和γ 光子接近100%的能量沉积效率不同,X 光的能量沉积效率只有14%～20%, 这是由于X 光的辐射自由程非常短, 大部分辐射能量会随着小行星表面物质的吸收再发射过程返照到太空当中。 相比于烧蚀物质产生的动量变化, 由于返照的光压较小, 以及作用在小行星上的时间极短,因此返照光子对小行星动量的影响可以忽略, 在本文中不予考虑。

表2 核爆X 射线对小行星产生的偏转速度Table 2 The velocity change of asteroid caused by X-rays of nuclear explosion

相比中子和γ 光子, X 光沉积能量的动能转化效率相对较高, 可以达到2%。 尽管如此, 由于X 光的能量沉积效率较低, 相同出壳能量条件下, X 光转化的动能依然较小。 例如, 核爆出壳X 光能量500kt TNT 炸药当量时, 沉积能量转化的动能约0.57kt; 相同出壳能量的中子和γ 光子,沉积能量转化的动能约1.99kt。 相比于能量的差异, 不同效应场对小行星偏转速度的差异更大。对于出壳能量500kt TNT 炸药当量的情况, 中子和γ 光子对小行星造成的偏转速度约0.34cm/s,而相同能量的X 光对小行星造成的偏转速度只有0.039 cm/s, 约偏小一个量级。 这是由于: 一方面, X 光转化的动能只有中子和γ 光子的1/4 左右; 另一方面, X 光的自由程较短, 导致小行星向外抛射物质的深度较浅, 抛出的物质质量Δm比中子和γ 光子的作用偏小约2 个量级, 而动量增量因此, 相同出壳能量条件下, X 光对小行星造成的速度偏转效应远小于中子和γ 光子对小行星的作用。

3.3 核爆多物理场对小行星的偏转效应

前面两小节分别讨论了中子-γ 光子联合作用, 以及X 光对小行星的偏转效应。 实际上, 核爆对小行星的偏转效应是多物理场耦合作用的结果。 本小节将考察不同核爆威力条件下, 中子、γ 光子、 X 光同时作用在小行星上, 对其偏转速度的影响。

在本文数值模拟过程中, 中子、 γ 光子、 X光都采用了接近太空核爆环境的高空核爆出壳能谱。 但是三者出壳的能量份额与核装置具体结构相关, 在考察多物理场耦合作用效果时, 表3 中给出了中子和γ 光子出壳能量份额分别占10%、30%、 50%的模拟结果。 其中,Ytot为核装置释放的总能量,Yn+γ为出壳中子和γ 光子的能量,YX-ray为出壳X 光的能量, 表3 其余各列物理含义与表1 相同。

表3 核爆中子、γ 光子和X 射线的共同作用对小行星产生的偏转速度Table 3 The velocity change of asteroid caused by the coupled interactions of neutrons, photons and X-rays

从表3 中可以看出, 由于中子和γ 光子的能量沉积效率较高, 因此随着中子和γ 光子能量份额的提升, 多物理场在小行星表面的能量沉积效率明显增加, 尽管动能转化效率不断下降, 但转化的动能依然随中子和γ 光子能量份额的增加而增加, 由此造成的速度偏转效果不断增大。 结合表1、 表2 和表3 可以看出, 在本文给定模型条件下, 核爆多物理场对小行星的偏转效果主要来自于中子和γ 光子。 对于核装置放能1000kt TNT炸药当量的情况, 中子和γ 光子能量份额占10%(即100kt TNT 当量) 时: 中子和γ 光子单独作用对小行星造成的偏转速度约0.14 cm/s (见表1);考虑多物理场耦合作用后, 对小行星造成的偏转速度约0.15cm/s (见表3), 90%能量份额的X光使偏转速度增加了约7%。 中子和γ 光子能量份额占50%时, X 光对偏转速度增加的影响约1%。

总的来说, 核爆威力越大, 对小行星的偏转速度越大; 中子和γ 光子的能量份额越高, 对小行星的偏转速度越大。

4 总结

核爆防御是目前唯一可应对大尺寸、 短预警时间小行星撞击地球的技术手段, 本文针对1km直径的硅岩质小行星开展了核爆多物理场单独作用及耦合作用对小行星偏转效应的影响研究。 核装置在距离小行星表面23m 处对峙爆炸, 受到小行星对核装置立体张角的限制, 只有约35%的能量可以作用到小行星表面。 理论分析表明, 中子-γ光子联合作用与X 光的作用机制不同: 中子和γ光子的自由程相对较长, 对小行星的加热深度更深, 能量沉积效率接近100%, 大部分沉积能量转化为物质内能, 约1%的能量转化为物质动能,小行星表面向外抛射物质的深度约10 ～100 cm,抛射物质的最大速度约9cm/μs; 而X 光的自由程比较短, 只能加热小行星浅表层的物质, 受到小行星表面物质对辐射吸收再发射物理过程的影响, 大部分辐射能返照回到太空, 能量沉积效率相对较低, 只有约20%, 尽管小行星受X 光作用向外抛射物质的速度非常快, 最大约80cm/μs,但是由于向外抛射物质的深度很浅(约1cm),导致向外抛射物质的质量及动量较小, 因此X 光造成的小行星偏转速度相对较小。

针对核爆中子、 γ 光子、 X 光多物理场耦合相互作用的研究表明, 在本文给定计算模型及近似条件下, 核装置放能1000kt TNT 炸药当量时,90%能量份额的X 光会使小行星的偏转速度增加约7%, 中子和γ 光子对小行星的偏转起主要作用。

本文针对一个特定的核爆防御小行星场景,分析了核爆多物理场对小行星的偏转效应。 事实上, 影响核爆对小行星轨道偏转效果的因素有很多, 如爆炸距离、 核爆威力、 X 光能量份额、 中子能量份额及能谱分布、 γ 光子能量份额及能谱分布等, 而且还与小行星的材质、 密度、 表面形状等相关。 另外, 本文模拟过程中尚未考虑核装置碎片携带动量撞击小行星的影响。 总之, 未来需要针对核爆防御技术, 开展更加精细的物理建模, 细致研究核爆多物理场耦合作用造成小行星偏转的物理机制。