孔洞缺陷对B炸药性能影响的理论计算

2019-09-17王玉玲

原子与分子物理学报 2019年4期

苗爽，王涛，王玉玲

(火箭军工程大学，西安 710025)

1 引言

B炸药是一种典型的熔铸炸药，由黑索金(RDX)和梯恩梯(TNT)按照一定的质量比例熔铸而成.B炸药兼顾了RDX的高能量和TNT良好的熔铸成型的特点，具有良好的爆轰性能和热安定性，是常规武器广泛采用的主装药.国内对B炸药开展了广泛的研究以改进B炸药的性能.牛国涛等[1]人研究了老化对TNT/RDX熔铸炸药性能的影响.聂少云[2]和贺传兰[3]等分别研究了不同聚合物添加剂对B炸药安全和爆轰性能的影响.但是对于孔洞缺陷对B炸药性能影响的研究比较少.为进一步深化对B炸药性能的了解，研究孔洞缺陷对B炸药性能影响是有必要的.

分子动力学方法广泛应用于分析含能材料的感度、力学性能和分子结构等性能.杭贵云[4, 5]等人成功采用分子动力学的方法研究了RDX和TNT不同质量配比以及温度对B炸药的力学性能和感度的影响.本文，建立了空位缺陷浓度为15%，孔洞大小不同的3种B炸药缺陷模型，采用分子动力学方法探究孔洞缺陷及其大小对B炸药感度、结合能、爆轰性能和力学性能的影响.

2 计算模型与计算方法

2.1 模型的建立

B炸药是熔铸炸药，属于非晶体结构，没有特定的晶格参数[6].因此，本文利用Material Studio(MS)软件中的Visualizer模块，采用无定形晶胞建模方法，根据65:35的质量分数比，建立了RDX和TNT分子数分别为114和60的“完美”型B炸药模型，标记为Model-Ⅰ.在“完美”型B炸药模型的基础上，通过剪切掉17个RDX分子和9个TNT分子建立B炸药第一种孔洞模型，标记为Model-Ⅱ；通过分别剪切掉10个RDX分子和5个TNT分子以及7个RDX分子和4个TNT分子建立含2个孔洞的第二种B炸药缺陷模型，标记为Model-Ⅲ；通过分别剪切掉14个RDX分子和7个TNT分子以及3个RDX分子和2个TNT分子建立含2个孔洞的第三种B炸药缺陷模型，标记为Model-Ⅳ.RDX、TNT、分子如图1所示，B炸药“完美”晶胞模型和缺陷晶胞模型如图2所示.

模型建立之后，采用Smart算法对初始模型

进行能量最小化，优化收敛精度设置为0.001 kcal/mol/Å.而后，对优化后的晶胞模型进行了模拟退火，退火精度设置为“Fine”.退火后，不同B炸药晶胞模型的晶格参数如表1所示.

图1 B炸药的分子模型Fig.1 Molecular models of Composition B

图2 B炸药及其缺陷模型Fig.2 The defective models of Composition B

Table 1 The lattice parameters of different Composition B models

lattice parametersa(Å)b(Å)c(Å)α(°)β(°)γ(°)Model-Ⅰ33.633.433.889.390.490.0Model-Ⅱ33.532.933.189.792.088.7Model-Ⅲ31.733.632.384.788.391.1Model-Ⅳ33.032.932.990.189.689.7

2.2 分子动力学模拟

Compass[7]力场对凝聚相材料性能的仿真计算具有出色表现，全部的B炸药模型均是在Compass力场中进行的分子动力学(MD)模拟.

MD模拟选择在恒压恒温(NPT)系综下进行，模拟的温度设置为298 K，压力设置为1.01×10-4GPa，精度设置为“fine”.模拟过程中分采用Nose[8]和Berendsen[9]方法对温度和压力进行控制，范德华(vdW)和静电作用(Coulomb)的加和方法分别采用Atom-based[10]和Ewald[11]，截断半径取15.5 Å，并对截断尾部进行矫正.原子运动的初始速度由Boltzmann分布确定，并采用Verlet方法求解牛顿运动方程的积分.对所建模型进行了100 ps的MD模拟，前30 ps对模型体系进行平衡，后70 ps用于统计分析确定能量、力学参数和其它参数.每0.1 ps取样一次，共获得700帧轨迹.

