徐逸敏 高 贫 王桂香
①南京理工大学化学与化工学院(江苏南京,210094)
②国家民用爆破器材质量监督检验中心(江苏南京,210094)
随着科学技术的不断进步,武器装备也在不断更新换代。含能材料作为武器的主要能量来源也被赋予了更高要求和标准。因此,高能量密度材料(HEDM)和高能量密度化合物(HEDC)成为各个国家的研究重点[1-4]。根据Kamlet和Jacobs总结得到的K-J公式[5],晶体密度对爆速和爆压等性能有重要影响,晶体密度越大,爆轰性能越好。故增大含能化合物的晶体密度是提高爆轰性能的重要手段。
计算含能化合物晶体密度的主要方法中,应用较广的有摩尔体积基团加和法、摩尔折射度基团加和法、基于物质紧密堆积理论的结晶化学法等。然而,这些方法均存在一些不足之处。为此,人们试图寻找一种能既简便、快速又较准确地估算炸药晶体密度的新方法。笔者研究团队曾提出一种基于量子化学计算分子摩尔体积V,进而求化合物理论密度的简易方法:即先采用密度泛函理论对分子进行优化,再用Monte-Carlo统计方法求V,进一步计算得到理论密度ρc,并将其近似作为晶体密度用于计算爆轰性能。该方法被广泛应用于预测各类含能化合物的密度,在寻求HEDC的过程中发挥了较好的作用[6-10]。而对于不同类型的化合物,选择的密度泛函方法和基组不一定相同。文献[6-7,9-10]对氮杂硝胺、多硝基芳烃、硝酸酯类、环状氟化物等有较详细的研究。
金刚烷类化合物具有无(小)张力、高对称刚性结构,有稳定性好、密度大、热值高等优点,是HEDC的理想目标物。但目前已合成的金刚烷类化合物较少。因此,选择已知实验晶体密度ρe的金刚烷类化合物,运用3种密度泛函理论方法和5种基组预测其ρc,通过比较分析,筛选适用于预测金刚烷类化合物晶体密度的方法和基组,为后续研究新型含能金刚烷衍生物的结构与性能提供理论参考。
运用Gaussian 09程序包中的密度泛函理论(B3LYP、M06-2X和B3PW91)方法[11-15]和6-31G、6-31G∗、6-31G∗∗、6-311G∗和6-31+G∗∗基组[16-20],对图1中15种金刚烷衍生物结构进行优化。经过振动分析无虚频,说明得到的是各自势能面上极小值对应的稳定构型。在此基础上,运用Monte-Carlo方法计算化合物的V,并结合分子摩尔质量M求得其理论密度:ρc=M/V。爆速和爆压通过K-J公式[5]计算:
图1 15种金刚烷衍生物的分子结构Fig.1 Molecular structures of 15 adamantine derivatives
式中:D为爆速,km/s;p为爆压,GPa;ρ为炸药密度,g/cm3;N为单位质量炸药爆轰时生成的气体产物的物质的量是气体产物的摩尔质量,g/mol;Q为单位质量炸药的最大爆热,J/g。
表1给出了15种金刚烷衍生物在B3LYP方法和不同基组水平下的ρc和相应的ρe,以便比较。表1中,括号内的数据为偏差D;Daa为平均绝对偏差;Drms是均方根偏差。由表1可见,5种基组计算得到的ρc与ρe符合较好,6-311G∗和6-31+G∗∗求得的ρc更好些。表2为不同方法和基组计算得到的ρe与ρc之间的关系、标准偏差Ds及相关系数R。B3LYP方法5种基组的R和Ds说明:ρc与ρe线性关系均较好;尤其是6-311G∗和6-31+G∗∗,不仅Daa和Drms小,R和Ds也好于其他基组。
表1 B3LYP方法结合不同基组求得的理论密度和实验密度Tab.1 Theoretical density and experimental density obtained by B3LYP method combined with different base groups g/cm3
表2 不同方法求得的理论密度与实验密度间的关系(ρc=a+bρe)Tab.2 Relationship between theoretical density and experimental density obtained by different methods(ρc=a+bρe)
M06-2X和B3PW91结合5种基组计算得到的ρc与ρe之间的偏差均较小。M06-2X方法的Daa分别为0.06、0.07、0.07、0.05、0.05 g/cm3,Drms分别为0.08、0.08、0.08、0.06、0.06 g/cm3;B3PW91方法的Daa对应为0.06、0.07、0.07、0.06、0.05 g/cm3,Drms分别为0.07、0.08、0.08、0.07、0.06 g/cm3,略高于相同基组水平下B3LYP的结果。M06-2X方法6-311G∗和6-31+G∗∗的Daa(0.05 g/cm3)最小;B3PW91方法6-31+G∗∗的Daa(0.05 g/cm3)最小。
M06-2X和B3PW91方法求得的ρc与ρe的R也都大于0.98,说明它们之间线性关系也较好。其中,M06-2X方法6-31G∗、6-31G∗∗、6-311G∗和6-31+G∗∗基组得到的结果(R=0.985 3~0.987 8,Ds=0.041 1~0.043 8)比6-31G基组更好些,也略好于B3PW91(R=0.983 5~0.984 1)和B3LYP方法(R=0.982 7~0.986 9)。
基组越大,计算所需机时和内存以及硬盘空间就越大。尽管3种方法结合5种基组求得的金刚烷类化合物的理论密度均接近于实验密度,且相关性也较好;但如果需要综合考虑计算时间和资源问题,例如后续研究很大的或大量的新型金刚烷含能衍生物时,则建议采用M06-2X/6-31G∗方法,不经校正即可得到较接近于实验密度的结果。这与以往对其他类型含能化合物得到的结论有所不同[6-7,9-10]。
静电势校正法是目前使用较多的一种计算密度的方法,该方法应用经验关系和静电势参数对理论密度做校正,获得密度ρ。但由于关系式中的系数是由B3PW91/6-31G∗∗计算结果拟合得到,因此,也只能在B3PW91/6-31G∗∗水平下使用[27-29]。作为预测密度的一种方法,也为比较和测试本文中的方法,在B3PW91/6-31G∗∗水平下优化15种化合物的结构,计算它们的理论密度和静电势参数;再应用经验公式求得密度,得到Daa为0.06 g/cm3,与M06-2X/6-31G∗的结果(Daa=0.07 g/cm3)非常接近。表明M06-2X/6-31G∗法可应用于金刚烷类化合物的密度预测。
为进一步验证预测密度的可靠性,选择其中6种密度高的金刚烷类化合物,由M06-2X/6-31G∗方法预测的理论密度,通过K-J公式[5]估算了它们的爆轰性能,结果见表3。
表3 由理论密度预测得到的爆轰性能Tab.3 Detonation performance predicted by theoretical density
金刚烷衍生物爆轰性能的实验研究目前还未见文献报道。用M06-2X/6-31G∗方法预测的理论密度求得化合物3#的爆速和爆压分别为8.67 km/s和34.64 GPa,与文献报道的估算结果[30]非常接近:表明采用M06-2X/6-31G∗方法预测得到的晶体密度以及由此估算的爆轰性能是较为可信的。根据HEDC的能量标准(ρ≈1.9 g/cm3,D≈9.0 km/s,p≈40.0 GPa)[31],从理论上判定这些化合物中7#和10#是潜在的HEDC目标物。
采用不同的密度泛函理论方法和基组对15种金刚烷衍生物的理论密度进行研究发现:B3LYP、B3PW91和M06-2X方法结合6-31G、6-31G∗、6-31G∗∗、6-311G∗和6-31+G∗∗基组求得的理论密度与实验密度符合良好,平均绝对偏差均较小,为0.04~0.07 g/cm3;且理论密度与实验密度之间存在良好的线性关系,线性相关系数均高于0.98;标准偏差均小于0.05:表明这些方法均能较准确地预测金刚烷类化合物的晶体密度。
建议使用M06-2X/6-31G∗方法,无需进行校正即可快速、准确地预测新型金刚烷类化合物的晶体密度,进而用于预估爆轰性能,筛选潜在的HEDC。根据预测结果,15种金刚烷化合物中,7#和10#达到HEDC的能量标准,值得进一步实验合成研究。