高温下气相CL-20热分解反应的分子动力学计算

张 力， 李朝伟

(山西警察学院 治安系，太原 030401)

1 引 言

含能材料在实际生产，使用，存储和运输过程中，受到各种外界刺激发生相变，不仅导致其感度和安全性发生变化，而且影响了极端条件下含能材料反应过程. 虽然人们采用仪器分析[1]、金刚石压腔[1]和冲击加载[1]以及烤燃[2]等研究含能材料在静态与冲击加载下的相变过程，但却很难直接获得不同相态含能材料的详细反应路径. 近年来，van Duin等[3]提出反应力场分子动力学计算方法(ReaxFF reactive molecular dynamics simulation，简称RMDS)，为从原子分子尺度研究含能材料的反应机理提供了新思路. 自2003年Strachan等[4]采用RMDS，研究冲击作用下RDX引发的化学反应以来，人们已经使用这种计算方法，对HMX[5]、RDX[6]、TATB[5]、PETN[7]等在高温或者冲击作用下的产物、反应路径及化学反应动力学过程进行大量研究.

2011年，Liu等[8]在van Duin等[4]开发的反应力场基础上，进一步考虑伦敦色散作用，发展了ReaxFF-lg力场，随后人们采用这种力场对包括硝基甲烷[9]、TNT[10]、HMX[11]等反应过程进行了模拟. 研究均表明密度显著影响了固相含能材料的反应路径及热分解产物的空间分布. 近期人们[12-15]采用RMDS研究了密度对固相CL-20和CL-20/TNT[16]共晶炸药的分解反应影响，结果表明密度降低了影响了其反应过程. 此外，Xue等[17]采用SCCF-DFTB计算方法对CL-20冲击作用下的反应过程进行计算，计算结果表明在弱冲击作用下，CL-20分子容易发生N-NO2键断裂生成硝基自由基，而在强冲击作用下，容易发生H原子的脱落反应. 但鲜有研究密度对非凝聚相CL-20热分解反应的影响.

本文采用RMDS，基于ReaxFF力场研究密度和温度对气相CL-20的势能、碎片数量随时间变化的影响，得到了中间产物和最终产物的基元反应路径以及反应速率常数，为认识极端条件下非凝聚相含能材料的冲击起爆、冲击点火和爆轰等瞬态过程奠定基础.

2 计算方法和模型

采用ADF软件的Reax Ams模块，建立了五种由48个ε-CL-20分子组成的气相CL-20结构，其晶格参数分别为70.4 Å × 70.4 Å × 70.4 Å，48.8 Å × 48.8 Å × 48.8 Å，41.2 Å × 41.2 Å × 41.2 Å，36.8 Å × 36.8 Å × 36.8 Å，33.9 Å × 33.9 Å × 33.9 Å，α、β和γ角均为90°，对应的气相CL-20体系密度分别为0.1 g·cm-3，0.3 g·cm-3，0.5 g·cm-3,0.7 g·cm-3,0.9 g·cm-3. 图1为0.1 g·cm-3时的计算模型，图1a和图1b分别为分子结构和超晶胞结构，灰色、蓝色、红色和白色原子分别表示碳、氮、氧和氢原子.

图1 ε相CL-20的分子结构(a)和低密度状态超晶胞结构(b)Fig. 1 The molecule (a) and supercell structure (b) of ε phase CL-20

首先采用牛顿梯度法进行能量最小化进行结构优化，收敛常数为10-7(kcal/mol)/Å. 然后基于ReaxFF力场和NVT系综，采用Nose-Hoover控温方法，控温参数为10 fs. RMDS采用两步模拟，第一步是常温下的平衡过程，对五种结构进行时间步长为0.1 fs，温度为298 K，共10 ps的平衡过程；第二阶段是热分解过程计算，在1500 K，2000 K，2500 K和3000 K下，采用时间步长为0.1 fs，共490 ps的RDMS. 根据势能变化曲线和热分解碎片数量在20 ps内的波动小于0.6%，即表明气相CL-20体系达到平衡状态. 采用一级反应速率方程拟合稳定分解产物随时间的变化过程.