3 结果与讨论

3.1 平衡判别和平衡结构

对原子运动轨迹的统计分析必须基于模拟体系达到平衡才有意义.模拟体系平衡的标志是温度和能量达到平衡，其判断标准是温度和能量随时间变化的波动幅度在5～10%.以B炸药“完美”晶胞模型的MD模拟为例，其温度、能量随时间变化的曲线图如图3所示.由图3可见，B炸药“完美”晶胞模型的温度和能量随时间的推移逐渐趋于平缓，温度和能量偏差较小，表明B炸药“完美”晶胞模型已达到平衡状态.

3.2 感度

在含能材料领域，感度是指含能材料在受到外界激励的作用下发生分解甚至爆炸的难易程度.作为含能材料安全性评估的重要内容，感度是影响炸药应用进程的重要环节之一.许多采用MD方法研究含能材料性能的研究文献[12-17]提出了评估含能材料感度的理论，其中杭贵云等人[4, 5]利用键连双原子作用能和内聚能密度作为判据成功对B炸药的感度进行了预测.因此，本文采用键连双原子作用能和内聚能密度作为判据，评估孔洞缺陷对B炸药感度性能的影响.

3.2.1引发键键连双原子作用能

引发键是指含能材料在受到外界刺激时，最先发生解离导致含能材料分解甚至爆炸的化学键.B炸药模型中，RDX的感度高于TNT且研究表明RDX的引发键为N-NO2基团中的N-N键[18].因此，选择N-N键作为B炸药及其孔洞模型的引发键.

B炸药的引发键的键连双原子作用能的定义为：

(a) 温度变化曲线

(b) 能量变化曲线图3 温度和能量随时间的变化曲线Fig.3 Temperature and energy curves versus time

(1)

式中，E1是指平衡状态下，炸药模型体系的总能量；E2是指平衡状态下，炸药模型体系内所有RDX分子中N原子被固定后的总能量；n是指炸药模型内所有N-N键的数量.

经MD模拟计算，得到了不同模型的引发键键连双原子作用能，如图4所示.从图4可以看出含孔洞缺陷的B炸药的引发键键连双原子作用能较“完美”型B炸药减小.其中，Model-Ⅰ下降得最大，Model-Ⅳ次之，Model-Ⅲ下降得最小.这表明孔洞缺陷会导致B炸药的引发键键连双原子作用能下降，并且孔洞越大，下降程度越大.引发键键连双原子作用能下降导致引发键对N原子的束缚能力下降，炸药的感度上升，安全性降低.

3.2.2内聚能密度

内聚能密度(CED)是指单位体积内1 mol凝聚态物质汽化为气态物质所需要的能量.计算公式为：

图4 不同模型的引发键键连双原子作用能Fig.4 The interaction energies of trigger bond of different models

(2)

式中，ΔHV是指摩尔蒸发热；RT是指汽化时所做的膨胀功；Vm是指摩尔体积.

内聚能密度是分子间相互作用力的综合反映，与含能材料的感度之间存在一定关联，含能材料的CED值越小，表明含能材料受热发生相变所需的能量越小，其对热越敏感.通过仿真计算得到了不同B炸药模型的CED值，如表2所示.

表2 不同B炸药模型的CED及其分量

Table 2 CED and its components of different Composition B models

Model(n)CED/kJ·cm-3vdW/kJ·cm-3Electrostatic/kJ·cm-3Ⅰ0.7210.3280.384Ⅱ0.6960.3170.370Ⅲ0.7130.3230.382Ⅳ0.7060.3220.376

