超高速动能撞击密实小行星的动量传递规律

刘文近，张庆明，龙仁荣，龚自正，任健康，任思远，武 强，宋光明，陈 川，张品亮

（1.北京理工大学 爆炸科学与技术国家重点试验室，北京 100081；2.北京卫星环境工程研究所，北京 100094）

引 言

虽然小行星和彗星撞击地球的概率低，但后果和危害极其严重[1]。动能撞击目前被认为是使对地球有潜在威胁小行星偏离与地球碰撞轨道的最有效方法[2]。2021年11月24日美国成功发射了“双小行星重定向测试”（Double Asteroid Redirection Test，DART）探测器，于2022年9月26日以超过6.2 km/s的速度撞击双小行星Didymos的次星Dimorphos，成功进行了人类历史上第1次全尺寸的小行星动能撞击偏转测试试验[3]。小行星是一个多样化的种群，在形状、密度、强度和结构方面有很大差异[4]。当通过动能撞击偏转技术改变小行星轨道时，小行星物质属性通常是未知的[5]。因此，需要通过理论分析、试验和数值模拟，分析撞击速度和小行星的物理特性对动量传递系数的影响规律，降低动能撞击偏转小行星的风险和不确定性。

Gault等[6]通过超高速撞击试验获得了6.3 km/s撞击玄武岩的动量传递系数。Hoerth等[7]撞击不同孔隙的岩石结果表明，动量传递系数随着岩石孔隙率的增加而降低。Chourey等[8]对多孔干砂、玻璃微球、月球风化层模拟物的撞击试验结果表明动量传递系数均随着撞击速度的增加而增加。George等[9]在美国国家航空航天局（National Aeronautics Space and Administration，NASA）艾姆斯垂直靶场对各种陨石和陨石模拟物进行了一系列超高速碰撞试验，结果表明含水陨石模拟物动量传递系数均显著大于同等孔隙度无水陨石的动量传递系数，动量传递系数随着陨石孔隙率的增加而降低。受试验条件的限制，试验速度大都低于6.5 km/s，而且大部分超高速撞击试验缺乏对小行星模拟材料性能的完整描述。试验室模型试验通常不能直接重现动能撞击小行星表面所涉及的大撞击器尺寸和极高撞击速度[10]。数值模拟可有效弥补模型试验在尺寸和速度的不足，得出大尺度和高速条件下的碰撞结果[11-12]。Raducan等[11]使用iSALE以双小行星重定向测试任务为例进行了2D数值模拟，发现动量传递效率系数主要受材料初始孔隙度和内摩擦系数的影响。Rainey等[12]使用CTH模拟了DART撞击均质小行星模型，分析了模型参数对动量传递效率系数、陨石坑直径和深度的影响，结果表明动量传递效率系数、陨石坑直径和深度主要由材料的强度和孔隙率决定，但相关模拟未与试验结果比较，验证模型参数的准确性。

因此，本文首先开展球形铝弹丸撞击玄武岩靶试验，获得不同撞击速度下成坑的结果和动量传递系数。利用AUTODYN软件光滑粒子流体动力学（Smoothed Particle Hydrodynamics，SPH）法进行数值模拟，将模拟结果与实验室规模试验数据进行了比较，验证了模型参数的准确性。在此基础上，模拟了最高10 km/s撞击速度对抛射物质量、速度分布和动量传递系数的影响。结合试验和数值模拟结果，获得动能撞击密实小行星的动量传递系数相似律，揭示撞击速度对动量传递系数的影响规律。

1 超高速撞击分析

1.1 玄武岩靶

玄武岩经水钻钻取成棒材，然后用磨床进行打磨得到如图1所示的圆柱靶，尺寸Φ123.3 mm×123.3 mm，同时制备了测量岩石力学特性的试验试件，单轴压缩试验试件尺寸Φ50 mm×100 mm，巴西圆盘劈裂试验试件尺寸Φ50 mm×50 mm。在应变速率为1 × 10–4/s的条件下，采用单轴压缩和巴西圆盘试验测定了试件的抗压和抗拉强度。测得玄武岩单轴抗压和抗拉强度的平均值分别为146.0 MPa和9.8 MPa。通过测量抗拉试件质量和体积，得到玄武岩的密度2.876 g/cm3。

图1 玄武岩靶和试件Fig.1 Basalt targets and specimens