C(t)=C∞{1-exp[-k3(t-ti)]}

(1)

式中，C∞为平衡状态时不同产物的数量，k3为热分解产物的生成速率常数，ti为产物开始出现的时间.

3 结果与讨论

3.1 势能和碎片数量随时间变化

图2是3000 K时，势能随时间和密度的变化曲线. 从图中可以看出，势能随着时间增加而上升，随后达到最大值，将此时所处的状态称为点火状态. 所有体系均在3000 K时发生点火，然后按照不同能量释放规律继续发生反应. 0.1 g·cm-3，0.3 g·cm-3，0.5 g·cm-3，0.7 g·cm-3和0.9 g·cm-3体系，分别在31.9 ps，17.1 ps，13.35 ps，12.45 ps，11.47 ps时发生点火，表明在气相状态下，CL-20炸药发生点火所需要时间随密度的增加而降低.

这是由于在低密度状态下，密度未显著影响气相CL-20体系结构，此时体系主要由分子之间长程相互作用占据主导作用. 随着密度的增加，长程相互作用能量增加，促进了CL-20分子的热分解反应，降低点火时间. 此外，在0.1 g·cm-3和0.3 g·cm-3时，势能尚未达到平衡状态，当密度增加到0.5 g·cm-3时，体系的势能达到平衡状态.

图2 3000 K时势能随密度和时间的变化曲线Fig.2 Time evolutions of potential energy of different densities at 3000 K

图3为3000 K时分解产物碎片数量随时间的变化曲线. 其变化过程分为三个阶段，第一阶段为第0 ps-31.7 ps，碎片数量随着CL-20发生分解反应而迅速增加，从最初的48个分子增长到500个左右的碎片，这是由于CL-20体系发生初始分解反应所形成的；第二阶段为31.7 ps-115.9 ps，碎片数量缓慢增加到650个左右的碎片，碎片数量的缓慢增加是由于CL-20初始分解产物继续发生反应所导致的；第三阶段时为115.9 ps-500 ps，此时CL-20热分解反应基本达到平衡状态，其碎片数量不会发生明显变化，在750个左右的平衡位置处波动.

图3 碎片数量随密度和时间的变化曲线Fig.3 Time evolutions of total number fragments of different densities at 3000 K

本文采用Döntgen等[19]提出的RMDS自动发现反应路径、生成物、分解产物和过渡态搜索方法，对不同体系的反应物/产物进行分析. 在3000 K和第10 ps-20 ps时，0.1 g·cm-3、0.3 g·cm-3、0.5 g·cm-3、0.7 g·cm-3、0.9g·cm-3时，反应物/产物的种类分别为184/181、231/261、250/274、320/373、396/443，表明低密度状态下气相CL-20发生分解反应生成的产物种类随着密度的增加而增大. 表1为气相CL-20在不同密度下主要反应物/产物及反应路径. 在不同密度下，NO2、HONO、HCN、NO等在反应物和产物中大量出现，但其反应路径却有所差异，以HONO为例，NO/HO、NO2/HCN、NO2/H均能发生反应生成HONO.

表1 3000 K时主要反应物和产物及反应路径

3.2 初始反应路径及产物

图4为不同温度下，CL-20超晶胞发生分解反应，生成硝基自由基的数量随密度和时间的变化曲线. 在图4d中的3000 K时，0.1 g·cm-3，0.3 g·cm-3，0.5g·cm-3，0.7 g·cm-3,0.9 g·cm-3的体系中，硝基自由基的峰值个数/所处时刻分别为240/18.95 ps，237/14.5 ps，205/12.5 ps，215/12.95 ps，204/12.0 ps，表明平均每摩尔CL-20生成硝基自由基数量分别为5 mol，4.94 mol，4.27 mol，4.48 mol和4.25 mol，由于CL-20分子包括2个与六元环联结的硝基官能团和4个与五元环联结的官能团，因此CL-20分子中的五元环和六元环中的硝基官能团均发生断裂，二者在反应过程中没有明显的区别.

