文章快速检索     高级检索
  含能材料  2018, Vol. 26 Issue (1): 34-45.  DOI: 10.11943/j.issn.1006-9941.2018.01.004
0

引用本文  

何飘, 杨俊清, 李彤, 张建国. 含能材料量子化学计算方法综述[J]. 含能材料, 2018, 26(1): 34-45. DOI: 10.11943/j.issn.1006-9941.2018.01.004.
HE Piao, YANG Jun-qing, LI Tong, ZHANG Jian-guo. Overview on the Quantum Chemical Methods for Energetic Materials[J]. Chinese Journal of Energetic Materials, 2018, 26(1): 34-45. DOI: 10.11943/j.issn.1006-9941.2018.01.004.

基金项目

国家自然科学基金资助(10776002)

作者简介

何飘(1989-),女,博士研究生,含能材料的理论计算研究。e-mail: hepiaodj@163.com

通信联系人

张建国(1974-),男,教授,博士生导师,主要从事含能材料的理论与应用研究。e-mail: zjgbit@bit.edu.cn

文章历史

收稿日期:2017-10-26
修回日期:2017-11-03
含能材料量子化学计算方法综述
何飘, 杨俊清, 李彤, 张建国     
北京理工大学爆炸科学与技术国家重点实验室, 北京 100081
摘要:概述了量子化学基础理论, 详细综述了含能材料关键参数(密度、生成热、爆热、爆速、爆压和撞击感度)的计算方法, 并比较这些方法的特点和适用范围。介绍了CHEETA、EXPLO5等计算软件在含能材料领域的应用。最后, 为满足新一代材料高能稳定与绿色环保的综合要求, 设计了20种新型高氮含能分子, 运用上述量子化学方法估算了其物化和含能性质, 并提出了新型含能化合物的设计原则: (1)高氮分子骨架; (2)零氧平衡; (3)基团调节修饰, 以期促进含能材料的发展。附参考文献76篇。
关键词含能材料     理论计算     量子化学     性能预估     分子设计    
Overview on the Quantum Chemical Methods for Energetic Materials
HE Piao, YANG Jun-qing, LI Tong, ZHANG Jian-guo     
State Key Laboratory of Explosion Science and Technology, Beijing Institute of Technology, Beijing 100081, China
Abstract: We summarized detailedly in this review the basic theory of quantum chemistry, and discussed the methods for calculating key parameters of energetic materials (density, heat of formation, detonation heat, detonation velocity, detonation pressure and impact sensitivity). We also compared the characteristics and applicable scope of these methods. In addition, the application of CHEETA, EXPLO5 and other computing software in the field of energetic materials were briefed. Finally, in order to satisfy the comprehensive requirements for new generation materials, high energy, good stability and environmental friendliness, we designed 20 new high nitrogen molecules, and estimated their physicochemical properties and energetic parameters using the quantum chemistry methods as mentioned above. The design concept and principle for the new energetic compounds was put forward. It was expected that it would promote the development of energetic materials. Seventy-six references were contained.
Key words: energetic materials    theoretical computation    quantum chemistry    properties prediction    molecular design    
1 引言

含能材料主要包括炸药、推进剂和烟火药剂[1]。在民用、国防、航天及空间站等领域已有广泛应用。含能材料的发展, 经历了以下几个重要的阶段。第一个阶段是传统CHON化合物, 如三硝基甲苯(TNT), 黑索今(RDX), 奥克托今(HMX), 四硝酸季戊四醇酯(PETN)等, 分子中同时含有氧化组分和燃料组分, 这类化合物的能量主要源于分子内C链骨架的氧化。第二阶段是笼型化合物(如六硝基六氮杂异伍兹烷(CL-20)、二硝基四氧杂二氮杂四环十二烷(TEX)、八硝基立方烷(ONC)等), 分子的笼状结构有效增加了密度, 显著提高了其性能。目前, 含能材料的研究已进入高能量密度材料(HEDM, high energy density materials)的新阶段。高氮化合物[2]是高能量密度材料中非常重要的一类化合物, 分子中含有大量的NN键和NC键, 分解产物主要为N2并放出大量热, 由于不同类型键的平均键能的差异(N—N键160 kJ·mol-1, N=N键418 kJ·mol-1, 三键954 kJ·mol-1)[3], 使得这些化合物具有高生成热。新一代含能材料必须满足以下三个方面的要求[1]:一是具有高效可靠的能量性质; 二是对破坏性刺激不敏感, 以确保安全使用和能量释放可控; 三是具有较高的热稳定性和水解稳定性、良好的相容性及生态环境友好性。

随着计算机技术的快速发展, 计算化学方法和模拟技术对材料科学领域产生了深刻影响, 通过对材料组分和材料结构的模型化计算, 可实现对材料设计和制备过程的定量描述, 建立结构与性能之间的关系, 并最终按指定目标设计新材料。一般来说, 通过电子及原子层次的计算实现对材料本质的认识, 是材料设计的基础[4]。计算研究方法主要包括量子化学(QC)和分子动力学模拟(MD)。其中, 量子化学的研究对象主要为平衡态、由单分子或几个分子组成的体系、以及周期性体系(晶胞或超晶胞)。近年来, 量子化学在含能材料研究中发挥着关键作用, 主要有三个方面的优势: (1)理论预测物化参数和能量性质, 减少了大规模实验测试(爆炸性能和感度测试)带来的高成本及潜在危险; (2)基于分子水平组合, 引入理想的分子骨架和官能团, 设计丰富多样的新含能化合物, 可高效选择有潜力的候选物用以实验合成; (3)通过研究含能材料合成或热分解过程中各基元反应历程, 建立反应物、过渡态(中间体)和产物之间的联系, 从而一定程度上解释反应机理。

结构-性质关系研究, 作为微(介)观构型与宏观性能的纽带, 对含能材料的发展具有重要意义。大量研究[5-8]表明, 含能化合物的关键属性(如密度、能量、熔点、热稳定性、感度、毒性等)都与其结构有着密切联系。根据含能材料的结构预测其物化性质的方法有基团贡献法和定量结构-性质/活性关系法(QSPR/QSAR)。基团贡献法[9], 是基于化合物分子中基团的数量和类型以及与相近分子的联系, 从而研究基团加和对其性质产生的影响。QSPR/QSAR法[10], 则是将物理化学性质关联到一系列的分子描述符, 通过数学和统计学的方法建立经验和理论模型。随着计算机软件和编码技术的提高, 使得基于化学结构预测含能化合物的性能成为可能。例如, 利用分子结构式、密度和生成热等参数, 借助计算软件或程序代码, 估算爆炸参数(爆压、爆速和比冲)以评价新型含能材料的综合性能。

本文主要综述量子化学计算方法在含能材料研究中的应用, 总结了含能化合物重要性能参数(密度、生成热、爆热、爆速、爆压和感度)的估算方法, 简要介绍了CHEETA、EXPLO5等计算软件的应用。最后, 本文设计了20种新型含能化合物, 并估算其物化与含能性质, 提出了含能化合物的设计原则。

2 量子化学方法

量子化学[11]是理论化学的一个分支, 集中于量子力学在物理模型和化学体系实验中的应用。研究体系为单个原子和分子的基态、激发态和在化学反应过程中产生的过渡态。量子化学方法主要有四类: (1)半经验方法, 包括AM1、MINDO/3和PM3, 利用实验参数求解薛定鄂方程。(2)从头算, 采用几个物理常数(光速、电子和核的质量、普朗克常数等)求解薛定鄂方程。其中最经典的是Hartree-Fock方法。(3)密度泛函理论方法(Density Functional Theory, DFT), 采用泛函对薛定鄂方程进行求解, 密度泛函包涵电子相关, 常用的有PW91、PBE、B3LYP、BLYP等方法。(4)Post-Hartree-Fock方法, 是在HF或自洽场(SCF)方法上改进的一组方法。相比HF, 它们增加了电子相关, 包括电子之间的相互排斥作用。常见的方法有: MP2、MP4、CI、G2、G4等。目前, 密度泛函理论(DFT)[12]方法得到广泛认可和使用。不同于波函数考虑电子自旋的复杂情况, DFT的优势在于随着电子数的增加电子密度维持一定的变量数目, 从而不依赖于体系大小。综合来看, 与Post-Hartree-Fock方法相比, DFT一方面显著降低了计算需求, 从而可以处理较大的多原子分子和大分子体系, 另一方面具有与MP2和CCSD(T)相当的准确性, 使它成为目前计算化学中最流行的方法之一。

3 密度的估算方法

