CL-20/HMX 共晶冲击起爆模拟*

2020-04-01李俊玲马兆侠

爆炸与冲击 2020年3期

刘海，李毅，李俊玲，马兆侠，陈鸿

（中国空气动力研究与发展中心超高速碰撞研究中心，四川绵阳 621000）

CL-20 是当前C/H/O/N 类含能材料中能量密度最高（约2.04 g/cm3）的炸药，氧平衡接近-11%，爆速和爆压分别为9.5 km/s 和43 GPa，但因较高的感度限制了在武器装备中的应用。2012 年，Bolton 等[1]成功制备了物质的量之比为2∶1 的CL-20/HMX 共晶，评估结果显示其爆速为9.484 km/s，优于β-HMX 和γ-CL-20，而感度低于多种晶型的CL-20 单晶。并且，CL-20/HMX 共晶燃烧时特征信号较低，因此在推进剂领域具有较大的应用潜力。2014 年，Sun 等[2]计算得到CL-20/HMX 共晶的多项物性参数，并解释了共晶感度低于ε-CL-20 单晶的原因。随后，Liu 等[3]分析了CL-20/HMX 在0～100 GPa 静水压条件下的物理、力学特性。高温热分解模拟显示CL-20/HMX 共晶及其单晶组分的初始分解路径均是N—N 键断裂。并且，热分解过程中CL-20 反应放热将传递至HMX，加速HMX 的分解[4]。Doblas 等[5]采用喷雾闪蒸技术制备了纳米尺度的CL-20/HMX 共晶，并研究了其结构特征以及热行为。对于单晶，单分子CL-20初始分解产物是稳定的芳香化合物1,5-二氢二咪唑[4,5-b:4′,5′-e]吡嗪[6]。Isayev 等[7]则给出了不同意见，其研究结果表明单分子CL-20 唯一初始反应路径是N—NO2键断裂。高温(2 000～3 500 K)条件下，CL-20 热分解的初始反应路径是N—NO2键断裂，最终反应产物是N2、CO2和H2O[8-9]。在冲击条件下，当冲击波速度为8、9 km/s 时，CL-20 单晶的初始分解路径是N—NO2键断裂；当速度提升至10、11 km/s时，H 原子脱落成为初始分解路径。N2、H2O、CO2为主要稳定产物[10]。当冲击波速度为8～10 km/s 时，HMX 冲击分解的主要初始产物是NO2，最终产物主要是N2、H2O、CO2、CO、H2、NH3[11]。而当冲击波速度提升至11 km/s 时，初始反应路径变化为C—H 键断裂[12-13]。

当前针对CL-20/HMX 共晶冲击行为的研究相对较少。直接冲击作用下CL-20/HMX 共晶冲击雨贡纽、爆炸参数、初始反应路径和稳定反应产物等物理、化学特征对于理解含能共晶起爆和爆轰特性具有现实意义。本文中采用非平衡分子动力学方法，同时结合在含能材料热及冲击行为研究中得到广泛应用的ReaxFF/lg 势[14]，模拟CL-20/HMX 共晶的冲击压缩过程，旨在获得不同冲击条件下含能材料中的冲击波传播时空特征和冲击诱导化学反应行为，揭示冲击起爆和爆炸参数，并与组份单晶开展对比分析，为共晶材料设计与工程应用提供数据支撑。

1 模型与模拟细节

物质的量之比为2∶1 的CL-20/HMX 单晶胞中含有4 个CL-20 分子、2 个HMX 分子[1]。构建20×4×4 的超晶胞，空间尺寸为32.69 nm×3.98 nm×4.86 nm。其中CL-20 和HMX 分子数量分别为1 280和640，共计64 000 个原子。首先通过不断迭代调整原子坐标以使系统能量最小化；随后采用等温等压NPT 系综将系统缓慢升至常温常压，此过程持续5 ps；接着将体系维持在常温常压下15 ps。温度和压力迭代的阻尼系数分别为300 fs 和1 000 fs，整个平衡过程采用周期边界条件，平衡后共晶的密度约为1.9 g/cm3。

