砾石堆结构超高速撞击溅射物特性

张鸿宇，迟润强，孙 淼，王 涵，庞宝君，张 熇

（1.哈尔滨工业大学 空间碎片高速撞击研究中心，哈尔滨 150001；2.北京空间飞行器总体设计部，北京 100094）

引 言

近地小行星（Near Earth Asteroid，NEA）指轨道近日点距离在1.3 AU以内的小行星，按其轨道可划分为4类：Atiras型、Atens型、Apollos型及Amors型[1]。由于与地球轨道存在交叠，Atens型与Apollos型近地小行星具有更高的撞击地球风险。现有研究将距地球最小的轨道距离（Minimum Orbit Intersection Distance，MOID）小于0.05 AU，自身直径大于140 m的小行星定义为对地球构成潜在威胁的近地小行星（Potentially Hazardous Asteroid，PHA）。若其与地球发生撞击，可在局部地区甚至是全球引发重大灾害[2]。截至2023年2月16日，已发现的NEA共31 338颗，其中PHA共2 330个[3]，且随着观测活动的推进，该数目仍在不断增大。

为应对潜在威胁近地小行星的撞击，提出了行星防御概念，并自1994年7月彗木撞击事件后开始得到国际社会的广泛关注，成立了国际小行星预警网（International Asteroid Warning Network，IAWN）、空间任务咨询小组（Space Mission Planning Advisory Group，SMPAG）等国际组织，并提出了核爆摧毁（Nuclear Explosion）[4]、动能撞击（Kinetic Impact）、引力牵引（Gravitational Traction）[5]、太阳光压（Solar Photon Pressure）[6]、“以石击石”[7]等系列主动防御手段。其中针对中长预警时间的危险小行星，采用动能撞击偏转小行星轨道被广泛认可[8-9]，并且在已开展的“深度撞击号”（Deep Impact）与“双小行星撞击转向试验”（Double Asteroid Redirection Test，DART）任务得到了部分验证。

1）Deep Impact以获取彗星9P/Tempel 1的物性特征为主要目标。飞跃目标彗星期间，由探测器发射铜制撞击器与进近传感器（Impact or Targeting Sensor，ITS），质量370 kg，撞击速度约10.3 km/s。对溅射物开展分析，发现其最高溅射速度达5 km/s，溅射物质量约1.4×105kg，彗星9P/Tempel 1位置变了约10 km[10-12]。

2）DART任务撞击目标为双小行星系统Didymos较小的Dimorphos并评估撞击偏转效果[13]。拍摄的图像显示，撞击区域为砾石堆组构，且在撞击后随即产生了大量溅射物并从Dimorphos逃逸。撞击后开展的观测显示，Dimorphos公转周期减缓了33 ± 1 min[14-15]，动量倍增因子β=[16]。该任务为人类首次针对地外天体开展的动能偏转效能评估试验，将很大程度推动该领域未来的发展。

动能撞击偏转任务中，撞击形成的溅射物是需重点关注的问题之一，体现在：①溅射物可产生动量倍增效应[17-18]，增大动量交换效率，目标小行星的动量变化是撞击航天器动量与撞击区域飞出的溅射物动量之和，探测器与远红外数据表明，大多数小行星表面覆盖着砾石堆状风化层[19-20]，且小行星表面重力极低，撞击后大部分溅射物将超过逃逸速度并从小行星飞离，该部分溅射物由于出射速度与撞击速度方向相反，且质量可达撞击航天器的数倍，转移至目标小行星的动量可能显著大于撞击航天器的动量；②溅射物观测数据可辅助撞击效果及小行星表面特性评估，通过拍摄撞击溅射物，可由图像预估溅射物的总质量及溅射速度[11]，辅助计算目标小行星运行轨道可能产生的变化，同时根据溅射物几何形貌及演化特性，利用地面建立的相关模型可推测其表面砾石粒径分布[21]。结合遥感探测获取的物质构成，即可较为全面地了解撞击区域小行星表面风化层物性特征。

本文根据砾石堆结构小行星表层风化层的特征开展建模，采用光滑粒子流体动力学（Smoothed Particle Hydrodynamics，SPH）方法分析其在撞击体超高速撞击下的动力学响应，重点分析砾石堆结构靶体中砾石直径、砾石质量占比对撞击溅射物特性的影响规律与机制，为未来中国将开展的小行星防御任务提供参考。

1 数值模拟