图4 硝基自由基的数量随密度、温度和时间的变化曲线Fig.4 Time evolutions of nitro radical fragments of different densities at various temperatures

2019年，Liu等[20]采用第一性原理分子动力学计算方法对ε相CL-20热分解反应进行计算，结果表明在3000 K时，平均每摩尔CL-20超晶胞结构有2个N-NO2键发生断裂生成硝基自由基. Okovytyy等[21]采用密度泛函理论对单分子CL-20的分解路径进行计算，表明六个不同的硝基官能团均能发生N-NO2键断裂生成硝基自由基，Wang等[14]计算结果表明固相CL-20中，72个CL-20分子生成197个硝基自由基. 与固相CL-20生成硝基自由基的数量相比，由于空间位阻效应的减小，使得低密度状态下体系生成硝基自由基的数量远大于固相CL-20发生分解反应生成硝基自由基数量.

Isayev等[22]采用第一性原理分子动力学计算方法，张力等[15]和Wang等[14]采用RSMD对高温下CL-20的反应路径进行研究，结果均表明其初始反应路径为N-NO2键断裂生成硝基自由基，且在高温下六个硝基官能团均能断裂生成硝基自由基. Liu等[20]的计算结果表明高温下固相CL-20初始分解反应的初始反应路径为五元环和六元环中的N-N键、C-N键和C-N键的断裂生成硝基自由基和C6H6N11O10，且固相CL-20在高温下的热分解反应为单分子反应路径. 这与本文计算得到的低密度状态下CL-20初始反应路径相同. 陈芳等[23]采用RSMD对高温下RDX的分解反应过程进行计算，结果表明RDX的初始反应路径中也存在N-NO2键断裂生成硝基自由基的过程.

3.3 中间反应路径及重要的中间产物

图5为0.9 g·cm-3时，气相CL-20在不同温度下的中间产物数量随时间的变化曲线，存在HNO3、NO3、N2O2、NO、HONO、CHN、CHN2、CHN2O2等中间产物. 在3000 K时，HONO、NO、CHN、N2O2等产物数量均出现了峰值，然后逐渐降低. Liu等[20]认为在CL-20分解反应中存在N2O5、N2O2、NO2、NO、HONO、HCN等产物，这与本文计算得到的低密度状态下CL-20分解反应的中间体及分解产物一致.

图5 0.9 g·cm-3的CL-20在不同温度下的中间产物随时间的变化曲线Fig.5 Time evolutions of 0.9 g·cm-3 CL-20 intermediate fragments at different temperatures

表2为3000 K时，N2O2和NO涉及的基元反应路径. 存在7种中间分解产物直接生成N2O2，N2O2与10种其他中间产物发生反应生成其他碎片的反应路径；NO的基元反应中存在6种中间体生成NO，NO与1种中间体生成其他产物的反应路径. 其中N2O2涉及到5个基元反应，硝基自由基参与四个基元反应，伴随着N-N键的断裂和生成以及C-H键断裂和N-H键的生成.

表2 3000 K时N2O2和NO的基元反应路径

表3为3000 K时，HCN和HONO的基元反应路径. 存在12种中间分解产物直接生成HCN，HCN与11种其他中间产物发生反应生成其他碎片的反应路径，涉及到10个基元反应；存在14种中间分解产物直接生成HONO，HONO与3种其他中间产物发生反应生成其他碎片的反应路径，涉及到8个基元反应.

表3 3000 K时HCN和HONO的基元反应路径

3.4 化学反应速率常数拟合

图6分别为3000 K时，气相CL-20最终分解产物CO、N2、CO2、H2O的数量随时间和密度的变化曲线. 在密度为0.1 g·cm-3时，CO的数量并没有出现达到平衡状态，而是持续增加，在密度为0.3 g·cm-3时，CO的数量在200 ps时达到平衡状态，在密度为0.5-0.9 g·cm-3时，CO的数量随时间的变化趋势相同，均为在不同时刻达到平衡.