从表2中可以看出孔洞缺陷的存在导致B炸药的内聚能密度下降.其中，Model-Ⅰ下降最大，Model-Ⅳ次之，Model-Ⅲ下降最小.这表明孔洞缺陷会导致B炸药的内聚能密度下降，并且孔洞越大，下降程度越大.内聚能密度下降意味着B炸药发生相变所需能量下降，感度上升，安全性降低.

利用引发键键连双原子作用能和内聚能密度为判据评估孔洞缺陷及其大小对B炸药感度性能的影响获得的结论一致.我们推测孔洞缺陷导致感度升高的原因是缺陷B炸药的自由体积增加，提高了炸药分子的活化程度，稳定性下降，感度增加，安全性降低.

3.3 结合能

熔铸炸药的稳定性可以通过各组分之间的结合能(Eb)来表征.结合能越大，说明熔铸炸药各组分之间的相互作用力越大，熔铸炸药的稳定性越好.

B炸药及其缺陷模型的结合能定义为：

Eb=-Einter=-[Etotal-(ERDX+ETNT)]

(3)

式中，Eb是结合能；Einter是B炸药中RDX与TNT之间的相互作用能；Etotal是平衡状态下，B炸药体系的总能量；ERDX是平衡状态下，B炸药体系中所有RDX的能量；ETNT是平衡状态下，B炸药体系中TNT的能量.经计算，得到了B炸药及其缺陷模型的结合能，如图5所示.

图5 不同模型的结合能Fig.5 The bind energies of different models

从图5中可以看出，含孔洞缺陷的B炸药模型的结合能较“完美”型B炸药模型急剧下降，孔洞缺陷模型之间的结合能变化幅度不大.其中，Model-Ⅰ下降得最大，Model-Ⅳ次之，Model-Ⅲ下降得最小.这表明孔洞缺陷的存在导致B炸药的稳定性下降，但是孔洞尺寸对其影响较为微弱.

3.4 爆轰性能

爆轰性能主要反映含能材料的威力与能量密度，主要通过氧平衡系数(OB)、爆速(D)、爆压(p)和爆热(Q)等参数进行表征.本文采用修正氮当量法[19]计算了炸药的爆速和爆压，采用质量加权法计算了炸药的爆热，以此对其能量特性进行评估.

对于C-H-O-N四种元素组成，化学式为CaHbOcNd的炸药，氧平衡系数计算公式如下：

(4)

式中，a、b、c分别为炸药分子中包含的C、H、O三种原子的数目；Mr为炸药的摩尔质量，g·mol-1.

根据修正氮当量理论，爆速(D)和爆压(p)的计算公式如下：

(5)

式中，D为炸药的爆速，m·s-1；p为炸药的爆压，GPa,；ρ为炸药的密度，g·cm-3；ΣNch为炸药的修正氮当量；pi为1 mol炸药爆炸时生成第i种爆轰产物的摩尔数；Npi为第i种爆轰产物的氮当量系数；BK为炸药分子中第K种化学键出现的次数；NBK为炸药分子中第K种化学键的氮当量系数；Gj为炸药分子中第j种基团出现的次数；NGj为炸药分子中第j种基团的氮当量系数.

关于修正氮当量法的具体计算过程以及pi，Npi，BK，NBK，Gj，NGj等参数的来源可以通过文献[20]获得详细信息.

混合炸药可以通过质量加权法计算得到爆热，其计算公式为：

Q=∑ωiQi

(6)

式中，ωi是混合炸药中第i种组分的质量分数；Qi是混合炸药中第i种组分的爆热，kJ·kg-1.根据文献[21]，RDX和TNT的爆热值分别为6025 kJ·kg-1和5066 kJ·kg-1.

通过MD模拟获得了不同B炸药模型的密度(ρ)，并根据公式(4)～(6)计算得到了氧平衡、爆速、爆压和爆热等爆轰参数，如表3所示.