本文应用AUTODYN对撞击体超高速撞击模拟砾石堆结构开展分析，模拟算法采用SPH，该方法可较好地处理超高速撞击过程中所发生的大变形问题，且可较好地反映溅射物形成后的运动。由于研究重点为砾石堆结构对溅射物形成的影响规律与机制，为提高数值计算效率，模型中撞击体直径300 mm，材质为Al 7075-T6，撞击速度6 km/s。靶体尺寸2 000 mm ×2 000 mm× 1 000 mm，由细砾与随机填充的砾石构成，数值模型未施加重力加速度。

1.1 砾石堆模型生成

砾石堆结构靶体模型由随机分布的砾石与细砾两部分构成，模型生成主要包含4个步骤，如图1所示。

1）大尺寸砾石随机分布几何模型建立。本文模拟砾石堆结构考虑2个主要变量：砾石粒径与砾石质量占比。砾石堆粒径设置3个梯度范围，即：小粒径50～100 mm、中粒径100～200 mm、大粒径200～300 mm，且各梯度范围砾石粒径采用靶体30%与50%的质量占比填充。

模型填充时，为模拟砾石分布的随机性，根据相关小行星表面地形的研究，表层砾石的尺寸分布概率密度函数为[22]

其中：d为砾石直径；α为幂指数，根据对小行星Itokawa的研究[23-24]，本文取为3.3；dmin为靶体中的最小砾石粒径，本文靶体细砾模型假设该值为20 mm。在小粒径、中粒径以及大粒径砾石直径范围内进一步细化粒径梯度，以100～200 mm中粒径砾石为例，设置粒径区间100～120 mm、120～140 mm、140～160 mm、160～180 mm以及180～200 mm，并以区间中值110、130、150、170和190 mm作为代表进行建模。根据式（1），对相应区间求取积分即可获得砾石的分布概率，结合靶体质量以及砾石质量占比，可计算获得砾石的具体数量。

应用离散元EDEM软件建立2 000 mm × 2 000 mm ×1 000 mm的砾石填充域，并根据砾石直径与对应数量随机填充。质量占比30%的大粒径、中粒径及小粒径砾石的填充结果如图2所示。填充完成后，输出砾石中心坐标及其半径大小。

图2 砾石随机填充模型Fig.2 Boulders random filling model

2）靶体几何模型建立与网格划分。建立尺寸大小2 000 mm × 2 000 mm × 1 000 mm的靶体几何模型并划分网格，网格大小20 mm，如图3所示。

图3 靶体几何模型建立与网格划分Fig.3 Target geometry and mesh model

3）大尺寸砾石几何模型映射于靶体有限元网格模型。根据输出的砾石中心坐标及半径大小，将所有砾石几何模型映射于靶体有限元网格模型中并完成单元替换。

4）砾石部分与细砾部分有限元网格模型导入AUTODYN。根据最小颗粒粒径尺寸大小，填充SPH粒子直径大小为20 mm。30%质量占比的大粒径、中粒径、小粒径靶体如图4所示。

图4 砾石30%质量占比靶体模型Fig.4 Target model with boulder mass accounting for 30%

1.2 材料模型及参数

1）撞击体材料模型。撞击体材料选择AUTODYN内置材料库中的Al 7076-T6，其强度模型为Steinberg-Guinan，状态方程为shock[25]。

2）细砾材料模型。细砾材料选择AUTODYN内置材料库中的SAND模型。该模型可较好地描述颗粒类靶体的压缩行为、内部应力波传播与衰减，其状态方程为Compaction，强度模型为MO Granular[26]。

3）砾石材料模型。小行星表面砾石的组成和力学特性因天体而异，如：S型小行星表面砾石与普通球粒陨石构成相似[27]，其密度与拉伸强度均高于构成与碳质球粒陨石相似的C型小行星[28]。虽然不同类型小行星表面砾石的力学性质有所不同，但由砾石与细砾构成的混合靶体对溅射物形成与演化的影响规律可能相似。因此，本文通过调研已开展的研究，砾石的状态方程选用Tillotson[29]，强度模型为Von Mises，失效模型为Hydro[30]。

2 数值模拟结果

2.1 溅射物形成过程与特性

溅射物形成与扩展主要经历撞击过程中对应的接触压缩（contact and compression）与开挖（excavation）2个阶段[31]。砾石堆结构靶体由于砾石随机分布，可在上述2阶段均对溅射物产生影响。

2.1.1 接触压缩阶段