密度是材料的基本物理性质之一。对于含能材料而言, 密度也是影响其能量属性的关键参数。根据Kamlet和Jacobs提出的半经验公式(K-J方程)[13], 爆压依赖于密度的平方, 爆速与密度成正比, 由此看出密度直接关联其爆轰性能, 是开发新材料需关注的要素之一。由于实验装药密度测试的复杂性, 通常用晶体密度或者理论密度代替。目前, 量子化学方法已经用于预测含能化合物的理论密度:基于第一性原理, 对体系的几何优化、分子轨道、热力学和静电势等进行研究, 再结合结构-性质关系(QSPR)参数, 推导出理论估算值。下面将从气态分子、凝聚相化合物、离子化合物三个类别进行介绍。

首先是气态分子, 基于化合物的优化构型, 用Monte-Carlo方法计算分子0.001e/bohr3等电子密度面[14]所包围的体积。经Monte-Carlo方法得到的分子体积具有随机性, 通常要取100次以上计算的平均值。理论密度(ρ)通过分子质量(M)与分子体积(V)之比求得, 计算公式如下:

$ \rho =\left( \frac{M}{{V\left( {0.001} \right)}}\right) $ (1)

基于此方程, 肖鹤鸣等人[15]对硝胺类化合物的密度进行了研究, Rice等[16]考察了180种CHNO化合物, 主要包括硝基芳烃、硝酸盐、硝酸盐酯和硝基化合物等, 发现在密度泛函B3LYP/6-31G(d, p)水平下, 计算的理论密度与实测值吻合较好。

由于含能化合物多为凝聚相, 因此需要进一步修正密度公式。Politzer等人[17]引入分子表面静电势概念来解释晶格中分子间的相互作用, 提出改进后的固体密度计算公式如下:

$ \rho = \alpha \left( \frac{M}{{V\left( {0.001} \right)}}\right) + \beta (\upsilon \sigma _{{\rm{Tot}}}^2) + \gamma $ (2)

式中, 第一部分为0.001 au电子密度等值面包围的分子体积求得的气态分子密度, 第二部分为静电相互作用, 第三部分为拟合校正。υ表示分子表面正负静电势平衡常数, σTot2表示总表面静电势的方差; α, βγ为常见含能化合物实验密度拟合参数。Politzer等人在密度泛函B3PW91/6-31G(d, p)水平下计算了36种常见含能化合物的理论密度, 与B3LYP方法相比, B3PW91方法得到的分子体积略小, 因而估算的理论密度高出约1%。

对于离子化合物MpXq(M表示阳离子, X表示阴离子), 关键在于计算体积。总体积(V)为离子晶体中阳离子体积(VM+)与阴离子体积(VX-)的总和, 表示如下:

$ V = pV_M^ + + qV_X^ - $ (3-a)

如果离子盐含有强氢键作用, 则需要校正氢键对体积的影响, 其计算公式为如下(其中N表示离子盐中氢键数目):

$ {V_{{\rm{correted}}}} = {V_{{\rm{uncorreted}}}} - [0.6763 + 0.9418 \times N] $ (3-b)

Rice等[16]采用B3LYP/6-31G**理论计算方法, 预测了71种离子化合物的密度, 结果表明, 校正前后的均方根偏差分别为5%和1.3%, 这种方法也试用于含有质子或甲基化四唑阳离子的含能离子化合物[18]。但是, 孤立的气相体系没有考虑阴阳离子间的静电相互作用, 使得体积比实际离子晶体的体积要大, 而估算的密度比实际值要小。基于离子表面静电势, Politzer等[19]提出了计算离子化合物密度的新模型:

$ \rho = \alpha \left( {\frac{M}{{{V_{\rm{m}}}}}} \right) + \beta \left( {\frac{{V_{\mathit{s}}^ + }}{{A_{\rm{s}}^ + }}} \right) + \beta \left( {\frac{{\overline V _{\mathit{s}}^ - }}{{A_{\rm{s}}^ - }}} \right) + \delta $ (3-c)

式中, As+表示阳离子表面包含正的静电势部分, Vs+则为该部分静电势的平均值; As-Vs-则代表阴离子相应的描述。α, β, γδ为实验拟合系数。Politzer等人预测了25种离子化合物的密度, 研究表明, 在B3PW91/6-31G(d, p)理论水平下, 计算结果的平均绝对误差为0.033 g·cm-3, 均方根误差为0.040 g·cm-3

4 生成热估算方法

爆轰性能与爆热密切相关, 而爆热可由爆轰反应过程中各物质的生成热(HOF)直接求得, 因此, 生成热是评价含能材料能量性质的重要参数。基于量子化学理论, 可以得到体系的总能量, 通过选用适当的方法, 从而准确计算含能材料分子的生成热。目前有三种方法用于计算分子的气相生成热, 分别是原子化能、原子当量、等键反应。

原子化能方法(Atomization Energies)[20], 利用已知的孤立原子的生成热, 通过计算原子化反应能来预测分子的气相生成热。

计算公式如下:

$ {\Delta _{\rm{f}}}{H^0}\left( {{A_x}{B_y}, 0\;{\rm{K}}} \right) = x{\Delta _{\rm{f}}}{H^0}\left( {A, 0\;{\rm{K}}} \right) + {\rm{ }}y{\Delta _f}{H^0}\left( {\mathit{B}, 0\;{\rm{K}}} \right) - \sum {D_0} $ (4)

式中, ΔfH0 (A, 0 K)和ΔfH0 (B, 0 K)分别表示原子A和B在0 K下的生成热, D0表示分子能量与原子求和总能量的差值。在高精度的G2[21]计算水平下, 这种方法可以准确计算各种有机和无机分子的生成热。但G2需要计算昂贵的电子相关处理, 因此对于大分子体系不太适用。DFT(B3LYP[22-23])预测结果的精度虽然不如G2, 但在合理的误差范围内, 而且计算相对便宜。BAC-MP4[24]则更准确, 这种方法通过实验的键加和, 矫正理论水平上的误差, Melius运用此方法研究90个分子, 计算结果的平均偏差为1.3 kcal·mol-1(5.43 kJ·mol-1)。

等键反应方法(Isodesmic Reaction)[25], 基于Hess定律[26], 通过设计等键反应来预测分子的气相生成热。标准反应焓变表示如下:

$ \Delta H_{{\rm{Reaction}}}^0 = {\Delta _{\rm{f}}}H^0 ({{\rm{Products}}})- {\Delta _{\rm{f}}}H^0({{\rm{Reactants}}}) $ (5)

式中, ΔfH0(Products)和ΔfH0(Reactants)分别表示生成物和反应物的生成热。由量子化学计算可以得到标准反应焓变, 已知其他组分的生成热数据, 就能求得某一组分的标准生成热。这种方法的关键在于设计合理的等键反应, 即反应前后的电子对和化学键的数目守恒, 从而减少在求解量子力学方程中近似处理电子相关时的误差。在B3LYP/6-31G*计算水平下, 预测得到RDX的气相生成热为52.8 kcal·mol-1(220.70 kJ·mol-1), 仅比实验值大[27]

原子当量方法(Atom Equivalent)[28], 与原子化能相似, 都是利用量子化学方法计算由分子形成原子时的能量变化, 并且结合已知原子的相关信息, 从而预测分子的生成热。原子当量方法计算公式如下:

$ \Delta {H_i} = {E_i} - {\sum _j}{n_j}{\varepsilon _j} $ (6)

式中, Ei是分子i的能量, εj记作一个“原子当量”, 定义为εj = (Ei -xj), 其中Ej是分子中原子j的能量, nj是分子i中原子j的数目, xj是对原子j的理论水平校正。Rice等[29]采用B3LYP/6-31G*理论水平计算研究了35个分子的生成热, 通过最小二乘法拟合方程确定了7种原子当量, 包含4种单键原子(C、H、N、O)和3种多重键原子(C′、N′、O′)。结合实验的生成热数值和量子化学方法计算的能量, 校正温度和零点能误差。这种方法一旦原子当量确定, 在处理电子相关或者输入实验数据时, 就不再要求高水平的计算[30]。Mole等[28]在B3LYP /6-31G*理论水平下, 利用原子当量法预测了23种碳氢化合物的生成热, 计算结果的均方根误差为1.72 kcal·mol-1(7.19 kJ·mol-1), 预测值与实验值的最大偏差为6.2 kcal·mol-1(25.92 kJ·mol-1)。为进一步提高计算精度, 在B3LYP /6-311G*理论水平进行研究, 计算值与实验的均方根误差仅有1 kcal·mol-1(4.18 kJ·mol-1), 最大偏差为2.34 kcal·mol-1(9.78 kJ·mol-1)。Habibollahzadeh等人[31]运用类似的方法预测了54种有机化合物的气相生成热, 计算水平以Becke[32]交换项和Perdew[33]相关项作为泛函, 以6-31G(d, p)作为基组。计算值与实验值的平均偏差为3 kcal·mol-1(12.54 kJ·mol-1)。

