今天偶见我硕士导师课题组发表在国内著名的物理化学学报上的一篇文章《适用于TATB, RDX, HMX含能材料的全原子力场的建立与验证》。

含能材料，就是俗称的炸药。就是能够在外界刺激下瞬时释放大量能量、热量和气体产物的化合物。这篇文章涉及到了三种含能材料：TATB，是著名的钝感炸药，用来添加到那些高热敏炸药中降低感度，但同时又不显著降低爆轰能量密度；RDX，就是知名的“黑索金”， 不仅能做心脏病药物、毒鼠药，而且做炸药时更是威力巨大，是洛克比空难中恐怖分子所使用的塞姆汀塑胶炸药的主要成分之一；HMX, 又名octogen，中文“奥可托金”， 比黑索今更加安全、高能，在维基中被称为“现今综合性能最好的炸药”，不仅曾用以起爆核武器，也用于固体火箭推进剂当中。推动嫦娥上天的功劳很可能也有HMX一份。不过HMX的合成成本及其昂贵，在分子结构上它是一个含有四个胺硝基基团的八元环。他们论文没有讨论的另外一个类似于HMX但是呈笼装的CL20，其用作火箭推进剂产生的能量比HMX多20%， 但是比HMX更加难以合成。RDX, HMX, CL20等被称为“高能量密度材料”（high energy density materials, HEDM），简称“高能材料”，也是我硕士论文研究的主题。

虽然我硕士期间做过一些探索，但是都没有达到发表要求而躺在了硕士大论文中。尽管硕士毕业已近7年，但我仍然对高能材料抱有特别的兴趣。因为含能材料代表了一类有趣的化合物，这类化合物常温常态下处在一种亚稳定存在的状态，一旦受到外界的刺激，比如加热、摩擦、机械冲击等，就会迅速解体生成大量气体，并释放大量热量、声光。千百年来，人们小心翼翼地搭配、存储、合成此类化合物，以满足工业开矿、航天推进、节日烟火以及军火工业的需求。我就一直非常的好奇，如何在化学上将原子组合成含有巨大能量的分子并存储为“亚稳态”？它们能瞬时释放大量能量的秘密是什么？常温常压保持一定稳定性的关键又在哪？如何设计性能更好的含能化合物？

因为好奇，虽然不再涉足含能材料的研究项目，但一直都在关注这方面的报道。由于这类研究多属于国防和企业的核心秘密，公开报道的主要是基础理论研究，比如含热、机械稳定性能的测量、冲击响应的分子机理、起爆机制等问题。美国有一份期刊，《含能材料》(Energetic Materials)，刊载了很多美国和欧洲的研究进展，中国的绵阳九院和南京大学也时有文章发表。不过由于军火工业的争议性，我所了解到的此类最好的论文都发表在物理化学类期刊如JPC， JCP， Mol.Phy.等。SCIENCE有时候也关注一下关于含能材料冲击响应的理论模拟进展。

高能材料不仅成本高，而且爆轰实验也很危险，爆炸产物和剩余物对环境也产生危害，因此，分子模拟研究就成为一种非常经济可行的研究方式，特别是这些材料在高温高压下的稳定性能、爆轰机理两个方面。前者可以使用包括经典分子动力学、第一原理分子动力学等方法模拟，后者则可以使用量子化学进行气态分子计算，或者第一原理或者半经验方法凝聚态计算模拟。在这些模拟研究中，加温加压不过是改变一下参数，再也不用冒生命危险了。

在师弟这篇论文中，他们以自己从第一原理计算拟合得到的分子力学力场为基础，使用分子动力学方法，模拟了三种材料的晶体结构，在不同压力下固体的体积变化曲线(状态方程)，并计算了衡量固体材料体积-压力响应性质的弹性模量。这是一个关于分子模拟的极好的研究训练，同时也勾起了我以往对这些问题的思考，在这里浅显的讨论一下，希望能为自己或者他人的进一步研究提供一些参考。