该阶段撞击体与靶体表面接触并压缩其向次表层运动。相较于撞击均质颗粒靶体，砾石堆结构靶体表面砾石高于靶体平面，撞击体先与砾石发生碰撞并破碎，少部分撞击体、砾石碎片构成初始高速飞出的溅射物。撞击体主体碎片沿撞击速度方向继续侵入靶体内部，同时压缩细砾进入次表层形成瞬时撞击坑。砾石由于强度较高，中心撞击点外侧砾石受冲击未完全破碎，其对撞击体碎片的继续侵入、细砾的后继运动具有限制作用，导致瞬时撞击坑在砾石分布较少的区域呈近圆形，而在砾石分布较多的区域表现为不规则外轮廓。

砾石的限制作用可能与砾石粒径与质量占比相关，如图5所示，随砾石粒径的减小，砾石的限制作用减弱，形成的瞬时撞击坑趋向于均质颗粒靶体超高速撞击的圆形。而当砾石粒径相同时，增大其质量占比，砾石的限制作用增强，所形成的瞬时撞击坑轮廓不规则度增强，如图6所示。此外，如图6（a）、（b）所示，接触压缩阶段撞击体与靶体表面相交位置处可形成部分紧贴靶体表面的高速喷射物[32]，进而撞击其附近突出靶体表面的砾石，形成次生溅射物，导致初始溅射质量进一步增大〔为方便辨识，引入图6（b）的应力云图进行展示〕。

图6 撞击体与靶体相交处喷射物撞击附近砾石形成的次生碎片Fig.6 Secondary debris formed by boulders near the intersection of the impactor and the target due to jet impact

2.1.2 开挖阶段

该阶段包含溅射幕（ejecta curtain）的形成与演化，也包含撞击坑的继续运动扩展。砾石堆结构靶体撞击后的开挖阶段具有2大特征。

1）撞击坑非对称扩展。已开展的均匀粒径靶体超高速撞击试验，撞击体碎片侵入靶体后形成近半球形不断扩展的空腔，反映于撞击坑边缘则呈对称的圆形扩展。然而，砾石的随机分布形成了大量间断面，增加了靶体的非连续性，其限制作用进一步影响次表层砾石碎片及细砾沿开挖流场的运动，形成非对称扩展的撞击坑。

2）形成射线形溅射物，即砾石碎片与细砾在局部发生汇聚并以高速离开靶体表面。Shualov[33]、Sabuwala[34]等研究表明天体表面存在早期形成的撞击坑、沟壑时，能量将在撞击坑底与沟壑内部发生汇聚并驱使靶体材料呈射线状溅射。而在本文的研究中，溅射物中的射线部分形成于砾石的间隙区域。

砾石粒径及质量填充率对射线形溅射物形成的影响，如图7所示，相同质量占比30%时，中、小粒径砾石工况形成了明显的射线形溅射物，大粒径砾石工况溅射物在砾石位置处仅形成了空缺区；相同质量填充率50%时，所有粒径砾石分布的靶体中均形成了明显的射线形溅射物。

图7 开挖阶段撞击溅射物形貌Fig.7 Morphology of impact-induced ejecta during excavation

砾石粒径与质量占比可对形成的溅射物溅射角与射线数量产生影响。砾石粒径大小相同时，随砾石质量占比的增大，溅射物中射线部分愈加明显，射线数量与长度增大，并且射线位置处的溅射物具有更大的溅射角，如图7和图8所示。而砾石质量占比相同时，随砾石粒径的减小，溅射物射线部分减弱，整体溅射角增大，质量分布更加集中。

图8 开挖阶段溅射物溅射角对比Fig.8 Comparison of ejecta angle at excavation stage

2.2 射线形溅射物形成机制

已开展的均质颗粒材料撞击研究表明[35]，撞击后的靶体形成冲击波以半球形向深层传播，其强度随传播距离的增加不断衰减成为塑性波，最终转变为弹性波。冲击波通过后材料将处于运动状态。当沿半球形传播的冲击波接触靶体表面时，反射卸载波卸载受压缩的靶体材料，形成开挖流场并产生瞬时空腔，近表面的开挖流场驱使靶体脱离瞬时空腔形成溅射物。砾石堆结构靶体中存在大量砾石与细砾相结合的界面，该界面的存在可影响冲击波的传播，并影响开挖流场中靶体的运动，导致射线形溅射物的形成。

50%质量占比砾石靶体撞击过程的应力云图如图9所示，冲击波传播至细砾与砾石的界面处发生反射与透射。细砾作为一种颗粒类材料波阻抗较小，而砾石波阻抗较大，当冲击波传播至该界面时将反射压缩波，同时形成应力幅值大于入射波的透射波。若砾石间存在缝隙，相邻砾石界面附近的高应力区域发生叠加，形成瞬时高压区并与次表层及其他区域形成压力梯度，促使细砾及破碎的砾石碎片以更高速度离开靶体表面，形成射线状溅射物。