尽管量子力学和理论模型能较好地预测气相生成热, 但含能材料的标准态通常都是凝聚相(固态或者液态), 因此预测凝聚相的生成热尤为重要。根据Hess定律, 凝聚相的生成热可以由气相生成热和相变热(升华或者汽化)[25]来最终确定。计算公式如下:

$ \Delta H\left( {{\rm{Solid}}} \right) = \Delta H\left( {{\rm{Gas}}} \right) - \Delta H\left( {{\rm{Sublimation}}} \right) $ (7-a)
$ \Delta H\left( {{\rm{Liquid}}} \right) = \Delta H\left( {{\rm{Gas}}} \right) - \Delta H\left( {{\rm{Vaporization}}} \right) $ (7-b)

许多研究发现凝聚相性质(比如升华焓或汽化焓)与分子的静电势之间有密切关联。相变焓计算公式如下[34-35]:

$ \Delta H\left( {{\rm{Sublimation}}} \right) = a{\left( {SA} \right)^2} + b\sqrt {\sigma _{{\rm{Tot}}}^2\upsilon } + c $ (8-a)
$ \Delta H\left( {{\rm{Vaporization}}} \right) = a\sqrt {\left( {SA} \right)} + b\sqrt {\sigma _{{\rm{Tot}}}^2\upsilon } + c $ (8-b)

式中, SA为分子表面积, υ表示分子表面正负静电势平衡常数, σTot2表示总表面静电势的方差。Politzer等[36]采用上述方法预测了34种有机化合物的升华焓, 计算结果与实验值的标准偏差为2.5 kcal·mol-1(10.45 kJ·mol-1), 并且预测了41种化合物的汽化焓, 预测值与实验值的平均误差为0.6 kcal·mol-1(2.51 kJ·mol-1)。根据Hess定律, 由计算得到的升华焓或气化焓, 可以得到含能分子化合物的标准生成热。

对于离子化合物, 首先基于量子化学计算得到气相生成热, 再结合相变热(晶格能)由Born-Haber能量循环最终得到凝聚相生成热(如示意图Scheme 1)[37]

Scheme1 Born-Haber energy cycle

离子化合物MpXq的凝聚相生成热计算公式如下:

$ \begin{array}{l} \Delta H_{\rm{f}}^0\left( {{\rm{ionic}}\;{\rm{salt}}, {\rm{ }}298\ {\rm{ K}}} \right) = {\rm{ }}\Delta H_{\rm{f}}^0\left( {{\rm{cation}}, {\rm{ }}298\ {\rm{ K}}} \right) + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\Delta H_{\rm{f}}^0\left( {{\rm{anion}}, {\rm{ }}298\ {\rm{ K}}} \right) - \Delta {H_{\rm{L}}} \end{array} $ (9)

式中, ΔHL为晶格能, 计算公式如下[38]:

$ \Delta {H_{\rm{L}}} = {U_{{\rm{POT}}}} + [p(n{\rm{M}}/2-2) + q(n{\rm{X}}/2-2)]RT $ (10)

式中, nM和nX表示阳离子Mp+和阴离子Xq-的特性参数, 单原子离子为3, 线性多原子离子为5, 非线性多原子离子为6, UPOT表示晶格势能, 计算公式如下[39]:

$ {U_{{\rm{POT}}}} = \gamma ({\rho _{\rm{m}}}/{{\rm{M_m}}})1/3 + \delta $ (11)

式中, ρm为密度, Mm为分子量, 参数γδ可以查阅文献[39]得到。Gutowski等[40]结合等键反应和MP2完全基组计算水平, 研究了取代唑类(氨基、叠氮基、硝基、含甲基)离子化合物的生成热, 计算结果与实验数据的误差小于3 kcal·mol-1(12.54 kJ·mol-1)。Shreeve等[37]基于Born-Haber能量循环, 系统地研究了119种含能离子化合物的凝聚相生成热, 这种方法简单方便而且能够很好地吻合实验值, 因此可以用来筛选大量的新设计合成的含能离子盐。除了等键反应, Byrd和Rice[41]采用G3MPP2B3或结合原子当量法预测了25种含能离子盐的生成热, 研究结果表明, 与原子当量法相比, G3MP2B3的结果更接近实验值。

5 爆炸参数估算方法 5.1 爆热估算方法

准确测定爆炸分解产物是评价含能材料爆轰性能的首要问题。气体产物的平衡组成可以通过实验测量、热化学平衡推导或建立适当爆轰反应方程得到。在Chapman-Jouguet状态下, 爆炸气体产物的组成依赖于下面两个重要的平衡方程:

$ {\rm{2CO}} \leftrightarrow {\rm{C}}{{\rm{O}}_2} + {\rm{C}}, \Delta {H_0} = - 41.2\ {\rm{ kcal}}\left( { - 172.45\ {\rm{ kJ}}} \right) $ (12)
$ {{\rm{H}}_{\rm{2}}}{\rm{ + CO}} \leftrightarrow {\rm{C}}{{\rm{O}}_2} + {\rm{C}}, \Delta {H_0} = - 31.4\;{\rm{kcal}}\left( { - 131.43\;{\rm{kJ}}} \right) $ (13)

根据LeChatelier′s原理, 高压(高密度)会使得平衡向右移, 而高温使平衡向左移; 固体碳的数量并没有实质性地影响平衡; 假定N2、H2O、CO2、CO和H2为主要的气体产物, 平衡方程各自极端的情况导致了不同的爆炸产物组分。下面主要介绍常用的两种爆炸反应模型。

Kamlet和Jacobs[13]建立的模型中, 假定平衡都向右, 主要爆炸产物是N2、H2O和CO2, 且H2O优先CO2形成。对于给定组成的含能化合物CaHbNcOd, 根据含氧量分为以下三种情况:

(1) db/2, 含氧量不足以使H完全转化为H2O, 则有:

$ {{\rm{C}}_a}{{\rm{H}}_b}{{\rm{N}}_c}{{\rm{O}}_d} \to \frac{1}{2}c{{\rm{N}}_2} + d{{\rm{H}}_2}{\rm{O}} + \left( {\frac{{b - 2d}}{2}} \right){{\rm{H}}_2} + a{\rm{C}}(s) $ (14-a)

(2) 2a+b/2>db/2, 含氧量可以将H转化为H2O, 但不足以将C完全转化为CO2, 则有:

$ \begin{array}{l} {{\rm{C}}_a}{{\rm{H}}_b}{{\rm{N}}_c}{{\rm{O}}_d} \to \frac{1}{2}c{{\rm{N}}_2} + \frac{1}{2}b{{\rm{H}}_2}{\rm{O}} + \left( {\frac{{d - 2b}}{2}} \right){\rm{C}}{{\rm{O}}_2} + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\left[{a-\frac{{\left( {d-2b} \right)}}{2}} \right]{\rm{C}}(s) \end{array} $ (14-b)

(3) d≥2a+b/2, 含氧量可以使全部的H转化为H2O, 并且使全部的C转化为CO2, 剩余的O转化为O2, 则有:

$ {{\rm{C}}_a}{{\rm{H}}_b}{{\rm{N}}_c}{{\rm{O}}_d} \to \frac{1}{2}c{{\rm{N}}_2} + \frac{1}{2}b{{\rm{H}}_2}{\rm{O}} + a{\rm{C}}{{\rm{O}}_2} + \left[{\frac{{\left( {2d-b-4a} \right)}}{4}} \right]{{\rm{O}}_2} $ (14-c)

含能化合物的爆热, 被定义为爆炸反应焓变的负值, 即反应物与爆炸分解产物的生成热差值, 表示为:

$ {Q_{{\rm{det}}}} \cong - \Delta {H_0} = - \frac{{[\Delta {H_{\rm{f}}}({\rm{detonationproducts}})-\Delta {H_{\rm{f}}}({\rm{explosive}})]}}{{Mr}} $ (15)

因为N2(g)、O2(g)和C(s)的标准生成热为零, 反应中只需考虑H2O(g)和CO2(g)以及含能化合物的标准生成热, 就能估算出爆热。此外, 这里的爆热指的是反应焓变, 当考虑反应能变时, 二者计算的爆热相差1~15 cal·g-1(4.18~62.78 J·g-1)。