首先是关于晶体堆积力，Crystal packing force，它很重要，因为这种力是常温常压下炸药保持一定热稳定和机械稳定的关键，也是决定固体对温度、压力的体积响应性质的关键。有一个起爆理论称为“热点起爆”机理，认为爆轰起始于晶体颗粒之间的缝隙、晶体内的缺陷空隙处的分子的分解。显然晶体堆积力是分子在晶格附近运动的关键因素。

晶体堆积力实际上就是固体中的分子间作用力，分为排斥和吸引两部分。吸引主要由分子的偶极-偶极作用引起，在分子力学中用原子点电荷之间的静电作用表示； 瞬时偶极引起的部分称为色散作用，它和范德华排斥作用在分子力学中一起构成非静电作用部分。这部分的力场参数，是通过拟合晶体的密度和升华焓得到的。这是拟合非键作用参数的通常做法。不过，由于待拟合参数个数远多于可获得的性质，这个拟合结果多少有一些不确定性。我一直在思考的一个问题是，我们到底需要多大精确度的力场，以描述晶体堆积力？力场模型本身是一种粗糙的、经验的模型，它是否足以达到我们需要的精度？如果可以达到，那么我们需要怎样做才能改善现有拟合过程的一些缺点，达到最好的结果？

已经发表的关于二聚体、三聚体分子间作用的量子化学计算很多。不过，多数是讨论其最低能的作用构型和作用力大小。但是，我认为我们需要一个胺硝基体系之间作用的一些定性认识，即为了表现胺硝基体系之间的作用，我们需要考虑哪些物理作用模型？胺硝基基团的局部偶极大致是多少？在不同的化合环境下如何变化？极化对分子间作用势能贡献占多大比例？为了准确描述晶体的弹性模量和剪切模量，极化作用是否是不可或缺的？

虽然目前我还不清楚这些问题的答案，但是，作为分子间作用的一个关键，是原子点电荷，因为晶体分子间吸引力的主要部分是其分子偶极-偶极作用，而一个准确的偶极距离不开准确的原子点电荷。师弟这篇论文中比较了气相分子DFT计算的Mulliken和ESP电荷，并因为Mulliken电荷在不同化合物中更为一致所以选作拟合参考，但若能比较一下分子偶极距的大小，应该会更有说服力。

其次，胺硝基中存在以N-N键为轴心的二面角，是此类化合物的一个结构特征。我曾经发现这个二面角的旋转动力学对TNAZ的熔点有显著影响。虽然此二面角的旋转能垒很高，但其势能面底部平缓。因此，推断此二面角的旋转对晶体压力-体积响应可能有一定贡献。若要进一步改善力场键参数，不妨考虑这一点。

第三，关于原子类型的定义和力场参数迁移性，是一直都令我困惑的问题。原子类型划分的越多越细，意味着拟合结果会更好，但是由于未知数个数增加，参数拟合结果的稳定性变差，而且参数迁移性变差；原子类型划分的越粗越少呢，则每个化合物的描述变差。师弟的力场，对TATB/RDX的密度拟合结果都不超过2%，但对HMX的两种晶型的密度拟合结果偏大约4~5%，当然，在这完全在可接受的误差范围内（实验平均误差10%）。这里的问题是，对于原子类型的划分，我们是否能够确立一个标准？即应该怎样衡量原子之间的差别从而区分为不同的原子类型，即衡量什么性质？其次是差别多大需要区分为不同原子类型？

最后，分子力场方法面对的主要竞争还在于固体第一原理方法以及半经验量子化学方法的逐步普及。在这个问题上，分子力场方法有其优势，比如对固体宏观密度、升华焓等性质的描述。这一方面，二者的对比应该也很吸引人。

Ref: 金钊等, 适用于TATB, RDX, HMX含能材料的全原子力场的建立与验证，物理化学学报
doi: 10.3866/PKU.WHXB201402113