为直接测量撞击玄武岩的动量传递系数，玄武岩靶被挂在2根42 cm的尼龙绳上，悬挂在靶室中央，形成弹道摆。为防止岩石材料在受高速撞击时崩落破坏，在玄武岩的侧面缠绕几圈玻璃纤维增强塑料胶带，如图2所示。

图2 双吊绳弹道摆Fig.2 Double rope pendulum

1.2 弹丸和弹托

试验采用如图3（a）所示的直径为6 mm球形2024铝合金弹丸，用聚碳酸酯弹托发射，弹托强制分离系统将弹丸与弹托分离，然后用磁测速测量弹丸速度[13]。为测量无弹丸情况下弹托的速度，将弹托切分为上下两瓣，弹托中心开槽，槽内放置一块铝片后，将分瓣弹托组合在一起，形成内部包含铝片的弹托。弹托内部铝片用于触发磁测速，测量无弹丸情况下弹托的速度，如图3（b）所示。

图3 弹丸和弹托Fig.3 Projectile and sabot

1.3 试验系统和装置

试验加载系统如图4所示，二级轻气炮加速弹丸和弹托，弹丸和弹托出炮口后，用弹托分离器将弹丸与弹托分离[14]，通过磁测速获得弹丸速度信号。为测量无弹丸情况下的弹托速度，将磁测速装置安装在弹托分离器前面，弹托首先通过磁测速装置，弹托内铝片触发磁测速，获得弹托速度信号。速度信号触发高速摄像，高速摄像记录弹道摆的整个摆动过程。高速摄像的拍摄帧率为1万帧/s，靶室内高速摄像机的对面设置了1个LED灯补光。靶室内环境的空气压力保持在100 Pa左右。当弹丸碰撞靶板后，轻气炮高压气室中的氮气和一级段泵管中的氢气会尾随弹丸进入靶室。为最大限度地减小轻气炮发射气体对目标的影响，在二级发射管与靶室间放置了一个圆形档板，档板中心开一个直径20 mm的圆孔，一方面可减小轻气炮发射气体对目标动量的影响；另一方面可减少进入靶室的黑烟，影响高速摄像的拍摄。

1.4 动量传递系数的计算和误差分析

1.4.1 计算方法

从高速摄像记录弹道摆运动过程图中找出弹道摆的初始位置和最大位置处如图5所示，将2幅图片导入Adobe Photoshop软件，测量得到弹道摆最大水平位移xmax。根据弹道摆的最大水平位移xmax和几何关系可算出弹道摆的最大摆角θmax和弹道摆上升高度h。

图5 撞击误差分析Fig.5 Impact error analysis

根据能量守恒，计算弹丸撞击靶板后靶板获得的速度增量v为

为排除发射气体对动量传递的影响，开展相同工况下不装载弹丸发射试验，标定发射气流引起弹道摆水平位移xbd和最大摆角 θbd。弹丸引起弹道摆角度的实际变化 θ等于弹丸撞击弹道摆的最大摆角θmax减去不装载弹丸发射试验最大摆角 θbd，实际动量传递系数的大小为

其中：mp为弹丸质量；U为根据磁测速信号计算的弹丸速度；M为撞击后靶和吊绳的质量；g为重力加速度。

1.4.2 误差分析

在靶板质心建立Oxyz坐标系，如图6所示，坐标原点位于靶板质心，弹丸速度和x轴的夹角记为 α，撞击角度对动量传递系数带来的误差 ∆1为

图6 超高速撞击动量传递测试系统示意图Fig.6 Schematic figures of the hypervelocity impact momentum transfer test system

撞击角为10°时，撞击角对动量传递系数带来的误差1.6%。实际试验中弹道摆垂直悬挂在靶室中，弹丸出炮口后均沿水平方向运动。因此，弹丸均垂直撞击靶板表面α ≈0，因此撞击角度带来的误差约为0。

超高速撞击动量传递测试系统如图6所示，当撞击点在靶板表面中轴线z′上下移动，如图6（b）所示，由于撞击点改变对动量系数无影响，当撞击点在靶板表面中轴线左右移动时，撞击点距z′轴线距离为e，撞击后弹道摆向后运动过程中同时发生旋转，因此试验前需要进行激光对中。由图5可知，弹道摆未发生偏转，试验带来误差可忽略不计。

1.5 结果分析