对34种CHNO炸药的热化学进行计算, 发现94%的气体产物为CO、H2O、H2、N2和CO2[42]。Keshavarz[43]建立的模型中, 对于含能化合物CHNOFCl, 假定所有N都形成N2, 所有F都形成HF, 所有Cl都形成HCl; 氧原子优先转化为H2O, 而碳原子则优先转化为CO而不是CO2。对给定组成的含能化合物CaHbNcOdFeClf, 分为以下四种情况:

(1) 当da时, 则有:

$ \begin{array}{l} {{\rm{C}}_a}{{\rm{H}}_b}{{\rm{N}}_c}{{\rm{O}}_d}{{\rm{F}}_e}{\rm{C}}{{\rm{l}}_f} \to e{\rm{HF}} + f{\rm{HCl}} + \frac{c}{2}{{\rm{N}}_2} + d{\rm{CO}} + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;(a - d){\rm{C}}\left( s \right) + \left( {\frac{{b - e - f}}{2}} \right){{\rm{H}}_2} \end{array} $ (16-a)

(2) 当da且(b-e-f)/2>d-a时, 则有:

$ \begin{array}{l} {{\rm{C}}_a}{{\rm{H}}_b}{{\rm{N}}_c}{{\rm{O}}_d}{{\rm{F}}_e}{\rm{C}}{{\rm{l}}_f} \to e{\rm{HF}} + f{\rm{HCl}} + \frac{c}{2}{{\rm{N}}_2} + a{\rm{CO}} + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;(\mathit{d} - \mathit{a}){{\rm{H}}_{\rm{2}}}{\rm{O}} + \left( {\frac{{b - e - f}}{2} - d + a} \right){{\rm{H}}_2} \end{array} $ (16-b)

(3) 当da+(b-e-f)/2且d≤2a+(b-e-f)/2时, 则有:

$ \begin{array}{l} {{\rm{C}}_a}{{\rm{H}}_b}{{\rm{N}}_c}{{\rm{O}}_d}{{\rm{F}}_e}{\rm{C}}{{\rm{l}}_f} \to e{\rm{HF}} + f{\rm{HCl}} + \frac{c}{2}{{\rm{N}}_2} + \left( {\frac{{b - e - f}}{2}} \right){{\rm{H}}_{\rm{2}}}{\rm{O}} + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\left( {2a - d + \frac{{b - e - f}}{2}} \right){\rm{CO}} + \left( {d - a - \frac{{b - e - f}}{2}} \right){\rm{C}}{{\rm{O}}_2} \end{array} $ (16-c)

(4) 当d>2a+(b-e-f)/2时, 则有:

$ \begin{array}{l} {{\rm{C}}_a}{{\rm{H}}_b}{{\rm{N}}_c}{{\rm{O}}_d}{{\rm{F}}_e}{\rm{C}}{{\rm{l}}_f} \to e{\rm{HF}} + f{\rm{HCl}} + \frac{c}{2}{{\rm{N}}_2} + \left( {\frac{{b - e - f}}{2}} \right){{\rm{H}}_{\rm{2}}}{\rm{O}} + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\mathit{a}{\rm{C}}{{\rm{O}}_{\rm{2}}} + \left[ {\frac{{\left( {2d - b + e + f} \right)}}{4} - a} \right]{{\rm{O}}_2} \end{array} $ (16-d)

对比Kamlet模型和Keshavarz模型这两种模型, 如果考虑气态水H2O(g)为分解产物, Kamlet计算结果与实验值的均方根误差为1.006 kJ·g-1, 而Keshavarz模型为0.954 kJ·g-1; 如果考虑液态水H2O(l)为分解产物, Kamlet模型和Keshavarz模型的均方根误差分别为1.049 kJ·g-1和1.364 kJ·g-1, 由此可见, 模型的选取对计算结果会造成一定的影响。近年来Kamlet模型受到含能材料研究者更多的青睐, 尤其是对于理想的CHNO化合物, 其预测结果相对可靠。

5.2 爆速爆压估算方法

爆速爆压是含能材料的重要能量属性。为了描述产物状态和与组成有关的热力学变量, 选择合适的状态方程非常重要。爆炸反应体系同时含有气体和固体产物, 因此需要合适的气体状态方程以及合理的混合规则。下面简要介绍含能材料研究的三个状态方程。Kistiakowsky-Wilson[44](K-W)方程, 是硬球分子的模型和固态模型的折衷。采用K-W方程和体积元参数, 研究CHNO组成的含能体系, 预测得到的爆轰温度相当低; 在高装载密度条件下, 主要的碳氧产物是CO2而不是CO, 估算得到的爆轰压力和爆轰速度都接近于实验值。K-W方程如下:

$ p = RT(1 + X{{\rm{e}}^X})/V + f(V)\;\;\;\;\;\;\;\;X = b/V $ (17-a)

Fickett和Cowan[45]通过降低对函数f (V)的依赖, 增加可调参数β, 使b成为温度的函数, 得到改进K-W方程如下:

$ p{V_{\rm{g}}}/RT = 1 + X{{\rm{e}}^{\beta X}}\;\;\;\;\;X = k\sum {x_i}{k_i}/{V_{{g}}} $ (17-b)

对于大多数CHNO炸药, 用HM代替K-W方程式中的∑xiki, H是常数13.76, M是平均气体分子量, 由M/Vg=ρg(ρg为气体产品的密度), 得到变形的K-W方程如下:

$ p = (RT{\rho _{\rm{g}}}/M)(1 + X{{\rm{e}}^{\beta X}})\;\;\;\;X = kH{\rho _{\rm{g}}}/{(T + \theta )^\alpha } $ (17-c)

引入两个重要的参量NG, N表示每克炸药产生的气体爆轰产物的摩尔数, G表示炸药转变为气体的重量分数, 即有定义, NM=G; 再令kH=A, βkH=B, 则有K-W方程的一般形式:

$ p = \frac{{NRT{\rho _{{g}}}}}{G}\left\{ {1 + \frac{{\left( {A{\rho _{{g}}}} \right)}}{{{{\left( {T + \theta } \right)}^\alpha }}}{\rm{exp}}\left[ {\frac{{\left( {B{\rho _g}} \right)}}{{{{\left( {T + \theta } \right)}^\alpha }}}} \right]} \right\} $ (17-d)

该方程中, ki项被消除, 因而更容易评估影响p的各种因素, 由此可以看出, 爆炸性质的计算不仅取决于爆炸产物的组成, 还同时依赖于平衡状态方程。

Kamlet和Jacobs[13]提出采用相对简单的K-J经验方程来预测CHNO体系的爆压和爆速, 这些关键的爆炸参数只依赖于每单位重量的炸药产生的爆炸性气体的摩尔数、气体的平均分子量、爆炸反应的化学能以及装药密度。K-J方程表示如下:

$ p = K{\rho _{02}}\mathit{\Phi }, \;\;K = 1.558, \;\;\;\mathit{\Phi } = N{M^{1/2}}{Q^{1/2}} $ (18)
$ D = A{\mathit{\Phi }^{1/2}}(1 + B{\rho _0}), \;\;A = 1.01, \;\;\;B = 1.30 $ (19)

式中, p表示爆压, GPa; D表示爆速, km·s-1; N表示每克炸药产生气体的摩尔数; Q表示爆热, cal·g-1; ρ0表示炸药的装药密度, g·cm-3, 实际计算过程中, ρ0往往用晶体密度或理论密度代替。K-J方程[46-48]适用于初始密度大于1.0 g·cm-3的CHNO体系, NMQ值的估计是以H2O-CO2分解产物假定为基础的, 因此计算只需输入炸药的元素组成、生成热和载荷密度信息, 即可得到相应的爆压p和爆速D。K-J方程计算得到的p, 与基于K-W方程的RUBY的计算机代码所预测的数值相当。K-J方程计算得到的D, 与实验值的平均绝对误差为1%。

另一种估算爆轰压力的经验模型是C-J方程, 实验表明C-J爆轰压力与载荷密度的平方基本成比例[49]。固态密度作为炸药最重要的特性之一, 决定了炸药的性能。C-J经验方程[43]表示如下:

$ {p_{{\rm{CJ}}}} = 1.588\alpha {(M{Q_{{\rm{det}}}})^{1/2}}{\rho _{02}} - 11.17 $ (20)