根据计算公式(4)～(6)，可以看出影响爆速和爆压有两个因素，密度和B炸药中各组分的组成；氧平衡和爆热主要受B炸药中各组分组成的影响.从表3可以看出，孔洞缺陷对B炸药的爆轰性能影响较小，其中氧平衡和爆热的变化几乎可以忽略不计.受密度减小的影响，爆速(D)和爆压(p)略微减小，毁伤能力下降.分析B炸药孔洞缺陷模型MD模拟的轨迹文件发现，孔洞缺陷对B炸药密度影响微弱的原因在于随着MD模拟时间的推移，B炸药的孔洞会逐渐塌陷收缩.

表3 不同PETN超晶胞模型的爆轰参数

Table 3 The detonation properties of different PETN supercell models

Model(n)Model-ⅠModel-ⅡModel-ⅢModel-ⅣOB(%)-40.0-39.9-39.9-39.9ρ(g·cm-3)1.6091.6051.5991.599D(m·s-1)7676.487661.787642.637642.63p(GPa)24.9824.8424.6624.66Q(kJ·kg-1)5689.355689.735689.735689.73

3.5 力学性能

炸药的力学参数主要包括弹性模量(E)、体积模量(K)、剪切模量(G)、体积模量与剪切模量的比值(K/G)以及柯西压(C12-C44)，分别表征炸药的抗变形能力、抗断裂强度、硬度和延展性.通过分析计算MD模拟后的平衡体系的原子轨迹获得炸药的弹性系数.根据广义虎克定律[22]，通过最小二乘法拟合弹性系数得出平均的拉伸应力应变，获得体积模量(K)和剪切模量(G).虎克定律、体积模量和剪切模量的计算公式如下所示：

σi=Cijεj

(7)

KR=[S11+S22+S33+2(S12+S23+S31)]-1

(8)

GR=15[4(S11+S22+S33)-4(S12+S23+S31)
+3(S44+S55+S66)]-1

(9)

式中，σ表示应力；ε表示应变；Cij(i，j=1，2，……，6)表示弹性系数矩阵；Sij表示柔量系数矩阵，等于Cij的逆矩阵，即S=C-1.

力学参数之间具有相互联系，关系式如下所示：

E=2G(1+μ)=3K(1-2μ)

(10)

根据公式(10)可以计算得到弹性模量(E)以及泊松比(μ)

(11)

(12)

从而求得不同B炸药模型的力学参数，如表4和图6所示.

从表4和图6中可以看出，孔洞缺陷导致B炸药的弹性模量(E)、体积模量(K)、剪切模量(G)和柯西压(C12-C44)降低，下降幅度分别为1.24、0.21和0.50 GPa；1.46、0.23和0.40 GPa；0.46、0.08和0.19 GPa以及0.44、0.09和0.10 GPa.体积模量与剪切模量的比值(K/G)几乎不变.其中，Model-Ⅰ下降得最大，Model-Ⅳ次之，Model-Ⅲ下降得最小.这表明，孔洞缺陷导致B炸药的刚性、断裂强度、硬度以及延展性变差，抗变形能力减弱，并且孔洞缺陷越大，B炸药的力学性能恶化得越严重.

表4 不同B炸药模型的弹性系数及力学参数

Table 4 Elasticity coefficients and mechanical parameters of different Composition B models

Model(n)Model-ⅠModel-ⅡModel-ⅢModel-ⅣC11/GPa10.995210.404111.122211.1205C22/GPa11.07858.942010.99469.4817C33/GPa10.74088.698610.758811.0595C44/GPa2.68012.57662.57312.3344C55/GPa2.58092.49242.24862.2494C66/GPa2.92702.29902.56052.5497C12/GPa5.95805.41505.76055.5086C13/GPa6.33206.02575.92365.8366C23/GPa6.27034.24735.80026.1201C15/GPa-0.0195-0.0567-0.1807-0.1465C25/GPa-0.0956-0.08210.1039-0.0996C35/GPa-0.0886-0.0641-0.06260.0586C46/GPa0.01370.24180.0339-0.0269E/GPa6.945.706.736.44μ0.350.350.350.35K/GPa7.776.317.547.37G/GPa2.572.112.492.38K/G3.022.993.033.10(C12-C44)/GPa3.282.843.193.17