原子分子系统的冲击模拟中，冲击加载方法主要有非平衡加载方法(动量镜[15-17]和活塞[18])，多尺度冲击技术(MSST)[19-20]以及分阶段压缩-膨胀方法[21]等。目前，主要使用前两类方法。活塞和动量镜加载方法的关键是在近真实碰撞条件下，针对较大模型体系，从波后粒子速度up出发，获得冲击波速度us以及沿空间路径的密度、温度等冲击波传播时空特征。而MSST 方法的核心是在冲击波关系式的约束下，针对较小模型体系施加冲击波条件us，获得沿时间路径的波后粒子速度up、密度、温度等动力学和热力学量变化。本文中，采用动量镜方法进行冲击加载。具体为CL-20/HMX 共晶块体材料以初始速度(粒子速度)up=-0.5～-4.0 km/s 穿过动量镜，动量镜反射原子动量以激发冲击波以速度us沿晶格矢量a 方向(x+方向)传播。沿冲击波传播方向采用收缩边界条件，垂直于冲击波传播方向采用周期边界条件。整个模拟过程采用的时间步长均为0.1 fs。CL-20/HMX单胞、超晶胞模型以及冲击波产生和传播原理如图1 所示。本文中，采用大尺度原子分子并行模拟器(large-scale atomic/molecular massivelyparallel simulator，LAMMPS)[22]开展计算。

图 1 CL-20/HMX 模型以及动量镜激发冲击波传播原理图Fig. 1 CL-20/HMX structure model and schematic diagram of the shock wave propagation induced by momentum mirror

2 结果与讨论

2.1 冲击波传播特征

沿冲击波传播方向对共晶材料进行一维分层，每层的厚度为0.5 nm，约含有1 000 个原子。对不同时刻每层中的密度和粒子速度进行统计分析形成x-t 图，如图2 所示。整体上，x-t 图主要分为三个区，Ⅰ区表示未压缩材料的初始状态，Ⅱ区表示冲击压缩区，Ⅰ、Ⅱ区界面表示冲击波阵面轨迹，而Ⅲ区表示冲击波到达右端自由面后反射形成的稀疏拉伸区。随着up的增大，冲击波到达右端自由面的时间逐渐缩短。当up=0.5 km/s，冲击压缩区密度最大为2.5 g/cm3。而当up提升至2.5 km/s 和4 km/s 时，冲击压缩区密度最大值分别为3.5 g/cm3和4 g/cm3。对于粒子速度，冲击波到达右端自由面后，反射形成的稀疏波进入含能材料并拉伸粒子向右运动，并且在拉伸初期，右端自由面粒子速度约为2up，随后粒子加速运动。另外，当up=0.5 km/s 时，Ⅰ和Ⅱ区界面有变形，这表明低速情况下，冲击波传播速度不稳定，而在高速情况下，冲击波传播速度较为稳定。

图 2 冲击波传播引发粒子速度和密度变化的x-t 图Fig. 2 x-t diagrams of particle velocity and density induced by shock wave propagation

2.2 冲击雨贡纽

冲击波在凝聚介质中传播时，在相当宽的速度范围内，冲击波速度us和波后粒子速度up之间具有如下形式的线性关系：us=sup+c0，其中s 表示绝热体积常数，c0表示材料声速。依据冲击波到达右端自由面的时间以及含能块体的长度得到冲击波速度，并绘制CL-20/HMX 的冲击波速度-粒子速度关系如图3（a）所示。空心矩形是本模拟结果，绿色线段是依据试验获得的材料声速(2.77 km/s) 和绝热体积常数(2.19)绘制[23-24]。而实心矩形、空心圆和实心圆则分别表示前期获得的CL-20 模拟结果[10]、HMX 冲击雨贡纽试验数据[25]和模拟获得的HMX 数据[11]。