式中, pCJ为爆轰压力, GPa; α为每克炸药产生的爆炸气体摩尔数;M为气体产物的平均分子量;ρ0为装载密度。C-J方程假定CHNOFCl炸药的分解产物主要为HF、HCl、CO、N2、H2O、H2和CO2, 根据氧平衡确定爆炸反应方程, 需要首先计算出炸药重要的参量α, MQdet, 即可得到给定组成炸药的爆压。在炸药密度大于0.8 g·cm-3的情况下, 计算值和实验值比较一致。C-J方程还适用于CHNOFCl的混合体系[50-51], 计算值也能较好地符合实验结果。

6 撞击感度估算方法

评价高能材料的感度性能, 对含能材料的使用、运输和储藏等过程具有重要意义。撞击感度是其中最受关注的感度性能之一, 一般通过初始分解反应或声音探测, 这两种测量方法都与炸药的分子结构参数密切相关[52]。由于实验测试撞击感度的复杂性, 理论预测表现出较大的优势。预测撞击感度的方法有很多, 本文主要介绍:量子化学方法、半经验模型方法、QSPR神经网络和其他相关性质方法。

利用量子化学计算, Poltzer等人分析孤立分子等电子密度表面的静电势(ESP), 基于统计全局描述分子表面EPS上电荷分布, 从而建立了这些统计量与材料宏观性质的关系。这种关系的函数表达被称为一般交互特性函数(General Interaction Properties Function, GIPF)[34-35, 53-54], 其形式取决于相应的宏观属性, 表示如下:

$ {\rm{Properties}} = {\rm{f}}({\rm{SA}}, \Pi, \sigma _{{\rm{Tot}}}^2, \nu ) $ (21)

式中, SA、Π、σTot2, 、ν为基于特定电子密度等表面ESP的全局性质参量。其中, SA为等电子密度表面的分子表面积(通常为0.001 e·bohr-3等值面), Π为分子表面静电势V(r)的平均偏导数, σTot2V(r)在分子表面的总偏差, ν为平衡常数。关联撞击感度H50的GIPF方法有5种模型, 分别表示如下:

(1) 模型1(中点正电荷模型):

$ {H_{50}} = {y_0} + a{\rm{exp}}\left( { - b{V_{{\rm{Mid}}}}} \right) + c{V_{{\rm{Mid}}}} $ (22-a)

式中, a=18922.7503 cm, b=0.0879 kcal·mol-1, c=-0.3675 cm/kcal·mol-1, y0=63.6485 cm。利用模型1预测的撞击感度与实验结果的均方根误差为26.1 cm, 拟合的相关系数为0.96。但是这个模型似乎只适用于研究硝基芳烃, 对硝铵类炸药的预测结果偏差较大。

(2) 模型2(GIPF参量Vs+Vs-):

$ {H_{50}} = {a_1} + {a_2}{\rm{exp}}\left[{-\left( {{a_3}\left| {\overline V _s^ +-\left| {\overline V _s^-} \right|} \right|} \right)} \right] $ (22-b)

式中, a1=9.1949 cm, a2=803.4464 cm, a3=0.3663 kcal·mol-1, Vs+Vs-分别表示分子表面正负静电势。方程拟合系数为0.94, 预测结果的均方根误差为31.2 cm。该模型适用于敏感炸药, 但对于FOX-7和NTO钝感炸药的预测结果并不理想。

(3) 模型3(GIPF平衡常数ν):

$ {H_{50}} = {a_1} + {a_2}{\rm{exp}}\left( {{a_3}\nu } \right) $ (22-c)

式中, a1=29.3248 cm, a2=0.001386 cm, ν为平衡常数, 方程拟合系数为0.80, 预测结果的均方根误差为54.9 cm。该模型对非常敏感的炸药(H50<30 cm)以及非常不敏感炸药(H50>100 cm)的预测结果较差。

(4) 模型4(Q参量):

$ {H_{50}} = {a_1} + {a_2}{\rm{exp}}\left[{{\rm{-}}\left( {{a_3}\left( {Q-{a_4}} \right)} \right)} \right] $ (22-d)

式中, a1=27.8331 cm, a2=0.1135 cm, a3=11.0793 kcal·g-1(46.37 kJ·g-1), a4=1.6606 kcal·g-1(6.95 kJ·g-1), 预测结果的均方根误差为24.1 cm, 但是利用该模型预测的撞击感度偏高, 尤其对于3-硝基-1, 2, 4-三唑-5-酮和硝基胍炸药。

(5) 模型5(平衡常数ν与爆热Q):

$ {H_{50}} = {a_1}\mathit{exp}\left[{{a_2}\nu-{a_3}\left( {Q-{a_4}} \right)} \right] $ (22-e)

式中, a1= 1.3410 cm, a2=8.1389 cm, a3=6.7922 kcal·g-1(28.43 kJ·g-1), a4=1.4737 kcal·g-1(6.16 kJ·g-1), 方程拟合系数为0.95, 预测结果的均方根误差为28.1 cm, 该模型对除了CL-14之外的苯并氧化呋咱类炸药的预测结果偏低, 但对于一些酸性的炸药预测结果偏高。

半经验模型方法, Keshavarz等人[55]提出了CHNO炸药的元素组成与撞击感度之间的简单关系, 可适用于不同类别的含能材料。对于CaHbNcOd炸药, 撞击感度表示为:

$ {\rm{log}}\;{\mathit{H}_{50}} = \frac{{{A_1}a + {B_1}b + {C_1}c + {D_1}d}}{M} $ (23)

式中, abcd表示化学计量系数, M表示分子量, A1B1C1D1为实验拟合参数。

这个方程只需要输入元素组成信息即可求得h50, 提供了简单的估算感度性能的方法。为了预测不同类型炸药的撞击感度, 将上述方程分为以下三种:

(1) 多硝基芳香化合物

$ {\rm{log}}\;{\mathit{H}_{50}} = \frac{{11.76a + 61.72b + 26.89c + 11.48d}}{M} $ (24-a)

(2) 含有α-CH和α-N-CH(如tetryl)的多硝基芳香化合物, 以及硝铵

$ {\rm{log}}\;{\mathit{H}_{50}} = \frac{{47.33a + 23.50b + 2.357c - 1.105d}}{M} $ (24-b)

(3) 多硝基脂肪族化合物

$ {\rm{log}}\;{\mathit{H}_{50}} = \frac{{81.40a + 16.11b - 19.08c + 1.089d}}{M} $ (24-c)

这些半经验模型描述了CaHbNcOd炸药分子的组成与撞击感度的相关性, 而且预测的结果与实验值吻合较好, 其最大的优势在于非常简单方便。

近年来, QSPR神经网络模型已经开始广泛应用于材料性质评估。Wang等[56]采用QSPR方法以156种硝基化合物的分子结构作为描述模型, 并结合电子拓扑特征(ETSI方法), 预测了目标分子的撞击感度。然而, 由于混合数据来自不同的结构(含有C—NO2、N—NO2基团的芳香族或脂肪族), 预测的结果并不理想。Xu等[57]在QSPR研究中, 使用3D描述符对156个结构上不同的含能化合物, 进行多元线性回归, 从而建立了非线性模型(由ANN建立), 可以很好地用于评价新的含能材料。但是这些QSPR方法往往无法评估化学物理的起始参量。

此外, 关联撞击感度的其他相关性质方法, 还包括:第一性原理能带[58]、硝基最大电荷[59]、冲击力学性质[60]、放热反应速率常数[61]、晶格自由体积[62]、晶体质量和粒度[63]等等。这些参量都在一定程度上反映了炸药的感度性能。

7 计算程序简介

随着计算机的快速发展, 计算模拟和仿真技术在高能材料中的应用日益广泛, 并且在空间尺度和时间尺度方面的研究均不断深入。运用大规模并行计算解决材料结构和性质问题, 表现出简便节约的优势。在含能材料领域, 理论计算不仅能够有效预测爆轰参数, 还可以分析反应产物的化学平衡组成, 以及物理化学性质(原子组成、生成焓、密度)等。这些应用的基础在于:

(1) 理想爆炸的数学模型; (2)吉布斯描述的特征函数极值原理; (3)真实气体的热力学方程(反应产物)具有广泛的压力和温度。现在许多数值方法和程序都被用于热力学计算, 预测凝聚炸药的爆轰参数, 例如BKW Fortran[64]、ARPEGE[65]、Ruby[66]、TIGER[67]、CHEETAH[68]、EXPLO5[69]、ZMWNI[70]、EMDB[71]等。其中, CHEETAH和EXPLO5是含能材料研究者常用的两种软件。