表4为第500 ps和3000 K时，最终稳定产物的数量. N2、H2O、CO2的数量随着密度的增加而增大，在密度超过0.3 g·cm-3时N2的数量为270个左右，在密度超过0.5 g·cm-3时H2O的数量为130个左右，在密度超过0.5 g·cm-3时H2O的数量为140个左右；CO、H2、NO的数量随着密度的增加而降低，当密度超过0.5 g·cm-3时CO的数量为130个左右，当密度超过0.5 g·cm-3时NO的数量迅速降低到5个. Strachan等[6]采用RMDS表明在RDX的分解过程中，密度显著影响了CO和CO2的分子生成的数量，且CO2更容易在高密度状态下生成，而CO更容易在低密度状态下生成. 陈芳和程新路[24]采用RMDS，对冲击作用下RDX的分解产物进行研究，发现在不同冲击作用下RDX的分解产物包括N2，CO和CO2分子.

表4 3000 K时不同密度的体系最终生成产物的数量

图7为采用式(1)拟合得到3000 K时，最终稳定产物速率常数k3随密度的变化. 在五种密度下，k3(H2O)的值均大于k3(CO2)、k3(N2)和k3(CO).在0.1 g·cm-3时，k3(CO2)、k3(N2)和k3(CO)的值接近，当密度大于0.3 g·cm-3时，CO、N2和CO2的生产速率常数的大小顺序为k3(CO)>k3(N2)>k3(CO2)，表明密度未影响不同密度下H2O的生成速率常数，但影响CO、N2和CO2的反应速率常数. Wang等[14]基于ReaxFF-lg力场对密度高于1.92 g·cm-3时和3000 K时，最终稳定分解产物的反应速率常数进行拟合，为k3(H2O)>k3(N2)>k3(CO2)，这也与本文计算得到的结果一致.

图6 3000 K时不同密度下最终产物随时间的变化曲线Fig. 6 Final products of different densities at 3000 K

表5 3000 K时N2和H2O的基元反应路径

图7 3000 K时主要产物的生成速率常数k3随密度的变化Fig. 7 Final products formation rate constants of different densities at 3000 K

4 结 论

通过采用RMDS，对气相CL-20体系在高温下的分解反应进行研究，结果表明密度和温度影响气相CL-20热分解的反应路径和化学反应速率常数，主要结论如下：

一是分解产物碎片数量随时间分为三个阶段. 第一阶段为CL-20分子分解生成初始产物碎片，当密度高于0.3 g·cm-3时，体系发生点火导致碎片数量迅速增加；第二阶段是初始产物碎片继续发生反应，此时0.1 g·cm-3的气相CL-20也发生点火，其他密度下气相CL-20势能降低，碎片数量继续增加；第三阶段是在不同的温度下势能和碎片数量在平衡位置波动.

二是密度和温度没有影响第一阶段气相CL-20体系的初始反应路径，但密度影响了CL-20分子发生N-NO2键断裂生成硝基自由基的数量和时间；在第二阶段，气相CL-20体系的分解反应为C-C键、C-N键和N-N键的断裂与生成导致的质子转移、开环、闭环反应，形成了包括N2O2、NO、HONO和HCN等中间产物；第三阶段为N2，H2O，CO2等最终稳定产物的生成过程，且中间产物的基元反应的种类大于最终稳定产物的基元反应种类.

三是密度和温度对气相CL-20体系的N2，H2O，CO2等分解产物数量及反应动力学过程的影响不同. 当密度高于0.3 g·cm-3时，密度未影响N2、H2O、CO2等产物的数量，但却影响了CO的数量. 当密度高于0.3 g·cm-3时，N2，H2O，CO2，CO等产物的反应速率常数随着密度和温度的增加而增大.