不同撞击速度下成坑直径、深度和动量传递系数如表1所示。撞击速度从2.3 km/s增加到3.9 km/s，β从1.96增加到2.51。不同撞击速度下，弹丸撞击玄武岩靶成坑形貌基本一致。6 mm弹丸以3.9 km/s撞击玄武岩典型成坑结果如图7所示。玄武岩靶上形成的空腔由一个碗状中心小坑（图7中白色虚线圈内）和一片较大的浅层剥落区（图7中白色和红色虚线圈内）组成，与之前超高速撞击致密岩石成坑形貌一致[7]。碗状的中心坑是由冲击开挖流产生的，浅层剥落区是冲击波与稀疏波在自由表面的相互作用导致材料发生拉伸破坏形成的[15]。

表1 超高速撞击条件下的试验参数及结果Table 1 Experimental parameters in hypervelocity impact conditions and results

图7 撞击速度为3.9 km/s时，玄武岩的成坑结果Fig.7 Craters formed in basalt，with impact velocity of 3.9 km/s

2 数值模型

2.1 算法及模型

AUTODYN软件与其它仿真软件在计算陨石坑大小、撞击过程中峰值压力和动量增强系数大小进行了基准比较[16-18]，广泛用于模拟动能侵彻、空间碎片防护、陨石成坑等超高速撞击问题[19-22]。动能撞击偏转小行星涉及超高速和大变形的成坑过程，相对于基于网格的方法，SPH模拟超高速撞击成坑问题不但可以避免网格畸变，而且可以追踪在陨石坑形成过程中抛射出的粒子质量和速度。基于AUTODYN软件二次开发功能，可以导出所有SPH粒子在不同时间的空间位置、质量和速度信息。根据所有抛射物的质量和速度信息，将撞击方向上产生的所有抛射物粒子动量分量相加，动量传递系数β-1等于抛射物在撞击方向上的总动量除以抛射物初始动量[23]。

利用AUTODYN软件建立二维轴对称模型进行数值计算如图8所示，球形弹丸直径6 mm，采用航天器常用材料Al2024。靶厚度和直径均为123.3 mm的玄武岩，弹丸半径20个粒子，SPH粒子大小为0.15 mm[16]。

图8 简化计算模型Fig.8 Simplified calculation model

2.2 材料模型及参数

Al2024采用Shock状态方程、Johnson-Cook强度模型和最大拉应力失效，失效应力取2.6 GPa[24]，材料参数来自AUTODYN材料库，具体材料参数见表2和表3。

表2 Al2024和玄武岩的Shock状态方程参数[19，25]Table 2 Parameters of shock for Al 2024-T4 and basalt

表3 Al-2024 的Johnson-Cook本构模型和失效参数[26]Table 3 The Johnson-Cook model and failu parameters for Al-2024

玄武岩采用Shock状态方程，玄武岩采用Drucker-Prager强度模型，模型考虑了屈服强度Y与静水压力P之间的线性关系。

玄武岩初始屈服强度Y0取146 MPa，摩擦系数f取0.6[11]，玄武岩采用最大拉应力失效，失效应力取50 MPa[27]。

2.3 计算结果与试验结果比较

6 mm铝弹丸3.9 km/s速度撞击玄武岩的成坑结果如图9所示。从0.05～0.2 ms成坑深度不变，表明开坑完成。撞击0.2 ms后，层裂碎片与靶完全分离，速度不会再发生改变。层裂碎片速度较低，完全离开靶板表面需要较长，可根据导出的SPH粒子空间位置、质量和速度信息，计算反溅抛射物和层裂碎片的动量。因此，为节约计算时间，取不同撞击速度下撞击0.2 ms后，计算反溅抛射物和层裂碎片的动量。

图9 撞击后成坑结果，左0.05 ms，右0.2 msFig.9 Crater results，left 0.05 ms，right 0.2 ms

撞击后0.2 ms不同撞击速度形成的抛射物的质量–速度分布如图10所示。对抛射物质量–速度分布曲线积分，得到抛射物累积动量PI，动量传递系数β-1=PI/mU。水平横线上是层裂碎片的质量和动量，层裂大大提高了抛射物的质量，如图10（a）所示。但层裂碎片速度较低，贡献动量较小如图10（b）所示。在不同速度下，使用AUTODYN模拟预测的陨石坑直径和深度分别与试验相差4.7%～10.9%和3.7%～8%。数值模拟得到的动量传递系数与试验的最大误差为10.2%。通过与不同速度的试验结果对比，验证了动量传递系数的统计方法、所选算法和模型参数的合理性。在此基础上，模拟了5～10 km/s这4种不同撞击速度撞击玄武岩靶，研究撞击速度对抛射物质量、速度分布和动量传递系数的影响，结果如表4所示。随着撞击速度提高，弹丸和靶板可能会发生融化和气化，Shock状态方程不能直接计算得到温度或焓，只能间接估算，因此本文未考虑融化和气化对抛射物产生和质量–速度分布的影响。之前的数值模拟结果表明，在10 km/s以下，弹丸和靶板发生融化和气化的比例较低，Shock状态方程是适用的[27-28]。但在更高速度的情况下，需考虑熔化和气化问题，选择可以描述物质的凝聚态（固液）、膨胀态（气）和混合态的状态方程[29]，研究融化和气化产物对动量传递的影响。为防止高速下靶板被撞碎，将靶板尺寸扩大1倍。当陨石坑直径达到几米时，在陨石坑表面形成大面积的损伤，损伤区域的抛射物无法获得足够的速度离开小行星表面[30]。因此，将小型试验和模拟结果外推到小行星尺度需要修改试验和模拟结果，去除层裂机制产生的抛射物对动量传递的贡献。