CHEETAH[68]是美国领先开发的代码, 对高能材料的性能进行化学模拟, 使用范围从核武器到枪支推进器以及简易爆炸装置。CHEETAH能够可靠运行分子计算, 基于热力学性质和密度参数, 并转换这些参量用以计算炸药的性能(爆速和爆压)。目前正在使用的CHEETAH 6.0程序, 涵盖了近半个世纪的炸药实验数据, 以及几百上千个反应物和产物。通过成千上万的微处理器, 仿真结果能够揭示许多转瞬即逝的爆炸反应细节。例如, 在少于100 ps的范围, 被引爆的高能炸药表现出类似于金属的导电性质, 氢离子可以移动, 而其他元素仍然紧密相连。在爆轰过程中, 高能炸药组分不仅具有分子特点, 还表现出类似于高能离子的集合等离子体。CHEETAH还可以用于先进的流体动力学模拟, 研究表明高能材料的安全性依赖于高度和粉末状炸药的颗粒大小。

EXPLO5[69]是一种热化学计算机程序, 可以预测理想的高能炸药以及各种爆炸配方、推进剂和烟火混合物的性能, 因此成为高能材料合成、配方和数值模拟的重要工具。一方面, 可以应用化学平衡爆轰的稳态模型, 估算爆轰参数(爆炸产物组成、爆轰速度和压力、爆炸能量、爆温和爆热等); 另一方面, 可以应用定压燃烧条件下守恒方程, 评价理论火箭性能(特定的冲动, 推力系数、流速等); 以及应用恒容燃烧条件计算特定的能量(或力, 动力)。通过改进的White、Johnson和Dantzig的自由能小型化技术, 计算爆轰和燃烧产物的平衡组成。该程序使用Becker-Kistiakowsky-Wilson(BKW)和Jacobs-Cowperthwaite-Zwisler(JCZ3)分别作为爆炸产物气体状态方程, 和燃烧产物的理想气体维里方程, 以及冷凝产物的Murnaghan方程。EXPLO5程序计算沿爆轰产物冲击绝热曲线状态下的平衡组成和热力学参数、C-J态和C-J状态下爆轰参数, 以及沿等熵膨胀状态的热力学参数。该程序在膨胀等熵线上, 将相对体积压力数据拟合到Jones-Wilkin-Lee(JWL)状态方程, 具有非线性曲线拟合特点, 使得可以评估JWL系数, 以及机械做功的能量。

8 新型含能化合物设计

本课题组主要从事高氮类含能材料的理论与实验研究工作, 长期的研究发现: (1)多氮杂环(三唑、四唑、呋咱等)的分子稳定性好, 且分子中含氮量越高, 释放的能量越高, 这主要是氮氮键较高的键能在爆炸过程中转化的缘故。另外, 多氮杂环分子的CH比例较低, 分解产物主要为氮气, 可以最大限度地降低对环境的危害。(2)氧平衡为零的化合物分子往往具有较高的正生成热和爆炸性能。这主要是因为分子中的O刚好使得H全部转化成H2O, C全部转化成CO2, 使得体系能量最大限度释放。(3)用不同含能基团修饰化合物一方面可以极大丰富材料的种类, 另一方面可以改善结构和性能。根据经验, 设计了20种新型高氮含能化合物[72-74](见图 1), 基于量子化学基础, 运用密度泛函理论, 在B3LYP/6-311++G**计算水平下, 对其几何结构优化、单点能、振动模式、热力学性能、静电势、电子密度等进行了研究。设计的新化合物分子均为势能面上的能量最低点, 振动分析显示没有虚频, 表明所有化合物理论上具有稳定的几何构型。

图 1 新设计的20种高氮含能化合物分子结构 Fig.1 Molecular structures for newly designed 20 high-nitrogen energetic compounds

(1) 分子结构

新设计的20种分子(图 1)主要分为四类:第一类是叠氮关环产物(1~5), 均由三唑和四唑环构成, 各分子中所有的C原子和N原子均处在同一个平面上, 且键长分布均匀, 说明分子共轭性较大。分子有对称性, 具有稳定的几何构型。第二类是FOX类似物(6~7), 以FOX-7结构片段和四唑作为母环, 引入多个硝基以提高其能量性质, 并适当引入氨基以增加其稳定性。第三类是零氧平衡分子(8~10), 以四唑基和三嗪为基础, 构筑六个氮原子直接相连的高氮结构, 引入多个硝基以提高其能量性质。由于电子的排斥作用, 使得硝基以一定角度排列在环外, 且不与环共平面。第四类是N10系列化合物(11~20), 分子中含有10个直接相连的氮原子, 以偶氮四唑为基本骨架结构, 引入不同含能的基团, 分子表现出对称结构。

(2) 能量性质

新设计的20种目标化合物的能量性质, 包括密度、生成热、爆热、爆速、爆压和氧平衡, 见表 1, 为便于比较, 常见含能化合物的能量性质也列于表 1中。

表 1 新设计的20种化合物及常规炸药的能量性质 Tab.1 Energetic properties of newly designed compounds and common explosives

高氮含量使得这些新化合物(除了11)都具有较高的正生成热, 且均高于已有的常见含能化合物的生成热, 其中化合物10的生成热最高, 达1625.31 kJ·mol-1。大多数化合物的密度在1.8 g·cm-3以上, 其中化合物681112的密度接近2.0 g·cm-3, 与高能量密度炸药CL-20相当。近一半化合物的爆热高于1500 cal·g-1, 其中化合物71617具有非常显著的爆热且均高于CL-20, 化合物17的爆热最高, 达到1833 cal·g-1。大多数化合物的爆速在9 km·s-1左右, 其中化合物7121617爆速高于CL-20。大多数化合物的爆压与HMX相当, 其中化合物78111217的爆压达到40 GPa, 化合物17的爆压最高, 为42.14 GPa。综上所述, 设计的新型化合物均具有高的正生成热、高密度、高爆热和良好的爆炸性能, 有望成为潜在的高能量密度含能材料。

9 结论与展望

综述了含能材料计算方法, 概述了量子化学理论基础, 详细介绍了物理化学性质(密度、生成热等)计算方法, 比较分析了爆炸参数(爆热、爆速、爆压)和感度性能(撞击感度)的预测方法。另外设计了20种新型高氮含能化合物分子, 并预测了其性能, 筛选出了性能良好的潜在含能材料, 提出了新型含能化合物分子的设计思路: (1)高氮分子骨架; (2)零氧平衡; (3)基团调节修饰。

在高氮分子骨架上引入不同的含能基团使其尽量达到零氧平衡是设计高能、稳定、绿色环保的新型含能化合物的一种有效手段。未来含能材料的开发和研制中, 量子化学将发挥不可替代的作用, 尤其是对新型含能化合物的探索, 可预先对其分子设计、性能估算和感度评价等开展理论研究, 提高了研制效率, 缩短了研究周期。但必须承认的是, 理论计算不能代替实验, 只有理论与实验相结合, 才能科学完善促进含能材料长远发展。