从图中可以看出，当us＜9.38 km/s(对应up=3 km/s)时，冲击波速度-粒子速度基本呈线性关系，而当us＞9.38 km/s 时，冲击波速度-粒子速度发生偏折，此时含能材料爆轰后冲击波速度-粒子速度的间断特征与CL-20/HMX 共晶多尺度模拟结果[26]和动量镜冲击加载获得的RDX 模拟结果一致[15]。图中黑色线表示Bolton 等通过Cheetah 模拟获得CL-20/HMX 共晶爆速D=9.484 km/s[1]。含能材料us-up关系在发生爆轰后出现向下偏折，原因在于爆速是定值，即随着up的增加，us会略高于爆速，处于过驱爆轰状态，但随着传播而逐渐衰减为定常爆轰。对共晶未爆轰时的冲击波速度-粒子速度数据进行拟合得到us=1.724up+4.215。另外，同冲击波速度条件下，CL-20/HMX 共晶波后粒子速度小于多尺度模拟结果[26]，但整体上介于HMX 和CL-20 单晶之间。但与实验结果有一定的偏离。可能的原因是模拟对象是完美晶体，而真实含能材料多含有孔洞、杂质、颗粒边界条件等缺陷。图3（b）表示p-V/V0雨贡纽数据，实心矩形表示本模拟结果，空心圆和实心圆则分别表示前期获得的CL-20 模拟结果[10]。模拟结果显示同压力条件下，CL-20 的压缩程度高于CL-20/HMX 共晶。依据冲击波关系式以及冲击波速度-粒子速度关系，得到冲击压缩条件下凝聚介质p-V/V0雨贡纽的一般形式[27]：

式中：ρ0为材料初始密度。将拟合冲击波速度-粒子速度关系得到的c0和s 代入式（1）得到：

依据上述关系绘制理论解，如图3（b）中黑色曲线所示。模拟结果显示，理论解与数值模拟结果基本一致。而当共晶发生爆轰后，其p-V/V0数据发生偏离。据此计算得到CL-20/HMX 爆压约为45.6 GPa，介于多尺度模拟结果获得的36.75～47.43 GPa 之间[26]。

图 3 冲击波速度-粒子速度关系及p-V/V0 雨贡纽Fig. 3 Shock wave velocity-particle velocity relation and p-V/V0 Hugoniot

2.3 反应物衰减速率与产物分布

图4 表示共晶中CL-20 和HMX 分子数量的衰减曲线。整体上，随着up的增加，共晶中CL-20 和HMX 分子数量的衰减速率加快，并且含能分子冲击分解前存在诱导期。当up=0.5，1 km/s 时，共晶并未发生分解。当up提升至1.5 km/s 时，共晶中仅CL-20 发生分解。up≥2 km/s 后，共晶中CL-20 和HMX 均发生分解，且分子数量基本呈线性衰减。在冲击波到达右端自由面后，反应物分子数量基本不变。这表明冲击波反射形成的稀疏波拉伸材料不会导致反应物的分解。

图 4 共晶中CL-20 和HMX 分子衰减曲线Fig. 4 Decay curves of CL-20 and HMX molecules in co-crystal

图5（a）中，左y 轴表示CL-20 和HMX 诱导期τ，该诱导期表示冲击加载后至激发释能化学反应的时间[28-29]，右y 轴表示对冲击波到达右端自由面时，共晶中CL-20 和HMX 分子数量衰减曲线进行归一化后得到的平均衰减速率（decay rate）rd。整体上CL-20 和HMX 诱导期均随着up的增加而缩短。在各冲击条件下，HMX 的诱导时间均多于CL-20，且随着冲击速度的增加，CL-20 和HMX 的诱导时间逐渐接近。通过对两组数据进行比较，获得各冲击条件下CL-20 和HMX 衰减速率的接近程度，如表1 所示。CL-20 和HMX 分子数量的衰减速率随着up的增加而增加，并且各冲击速度条件下，共晶中CL-20 分子数量的衰减速率均大于HMX，但随着冲击波速度的增加，CL-20 和HMX 的衰减速率逐渐接近。这主要是由于含能材料爆轰后，CL-20 反应放热传递至HMX，增加了HMX 的反应速率，同时使得CL-20 分解速率降低，进而导致CL-20 和HMX 分子数量的衰减速率逐渐接近。在CL-20/HMX 热分解[4]与CL-20/TNT[30]冲击诱导分解中均发现热传递现象。