表4 数值模拟与试验结果的比较Table 4 Comparison between numerical simulation and experimental results

图10 不同速度下的抛射物的质量和速度分布Fig.10 Mass-velocity distribution of ejecta at different velocities

2.4 数值计算结果及分析

2.4.1 抛射物的质量与速度分布规律

速度大于v的累积抛射物质量M（>v）分布曲线如图11所示，抛射物的速度逐渐变小，抛射物累计质量分布趋于平缓，接近抛射物质的总质量。不同撞击速度下累计抛射物质量M（>v）分布规律相似，撞击速度越大，抛射物的累计质量越大。撞击速度从5 km/s增加到10 km/s，β-1从1.64增加到2.72（见表5）。抛射物的质–速度分布函数可简化为幂律分布[10]，图11中黑色直线斜率等于-3µ，µ是耦合参数，黑色直线斜率等于–1.65。

图11 速度大于 v抛射物的质量 M(v)分布规律Fig.11 Distribution of ejecta with velocity over v

2.4.2 动量传递相似律

相似律模型广泛用于将试验室尺度的试验和模拟结果外延至行星级别，强度区间的动量相似律[10]为

其中：δ、U、Y、ρ分别为弹丸密度、撞击速度、靶强度和密度；ν为耦合参数；n=3µ；Ks、Kvs为材料相关的常数。

其中：m为弹丸质量；Me为速度大于v的抛射物质量；v∗为抛射物的最小速度。

速度大于v的抛射物质量如图12所示，玄武岩动量相似关系常数Ks=0.29，Kvs=0.47。将弹丸和靶板的参数带入到式（7），得到撞击密实小行星的动量传递曲线（U单位km/s）：

图12 速度大于 v的抛射物质量Fig.12 Mass with ejection velocity v

在10 km/s撞击速度范围内模拟的结果（图13）和相似律结果的最大误差在10%之间（见表5）。因此，在10 km/s以内的速度范围内，动量增强系数 β随撞击速度的0.65次方的增大而增大，相律模型可以很好地预测动量增强系数。

图13 密实小行星动量增强系数 β与撞击速度的关系Fig.13 Relationship between momentum enhancement factor β and impact velocity

根据威布尔理论，脆性岩石的强度Y随着其质量长度xmax的增大而减小，Y=Y0(L0/L)1/i，其中i通常为2～3[31]，靶板质量长度为L0时靶板强度大小为Y0。将实验室靶板质量长度L0外推至近地小天体大小，需要考虑岩石强度的变化。因此，考虑近地小天体尺寸的动量传递因子 β为

其中：β0为试验室条件下动量传递系数的大小，由式(11)可知，对于µ=0.55，i=2，β-1随近地小天体尺寸的0.16次方增大而增大。将试验室岩石的尺寸12 cm增加到140 m，近地小天体岩石的强度从146.0 MPa降低到4.2 MPa，则撞击140 m密实小行星的动量传递曲线（见图13）（U单位m·s–1）为

3 结 论

本文利用二级轻气炮进行了6 mm铝弹丸以2～4 km/s超高速撞击玄武岩靶试验，获得了在不同速度下撞击玄武岩的动量传递系数。将模拟结果和试验结果进行了比较，验证了计算模型和参数的准确性。模拟和相似律模型结果均表明，密实小行星的动量增强系数明显依赖撞击速度。撞击岩石目标数值模拟结果与相似律模型计算结果吻合较好，相似律模型可以很好地预测动量增强结果，在10 km/s以内的速度范围内，动量增强系数 β随撞击速度的0.65次方的增大而增大。