参考文献
[1]
Gao H, Shreeve J M. Azole-based energetic salts[J]. Chemical Reviews, 2011, 111(11): 7377-7436. DOI:10.1021/cr200039c
[2]
Zhang Q, Shreeve J N M. Growing catenated nitrogen atom chains[J]. Angewandte Chemie International Edition, 2013, 52(34): 8792-8794. DOI:10.1002/anie.201303297
[3]
Talawar M B, Sivabalan R, Mukundan T, et al. Environmentally compatible next generation green energetic materials(GEMs)[J]. Journal of Hazardous Materials, 2009, 161(2): 589-607.
[4]
居学海, 叶财超, 徐司雨. 含能材料的量子化学计算与分子动力学模拟综述[J]. 火炸药学报, 2012, 35(2): 1-9.
JU Xue-hai, YE Cai-chao, XU Si-yu. Overview on quantum chemical computing and molecular dynamic simulations of energetic materials[J]. Chinese Journal of Explosives & Propellants, 2012, 35(2): 1-9.
[5]
王桂香, 肖鹤鸣, 居学海, 等. 含能材料的密度、爆速、爆压和静电感度的理论研究[J]. 化学学报, 2007, 65(6): 517-524.
WANG Gui-xiang, XIAO He-ming, JU Xue-hai, et al. Theoretical studies on densities, detonation velocities and pressures and electric spark sensitivities of energetic materials[J]. Acta Chimica Sinica, 2007, 65(6): 517-524.
[6]
严启龙, 宋振伟, 安亭, 等. 含能材料物理化学性能理论预估研究进展[J]. 火炸药学报, 2016, 39(5): 1-12.
YAN Qi-long, SONG Zhen-wei, AN-ting, et al. Research progress in theoratical prediction of physicochemical properties for energetic materials[J]. Chinese Journal of Explosives & Propellants, 2016, 39(5): 1-12.
[7]
曹霞, 向斌, 张朝阳. 炸药分子和晶体结构与其感度的关系[J]. 含能材料, 2012, 20(5): 643-649.
CAO Xia, XIANG Bin, ZHANG Chao-yang. Review on realationships between the molecular and crystal structure of explosives and their sensitivities[J]. Chinese Journal of Energetic Materials(Hanneng Cailiao), 2012, 20(5): 643-649.
[8]
葛忠学, 王伯周, 牛永洁, 等. 含能材料分子设计与性能预估研究[J]. 计算机与应用化学, 2007, 24(12): 1714-1716.
GE Zhong-xue, WANG Bo-zhou, NIU Yong-jie, et al. Molecular design and property prediction of energetic materials[J]. Computers and Applied Chemistry, 2007, 24(12): 1714-1716. DOI:10.3969/j.issn.1001-4160.2007.12.027
[9]
Gmehling J. Potential of group contribution methods for the prediction of phase equilibria and excess properties of complex mixtures[J]. Pure and Applied Chemistry, 2003, 75(7): 875-888.
[10]
Dearden J C, Cronin M T D, Kaiser K L E. How not to develop a quantitative structure-activity or structure-property relationship (QSAR/QSPR)[J]. SAR and QSAR in Environmental Research, 2009, 20(3-4): 241-266. DOI:10.1080/10629360902949567
[11]
Atkins P W, Friedman R. Molecular Quantum Mechanics (4th ed.)[M]. New York: Oxford University Press, 2005.
[12]
Parr R G, Yang W. Density-functional theory of atoms and molecules[M]. New York: Oxford University Press, 1989.
[13]
Kamlet M J, Jacobs S J. Chemistry of detonation. Ⅰ.A Simple method for calculating detonation properties of C—H—N—O explosives[J]. Journal of Chemical Physics, 1968, 48: 23-35. DOI:10.1063/1.1667908
[14]
Bader R F W, Carroll M T, Cheeseman J R, et al. Properties of atoms in molecules: atomic volumes[J]. Journal of the American Chemical Society, 1987, 109(26): 7968-7979. DOI:10.1021/ja00260a006
[15]
Qiu L, Xiao H, Gong X, et al. Crystal density predictions for nitramines based on quantum chemistry[J]. Journal of Hazardous Materials, 2007, 141(1): 280-288. DOI:10.1016/j.jhazmat.2006.06.135
[16]
Rice B M, Hare J J, Byrd E F C. Accurate predictions of crystal densities using quantum mechanical molecular volumes[J]. The Journal of Physical Chemistry A, 2007, 111(42): 10874-10879. DOI:10.1021/jp073117j
[17]
Politzer P, Martinez J, Murray J S, et al. An electrostatic interaction correction for improved crystal density prediction[J]. Molecular Physics, 2009, 107(19): 2095-2101. DOI:10.1080/00268970903156306
[18]
Zhang X, Zhu W, Wei T, et al. Densities, heats of formation, energetic properties, and thermodynamics of formation of energetic nitrogen-rich salts containing substituted protonated and methylated tetrazole cations: a computational study[J]. The Journal of Physical Chemistry C, 2010, 114(30): 13142-13152. DOI:10.1021/jp1050782
[19]
Politzer P, Martinez J, Murray J S, et al. An electrostatic correction for improved crystal density predictions of energetic ionic compounds[J]. Molecular Physics, 2010, 108(10): 1391-1396. DOI:10.1080/00268971003702221
[20]
Curtiss L A, Raghavachari K, Redfern P C, et al. Assessment of Gaussian-2 and density functional theories for the computation of enthalpies of formation[J]. The Journal of Chemical Physics, 1997, 106(3): 1063-1079. DOI:10.1063/1.473182
[21]
Curtiss L A, Raghavachari K, Trucks G W, et al. Gaussian-2 theory for molecular energies of first-and second-row compounds[J]. The Journal of Chemical Physics, 1991, 94(11): 7221-7230. DOI:10.1063/1.460205
[22]
Becke A D. Density-functional thermochemistry. Ⅲ. The role of exact exchange[J]. The Journal of Chemical Physics, 1993, 98(7): 5648-5652. DOI:10.1063/1.464913
[23]
Lee C, Yang W, Parr R G. Development of the Colle-Salvetti correlation-energy formula into a functional of the electron density[J]. Physical Review B, 1988, 37(2): 785-789. DOI:10.1103/PhysRevB.37.785
[24]
Melius C F. Chemistry and Physics of Energetic Materials[C]//Vol. 309 (S. N. Bulusu, Ed. ), Kluwer Academic Publishers, Dordrecht, 1990.
[25]
Atkins P W. Physical Chemistry[M]. Oxford: Oxford University Press, 1982.
[26]
Hehre W J, Radom L, Schleyer P V R, et al. Ab Initio Molecular Orbital Theory[M]. New York: John Wiley, 1986: 271-298.
[27]
Hariharan P C, Pople J A. The influence of polarization functions on molecular orbital hydrogenation energies[J]. Theoretica Chimica Acta, 1973, 28(3): 213-222. DOI:10.1007/BF00533485
[28]
Mole S J, Zhou X, Liu R. Density functional theory(DFT) study of enthalpy of formation. 1. consistency of DFT energies and atom equivalents for converting DFT energies into enthalpies of formation[J]. The Journal of Physical Chemistry, 1996, 100(35): 14665-14671. DOI:10.1021/jp960801h
[29]
Rice B M, Pai S V, Hare J. Predicting heats of formation of energetic materials using quantum mechanical calculations[J]. Combustion & Flame, 1999, 118(3): 445-458.
[30]
Dewar M J S, Storch D M. Development and use of quantum molecular models. 75. Comparative tests of theoretical procedures for studying chemical reactions[J]. Journal of the American Chemical Society, 1985, 107(13): 3898-3902. DOI:10.1021/ja00299a023
[31]
Habibollahzadeh D, Grice M E, Concha M C, et al. Nonlocal density functional calculation of gas phase heats of formation[J]. Journal of Computational Chemistry, 1995, 16(5): 654-658. DOI:10.1002/(ISSN)1096-987X
[32]
Becke A D. Density-functional exchange-energy approximation with correct asymptotic behavior[J]. Physical Review A, 1988, 38(6): 3098-3100. DOI:10.1103/PhysRevA.38.3098
[33]
Perdew J P. Density-functional approximation for the correlation energy of the inhomogeneous electron gas[J]. Physical Review B, 1986, 33(12): 8822-8824. DOI:10.1103/PhysRevB.33.8822
[34]
Politzer P, Murray J S, Brinck T, et al. Immunoanalysis of Agrochemicals[M]. Washington, DC: ACS Symp. Ser. 586, American Chemical Society, 1994.
[35]
Murray J S, Politzer P. Quantitative Treatment of Solute/Solvent Interactions, Theoretical and Computational Chemistry[M]. Amsterdam: Elsevier Scientific, 1994.
[36]
Politzer P, Murray J S, Edward Grice M, et al. Calculation of heats of sublimation and solid phase heats of formation[J]. Molecular Physics, 1997, 91(5): 923-928. DOI:10.1080/002689797171030
[37]
Gao H, Ye C, Piekarski C M, et al. Computational characterization of energetic salts[J]. The Journal of Physical Chemistry C, 2007, 111(28): 10718-10731. DOI:10.1021/jp070702b
[38]
Jenkins H D B, Roobottom H K, Passmore J, et al. Relationships among ionic lattice energies, molecular (formula unit) volumes, and thermochemical radii[J]. Inorganic Chemistry, 1999, 38(16): 3609-3620. DOI:10.1021/ic9812961
[39]
Jenkins H D, Tudela D, Glasser L. Lattice potential energy estimation for complex ionic salts from density measurements[J]. Inorganic chemistry, 2002, 41(9): 2364 DOI:10.1021/ic011216k
[40]
Gutowski K E, Rogers R D, Dixon D A. Accurate thermochemical properties for energetic materials applications. Ⅱ. Heats of formation of imidazolium-, 1, 2, 4-triazolium-, and tetrazolium-based energetic salts from isodesmic and lattice energy calculations[J]. The Journal of Physical Chemistry B, 2007, 111(18): 4788-4800. DOI:10.1021/jp066420d
[41]
Byrd E F C, Rice B M. A comparison of methods to predict solid phase heats of formation of molecular energetic salts[J]. The Journal of Physical Chemistry A, 2009, 113(1): 345-352. DOI:10.1021/jp807822e
[42]
Rice B M, Hare J. Predicting heats of detonation using quantum mechanical calculations[J]. Thermochimica Acta, 2002, 384(1): 377-391.
[43]
Keshavarz M H, Pouretedal H R. An empirical method for predicting detonation pressure of CHNOFCl explosives[J]. Thermochimica Acta, 2004, 414(2): 203-208. DOI:10.1016/j.tca.2003.11.019
[44]
Kistiakowsky G B, Wilson Jr E B. The hydrodynamic theory of detonation and shock waves[R]. OSRD Rept, 1941: 114.
[45]
Cowan R D, Fickett W. Calculation of the detonation properties of solid explosives with the Kistiakowsky-Wilson Equation of state[J]. The Journal of Chemical Physics, 1956, 24(5): 932-939. DOI:10.1063/1.1742718
[46]
Kamlet M J, Ablard J E. Chemistry of detonation. Ⅱ. Buffered Equilibria[J]. Journal of Chemical Physics, 1968(1): 36-42.
[47]
Kamlet M J, Dickinson C. Chemistry of detonations. Ⅲ. Evaluation of the simplified calculational method for Chapman-Jouguet detonation pressures on the basis of available experimental information[J]. Journal of Chemical Physics, 1968, 48(1): 43-50. DOI:10.1063/1.1667939
[48]
Kamlet M J, Hurwitz H. Chemistry of detonations. Ⅳ. Evaluation of a simple predictional method for detonation velocities of C—H—N—O explosives[J]. Journal of Chemical Physics, 1968, 48(8): 3685-3692. DOI:10.1063/1.1669671
[49]
Chirat R, Pittion-Rossillon G. A new equation of state for detonation products[J]. The Journal of Chemical Physics, 1981, 74(8): 4634-4642. DOI:10.1063/1.441653
[50]
Mader C L. Numerical modeling of explosives and propellants[J]. Numerical Modeling of Explosives & Propellants, 1998, 99(5): 1384-1392.
[51]
Hobbs M L, Baer M R. Calibrating the BKW-EOS with a large product species data base and measured C-J properties[J]. Chemical Explosives, 1992
[52]
Zeman S, Jungová M. Sensitivity and performance of energetic materials[J]. Propellants Explosives Pyrotechnics, 2016, 41(3): 426-451. DOI:10.1002/prep.201500351
[53]
Murray J S, Brinck T, Lane P, et al. Statistically-based interaction indices derived from molecular surface electrostatic potentials: a general interaction properties function(GIPF)[J]. Journal of Molecular Structure Theochem, 1994, 307(Supplement C): 55-64.
[54]
And B M R, Hare J J. A quantum mechanical investigation of the relation between impact sensitivity and the charge distribution in energetic molecules[J]. Journal of Physical Chemistry A, 2002, 106(9): 1770-1783. DOI:10.1021/jp012602q
[55]
Keshavarz M H, Pouretedal H R. Simple empirical method for prediction of impact sensitivity of selected class of explosives[J]. Journal of Hazardous Materials, 2005, 124(1-3): 27-33. DOI:10.1016/j.jhazmat.2005.05.009
[56]
Wang R, Jiang J, Pan Y, et al. Prediction of impact sensitivity of nitro energetic compounds by neural network based on electrotopological-state indices[J]. Journal of Hazardous Materials, 2009, 166(1): 155 DOI:10.1016/j.jhazmat.2008.11.005
[57]
Xu J, Zhu L, Fang D, et al. QSPR studies of impact sensitivity of nitro energetic compounds using three-dimensional descriptors[J]. Journal of Molecular Graphics & Modelling, 2012, 36(6): 10
[58]
Zhu W, Xiao H. First-principles band gap criterion for impact sensitivity of energetic crystals: a review[J]. Structural Chemistry, 2010, 21(3): 657-665. DOI:10.1007/s11224-010-9596-8
[59]
Zhang C. Review of the establishment of nitro group charge method and its applications[J]. Journal of Hazardous Materials, 2009, 161(1): 21-28. DOI:10.1016/j.jhazmat.2008.04.001
[60]
Mathieu D. Toward a physically based quantitative modeling of impact sensitivities[J]. Journal of Physical Chemistry A, 2013, 117(10): 2253-2259. DOI:10.1021/jp311677s
[61]
Mathieu D, Alaime T. Predicting impact sensitivities of nitro compounds on the basis of a semi-empirical rate constant[J]. Journal of Physical Chemistry A, 2014, 118(41): 9720-9726. DOI:10.1021/jp507057r
[62]
Svatopluk, Zeman, Marcela, et al. Impact sensitivity in respect of the crystal lattice free volume and the characteristics of plasticity of some nitramine explosives[J]. Chinese Journal of Energetic Materials, 2015, 23(12): 1186-1191.
[63]
Jungová M, Zeman S, Yan Q L. Recent advances in the study of the initiation of nitramines by impact using their N-15 NMR chemical shifts[J]. Central European Journal of Energetic Materials, 2014, 11(3): 285-294.
[64]
Mader C J. FORTRAN BKW: a code computing the detonation properties of explosives[R]. Los Alamos Science Laboratory, Report LA-3704, 1967.
[65]
Cheret R. The numerical study of the detonation products of an explosive substance[R]. French Commission of Atomic Energy, Report CEA-R-4122, 1971.
[66]
Levin H B, Sharples R E. Operator′s manual for RUBY[R]. Lawrence Livermore Laboratory, Report UCRL-6815, 1962.
[67]
Cowperthwaite M, Zwisler W H. Tiger computer program documentation[M]. Stanford Research Institute, Publication No. Z106, 1973.
[68]
Fried L E, Glaesemann K R, Howard W M, et al. CHEETAH 5. 0 Users Manual[M]. Manuscript UCRL-MA-117541 Rev 3. Lawrence Livermore National Laboratory. 2007.
[69]
Sućeska M. EXPLO5 V6. 01[M]. Zagreb: Croatia; Brodarski Institute, 2012.
[70]
Grys S, Trzciński W A. Calculation of combustion, explosion and detonation characteristics of energetic materials[J]. Central European Journal of Energetic Materials, 2010, 7(2): 97-113.
[71]
Keshavarz M H, Klapötke T M, Sućeska M. Energetic materials designing bench(EMDB), Version 1.0[J]. Propellants, Explosives, Pyrotechnics, 2017, 42(8): 854-856. DOI:10.1002/prep.v42.8
[72]
He P, Zhang J G, Wang K, et al. Extensive theoretical studies on two new members of the FOX-7 family: 5-(dinitromethylene)-1, 4-dinitramino-tetrazole and 1, 1′-dinitro-4, 4′-diamino-5, 5′-bitetrazole as energetic compounds[J]. Physical Chemistry Chemical Physics, 2015, 17(8): 5840-5848. DOI:10.1039/C4CP04883K
[73]
He P, Zhang J G, Wang K, et al. Combination multinitrogen with good oxygen balance: molecule and synthesis design of polynitro-substituted tetrazolotriazine-based energetic compounds[J]. Journal of Organic Chemistry, 2015, 80(11): 5643 DOI:10.1021/acs.joc.5b00545
[74]
He P, Zhang J G, Wu L, et al. Computational design and screening of promising energetic materials: novel azobis(tetrazoles) with ten catenated nitrogen atoms chain[J]. Journal of Physical Organic Chemistry, 2016
[75]
Politzer P, Murray J S. Some perspectives on estimating detonation properties of C, H, N, O compounds[J]. Central European Journal of Energetic Materials, 2011, 8(3): 209-220.
[76]
Song X, Li J, Hou H, et al. Extensive theoretical studies of a new energetic material: Tetrazino-tetrazine-tetraoxide(TTTO)[J]. Journal of Computational Chemistry, 2009, 30(12): 1816-1820. DOI:10.1002/jcc.v30:12
图文摘要

The basic theory of quantum chemistry was summarized and the calculation methods of key parameters of energetic materials were discussed in detail. The characteristics and applicable scope of these methods were compared. In addition, the application of CHEETA, EXPLO5 and other computing software in the field of energetic materials were briefly introduced. Finally, 20 kinds of new high nitrogen molecules were designed and their physicochemical properties and energetic parameters were estimated by using above-mentioned quantum chemistry methods.