图5（b）表示冲击波到达右端自由面时，共晶分解产生碎片的种类及其数量分布。随着up的增加，反应产物的类别和数量均逐渐增多。当up=1.5 km/s 时，共晶开始分解并产生NO2，对比图4 可知，此NO2来自于CL-20 分子中N-NO2键断裂。随着冲击波速度的增加，系统中产物类别逐渐增多。up=2.5 km/s，除NO2外，其他中间小分子产物NO、OH、HONO 出现，并且最终稳定产物N2形成。当up增至4 km/s时，系统中小分子产物NO2、NO、OH、HONO 和H 中，NO2和OH 的数量较多。而稳定产物为N2、CO2和H2O。并且，N2分子数量较多，H2O 分子数量居其次，CO2数量最少。另外，NO2是共晶在各冲击条件下分解的主要中间产物。CL-20 和HMX 单晶模拟结果显示，较高冲击速度条件下，N—NO2键断裂受到抑制，H 原子首先脱落[6-7]。对比多尺度模拟结果[26]发现，二者在初始反应路径、中间小分子产物类别和分布特征上均一致，但在最终反应产物上，多尺度模拟结果显示最终稳定产物为N2、CO2、H2O、CO 和H2。这主要是由于本模拟仅针对冲击压缩阶段开展分析，而多尺度模拟结果涵盖了冲击波扫过共晶材料后的压缩阶段以及部分化学反应区(膨胀阶段)中的稳定产物。因此可以推测本模拟中，压缩波到达右端自由面反射拉伸材料，体积膨胀后，压缩阶段形成的大分子团簇分解，CO 和H2分子将出现。在Li[17]的模拟结果中，证明了RDX 冲击压缩后拉伸膨胀阶段的新产物形成情况。

图 5 诱导时间和衰减速率及冲击到达右端自由面时主要产物的数量统计Fig. 5 Induction time and decay rate under different shock conditions and quantitative statistics ofthe main products when the shock wave reaches the free surface

表 1 各冲击条件下共晶中CL-20 和HMX 诱导时间和衰减速率比较Table 1 Comparison of induction time and decay rate of CL-20 and HMX in co-crystal under different shock conditions

当up=1.5 km/s (us=6.8 km/s) 时，冲击诱导共晶分解。采用兰金-雨贡纽关系P=ρ0usup，其中us、up分别表示冲击波激发共晶分解的冲击波速度和对应的粒子速度，计算得到冲击起爆压力为19.38 GPa，与多尺度模拟方法获得的结果(16.52 GPa)[26]基本一致，低于HMX 单晶的冲击起爆压力35.1 GPa[11]。

图6 直观展示了up=1.5 km/s 和4.0 km/s 时冲击波从左向右传播的序列图像。同一时刻，up=4.0 km/s对应的冲击波传播距离大于up=1.5 km/s 的情况。当up=4.0 km/s 时，在t=3.15 ps 时刻，含能材料完全压缩。而后反射形成的稀疏波入射至材料中，自右向左传播并拉伸已压缩材料，如图中t=3.5，4.0 ps 对应的图像所示。而当up=1.5 km/s，在4 ps 内，冲击波未到达右端自由面。另外，通过压缩区的疏密程度可以直观反映出up=4.0 km/s 对应的冲击压缩区的密度大于up=1.5 km/s 时的压缩程度。