图9 砾石50%质量占比靶体冲击波的透射与反射Fig.9 Transmission and reflection of shock wave from boulder target with 50% mass ratio

此外，砾石的介质波速大于细砾，故砾石中的透射波先于细砾中的冲击波到达砾石远端，如图10（a）所示。而该透射波到达砾石远端界面时，从波阻抗较大的介质向波阻抗较小的介质传播，将反射卸载波卸载砾石中的高应力区域，同时形成应力幅值小于入射波的透射波。保持相同的大砾石质量填充率，减小砾石粒径时，如图10（b）所示，砾石中透射波传播距离减小，其与冲击波到达砾石远端的时间差减小，冲击波波阵面在靶体表面逐渐呈现圆形扩展，撞击形成的撞击坑不规则度与溅射物射线部分减弱，如图5和图8所示。

图10 砾石中透射波的传播Fig.10 Propagation of transmitted wave in gravel

值得注意的是，随砾石直径增加以及砾石质量占比增大，冲击压缩区域内砾石间发生接触，促使冲击波沿砾石与砾石构成的通路快速传播，导致撞击体撞击所造成的影响范围增大，如图11所示。Deller等[36]数值模拟同样发现了该特性，而该特性可能促使远离撞击区域的靶体次表层砾石与细砾发生重排，小粒径细砾向距离表面更深的次表层运动[21,37]。

图11 靶体剖面冲击波的传播Fig.11 Propagation of shock wave on target profile

2.3 溅射物动量

统计数值模型计算时间范围内（0.6 ms）溅射物沿撞击速度反方向的动量大小，如图12所示，为大粒径、中粒径及小粒径砾石靶体溅射物动量随时间的变化图。由图12中的（a）和（b）可知，砾石质量占比30%与50%占比时，大粒径砾石（与弹丸直径大小相近似）靶体撞击所产生的溅射物具有最大的动量，中粒径砾石靶体次之，小粒径砾石靶体最小。

图12 砾石不同质量占比靶体的撞击溅射物动量Fig.12 Impact ejecta momentum of boulder with different mass ratios

然而，Ormö等[38]所做的试验及数值模拟结论显示靶体中砾石的存在将导致溅射物动量减小，削弱撞击过程的动量传递。结合本文与Ormö等工况设置，分析可能原因如下：①撞击速度差异，文献试验与仿真速度400 m/s，而本文撞击体撞击速度6 km/s，属于超高速撞击，可导致靶体中更多的砾石发生破碎并参与构成溅射物；②砾石质量占比差异，文献中砾石质量占比较低，而本文工况砾石质量占比为30%与50%，填充得更为紧密，撞击过程中冲击波的影响范围更大[36]，导致形成更多的溅射物；③砾石分布差异，文献中砾石分布于表面细砾下层，而本文中砾石凸出于靶体表面，撞击体先与砾石发生碰撞并使其破碎，形成更大质量的溅射物。

3 结 论

动能撞击偏转作为现阶段可实施性与成熟度最高的小行星防御策略，在面对部分小行星表面复杂的砾石堆结构时，获得最大的动量传递并利用溅射物观测数据评估撞击效果为其待解决的关键问题。本文针对撞击体超高速撞击不同粒径与质量占比的砾石堆靶体开展了数值模拟，获得结论如下：

1）撞击体超高速撞击砾石堆结构靶体，可形成射线形溅射物，位于射线部分的溅射物拥有较其它部分更大的溅射角度，且砾石质量占比越大，射线形溅射物越明显，射线数量与长度增大。

2）溅射物在砾石的限制与导向作用下形成射线形溅射物，可能与相邻砾石间的缝隙区域相关。撞击发生后，靶体中产生近半球形传播的冲击波，当其由波阻抗较小的细砾部分传至波阻抗较大的砾石部分时反射压缩波促使砾石间缝隙范围内压力增大并与其它部分形成压力梯度，导致溅射物以射线状飞出。

3）冲击波在砾石堆靶体中的传播，随砾石直径与质量占比增大，其扰动范围增大，并对远离撞击区域的靶体次表层产生影响。该特性可能引发砾石组构的重排，影响砾石堆结构沿深度方向粒径的分布。

4）本文建立的不同粒径与不同质量占比砾石堆靶体模型中，大粒径砾石（与弹丸直径大小相近似）靶体撞击产生的溅射物具有更大的沿撞击速度反方向动量，可产生更大的动量交换效率。