表1(Table 1)
2. 嘉兴学院生物与化学工程学院, 浙江 嘉兴 314001;
3. 中国航天科技集团公司四院四十二所, 湖北 襄阳 441003;
4. 中国工程物理研究院流体物理研究所, 四川 绵阳 621900
2. College of Biological, Chemical Sciences and Engineering, Jiaxing University, Jiaxing 314001, China;
3. The 42nd Institute of the Fourth Academy of China Aerospace Science and Technology Corporation, Xiangyang 441003, China;
4. Institute of Fliud Physics, China Academy of Engineering Physics, Mianyang 621900, China
实际使用的高能复合材料如固体推进剂, 主要通过多年实践经验积累, 依靠反复配方实验而获得, 这一直是能源材料领域研究和发展的主流[1-3]。但配方实验如果缺乏理论指导, 则将耗费大量人力、物力和财力, 研发周期也较长, 还存在不安全隐患, 尤其是无法预测未知配方及其性能, 故而难以适应和满足航天、国防和国民经济对新型固体推进剂日益增长的迫切需求。
近年来迅猛发展的分子动力学(MD), 可对庞大复杂体系结构和性能的模拟计算, 能为高能复合材料配方设计提供理论指导。通过MD模拟, 揭示研究高能材料的组成、结构与性能关系的研究工作虽已不少。但最初主要集中于单体炸药[4-8], 后来发展到双组分体系[9-13], 较少四组分体系的MD研究[14-16]也未涉及爆炸性能判别和热膨胀性能。因实际使用的固体推进剂均为许多组分构成的复杂体系, 而至今未见五组分以上高能体系的MD研究文献。提供多组分高能复合体系组成、结构及其与安全性、热膨胀性和力学性能, 具有重要的理论和实用价值。
本研究对4种不同配比(PEG/NG/BTTN)/AP/HMX五组分体系和5个不同温度下(PEG/NG/BTTN)/AP/HMX/Al六组分体系进行MD模拟研究。其中PEG为聚乙二醇, NG为硝化甘油, BTTN为1, 2, 4-丁三醇三硝酸酯, AP为高氯酸铵, HMX为环四亚甲基四硝胺。用文献[12-13, 17-19]建议和采用的易爆燃组分引发键的最大键长(Lmax)判据, 对不同配比五组分和不同温度六组分体系的感度相对大小进行判别; 还预测了六组分体系的热膨胀系数, 并以静态力学分析求得其弹性力学性能。
2 模型构建和模拟细节 2.1 (PEG/NG/BTTN)/AP/HMX五组分模型搭建与MD模拟按照所设计的组分及其含量和混合先后次序, 先将NG和BTTN以1:1质量比混合搭好模型, 再与PEG混合组成粘结剂。然后再依据(PEG/NG/BTTN)、AP和HMX质量比依次为2.5:4.1:1.4、2.5:3.9:1.9、2.5:3.5:2.3和2.5:2.9:2.9的配比要求, 分别构建4种五组分体系的初始模型, 它们依次含2572、2839、3222和3872个原子。
将上述搭建好的初始模型, 均匀地分别放入200Å×200Å×200Å的周期箱中。在自行建立的(含AP)MPCFF力场[13]和NVT系综下进行20 ps MD模拟, 以确保体系达到热力学平衡; 然后缓慢缩小周期箱体积, 再进行MD模拟, 直到达到新的热力学平衡; 重复此过程直到体系密度接近其理论值。取其最终构象, 温度设定为295 K, 在NPT系综下, 运行平衡后，再进行了1 ns的MD模拟研究。初始原子运动速度按Maxwell-Boltzman分布确定。在周期性边界条件和时间平均等效于系综平均等基本假设之上求解牛顿运动方程。积分采用Verlet方法, 时间步长为1 fs。模拟过程中温度采用Andersen[20]方法, 范德华(vdW)和静电作用(Coulomb)分别用atom-based[21]和Ewald[22]加和方法, 截断半径取9.5×10-10 m, 并进行截断尾部校正。利用平衡后的原子运动轨迹进行统计分析。
图 1给出MPCFF力场下经NVT-MD模拟所得4种配比的(PEG/NG/BTTN)/AP/HMX五组分混合体系的平衡结构。由该图可见, 5种分子互相混合在一起, 它们之间存在较强的相互作用。图中a、b、c和d分别对应的是粘结剂(PEG/NG/BTTN)与AP和HMX的质量比。
|
图 1 不同配比(PEG/NG/BTTN)/AP/HMX混合体系的平衡结构 Fig.1 Equilibrium structures of (PEG/NG/BTTN)/AP/HMX mixture with different proportions |
对于(PEG/NG/BTTN)/AP/HMX/Al六组分混合体系, 利用上面五组分(PEG/NG/BTTN)、AP、HMX质量比为2.5:2.9:2.9的平衡结构(图 1d), 添加金属铝Al的晶体, 其结构来源于MS程序自带的晶体数据库, 于是得到(PEG/NG/BTTN)、AP、HMX和Al的质量比为2.5:2.9:2.9:1.7的六组分初始模型, 其中共含4179个原子。模拟方法与细节与上述五组分体系相同, 分别进行5个温度(195，245，295，345，395 K)的NVT-MD模拟研究。作为示例, 图 2给出295 K下该六组分体系的MD模拟平衡结构。
|
图 2 295 K下(PEG/NG/BTTN)/AP/HMX/Al六组分体系的平衡结构 Fig.2 Equilibrium structure of (PEG/NG/BTTN)/AP/HMX/Al system at 295 K |
研究复合体系的安全性, 应抓住和关注其中易爆燃组分结构的变化[23-24]。在(PEG/NG/BTTN)/AP/HMX五组分混合体系中, 易爆燃组分为NG、BTTN和HMX。一般认为硝酸酯炸药(NG、BTTN)的引发键为O—NO2键, 硝胺炸药(HMX)的引发键为N—NO2键, 它们最易引发分解和起爆。表 1给出五组分混合体系中NG、BTTN和HMX引发键(分别为O—N和N—N键)的平均键长(Lave)和最大键长(Lmax)。为利于比较, 表 1中最后一行给出295 K该三者晶体的MD模拟所得Lave和Lmax。
|
表 1 295 K不同配比(PEG/NG/BTTN)/AP/HMX五组分体系的Lave和Lmax Tab.1 Lave and Lmax of (PEG/NG/BTTN)/AP/HMX five-component system with different mass ratio at 295 K |
由表 1可见, 与NG、BTTN和HMX晶体中MD模拟所得平均键长(Lave)和最大键长(Lmax)相比, 由于混合体系中多组分之间存在彼此相互作用, 故在(PEG/NG/BTTN) /AP/HMX五组分体系中, 各引发键的平均和最大键长均有所增大。但在多组分体系中Lave几近相等, 与配方组成几无关联; 而各Lmax随配比不同变化较为敏感。根据易爆燃组分引发键最大键长(Lmax)判据[12-13, 17-19]可推知, 配比为2.5:3.5:2.3五组分体系的感度将最大。有趣的是, 此时NG、BTTN和HMX的Lmax依次为1.450, 1.461, 1.540Å, 在4种配比下恰均为最大; 这就表明, 在2.5:3.5:2.3配比下, 这些引发键相对较弱, 易于断裂引发分解和起爆, 亦即导致该体系感度增大。由此可见, 这一配比不太安全, 应予回避。
MD模拟能提供各种体系在广义平衡结构下的键长统计分布, 包括在不同温度下的键长分布, 这是MD的优势所在。表 2给出5个不同温度下(PEG/NG/BTTN)/AP/HMX/Al六组分混合物中易爆燃物NG、BTTN和HMX的引发键Lave和Lmax。
|
|
表 2 不同温度下(PEG/NG/BTTN)/AP/HMX/Al六组分体系中易爆燃物引发键的平均键长(Lave)和最大键长(Lmax) Tab.2 Average trigger bond length (Lave)and maximum trigger bond length (Lmax) of the deflagrating components in (PEG/NG/BTTN)/AP/HMX/Al six-component system at different temperatures |
从表 2可见, 对于该六组分体系, 当温度从195 K升至395 K时, HMX、NG和BTTN的引发键平均键长Lave变化很小, 几近相等, 亦即200 K的温度差不足以改变其Lave值。另一方面, 三者各自的最大键长Lmax均随温度升高而单调地增大, 与感度随温度升高而增大的实验事实相符, 这就表明, 对于复杂多组分体系热和撞击感度的相对大小, 确实也可用引发键最大键长(Lmax)加以判别。
3.2 不同温度下六组分体系热膨胀系数物体受热膨胀现象在自然界普遍存在。温度变化引起的热变形严重影响加工或测量精度, 直接影响到生产和使用, 在高能材料领域颇受关注。度量物体热变形的重要参量是组成该物体材料的热膨胀系数。较为精确的热膨胀系数是将实用热膨胀系数以微分方法表示成温度的函数。由此作成表格或用图形表示, 以便在实际使用时能精确知道某一温度点的热膨胀系数[25]。
先前基于MD模拟结果只对单体炸药热膨胀参数进行过研究报导[26-27], 这里借助以上六组分体系超晶胞不同温度下MD模拟结果, 由公式(1)和(2)求得微分线膨胀系数(α)和微分体膨胀系数(β):
| $ \alpha = \frac{1}{{{L_n}}}\frac{{{\rm{d}}L}}{{{\rm{d}}T}} $ | (1) |
| $ \beta = \frac{1}{{{V_n}}}\frac{{{\rm{d}}V}}{{{\rm{d}}T}} $ | (2) |
式中, Ln、Vn分别为所测温度下超晶胞的平均晶胞边长和体积。在不追求高精度的情况下, 晶胞边长和体积与温度T的关系均可表示为:
| $ Y = A + BT $ | (3) |
对不同温度下晶胞参数和体积进行一元线性拟合即可求得A和B值。相应的B值即dL/dT或dV/dT。
图 3a和图 3b分别给出该六组分体系的体积和密度随温度的递变曲线。从图 3a和图 3b可见, 该六组分体系的体积随温度升高而增大, 密度则随温度升高而下降。395 K时的晶胞体积相对于195 K膨胀了约3.3%, 密度则相应地下降约3.1%;这主要归因于混合体系中的原子间平均距离随温度升高而增大, 温度越高, 相邻原子间平均距离越大, 以至晶胞参数增大, 晶体体积膨胀, 密度则下降。
|
图 3 六组分混合体系的体积-温度和密度-温度关系 Fig.3 Volume vs. temperature and density vs. temperature for six-component mixture |
根据热膨胀系数定义, 从晶胞体积与温度之间的关系, 拟合得到体膨胀系数和平均线膨胀系数, 如表 3所示。从表 3可见, 六组分体系的体膨胀系数和平均线膨胀系数均随温度升高缓慢下降, 表明混合体系的热膨胀能力随温度的升高有所减弱, 这从图 3a和图 3b也可推知。物质的热膨胀行为是其原子非简谐振动的直接结果。热膨胀大小与材料的化学组成、晶体结构和化学键强度等多种因素密切相关。热膨胀测量可求得Grüneisen和Rayleigh值[28-29], 借以揭示热膨胀系数与其它性能之间的关系。
|
|
表 3 不同温度下(PEG/NG/BTTN)/AP/HMX/Al六组分体系的热膨胀系数 Tab.3 Thermal expansion coefficients of (PEG/NG/BTTN)/AP/HMX/Al six-component system at different temperatures |
力学性能直接关系到固体推进剂的加工成型和实际使用, 是固体推进剂的最重要性能之一。表 4给出经静态力学分析[30-31]求得的(PEG/NG/BTTN)/AP/HMX/Al六组分体系295 K下的力学性能。
|
|
表 4 (PEG/NG/BTTN)/AP/HMX/Al六组分体系的力学性能 Tab.4 Mechanical properties of (PEG/NG/BTTN)/AP/HMX/Al six component system |
由表 4可见, 除对角元Cij和C12、C13、C23共9个矩阵元外, 未列入的弹性系数值等于或接近于零, 可见所模拟的(PEG/NG/BTTN)/AP/HMX/Al六组分体系近似为正交各向异性弹性体。弹性系数Cij可反映非均匀材料在各处的不同弹性效应。对(PEG/NG/BTTN)/AP/HMX/Al六组分体系而言, 其C11、C22很大, 说明产生同样的形变需承受很大的应力; 而其它的系数较小, 表明(PEG/NG/BTTN)/AP/HMX/Al六组分体系具有明显的各向异性。
表 4中的K/G表示体系的体模量与剪切模量的比值, 可用来关联体系的延展性[32], 即材料发生形变而不产生裂缝的能力, 其值越高则体系的延展性越好; C12—C44为柯西压, 亦可用来关联体系的延展性[33], 若其值为负, 表明材料显脆性; 若为正, 则表明材料延展性较好; 二者的差异可见文献[17]。从表 4还可见, (PEG/NG/BTTN)/AP/HMX/Al六组分体系的柯西压为正值, 表明其晶体具有一定的延展性。
4 结论通过(PEG/NG/BTTN)/AP/HMX五组分体系不同配比和(PEG/NG/BTTN)/AP/HMX/Al六组分体系不同温度的MD模拟, 得出如下主要结论:
(1) 在4种配比的五组分体系中, 当(PEG/NG/BTTN)/AP/HMX＝2.5:3.5:2.3时, 各易爆燃组分引发键的最大键长均为最大, 表明该体系的感度相对最大。
(2) (PEG/NG/BTTN)/AP/HMX/Al六组分体系中易爆燃组分引发键最大键长(Lmax)均随温度升高而单调递增, 与感度随温度升高的实验事实相一致。表明对于复杂多组分体系的热和撞击感度, 亦可用Lmax加以判别。
(3) 利用曲线拟合法, 首次预测了六组分高能体系各温度点的热膨胀系数, 发现热膨胀系数随温度升高缓慢下降, 表明混合体系的热膨胀能力随温度升高有所减弱。
(4) 由静态分析法求得六组分体系的弹性力学性能, 比较弹性系数, 可见其具有明显的各向异性; 从K/G值和柯西压(C12—C44)值较大, 说明其延展性较好。
| [1] |
郑剑, 侯林法, 杨仲雄. 高能固体推进剂技术回顾与展望[J].
固体火箭技术, 2001, 24(1): 28-34. ZHENG Jian, HOU Lin-fa, YANG Zhong-xiong. The progress of high energy propellants[J]. J Solid Rocket Technol, 2001, 24(1): 28-34. |
| [2] |
李上文, 赵凤起, 袁潮, 等. 国外固体推进剂研究与开发的趋势[J].
固体火箭技术, 2002, 25(2): 36-42. LI Shang-wen, ZHAO Feng-qi, YUAN Chao, et al. Tendency of research and development for overseas solid propellants[J]. J Solid Rocket Technol, 2002, 25(2): 36-42. |
| [3] |
庞爱民, 郑剑. 高能固体推进剂技术未来发展展望[J].
固体火箭技术, 2004, 27(4): 289-293. PANG Ai-min, ZHENG Jian. Prospect of the research and development of high energy solid propellant technology[J]. J Solid Rocket Technol, 2004, 27(4): 289-293. |
| [4] |
Sorescu D C, Rice B M, Thompson D L. Isothermal-isobaric molecular dynamics simulations of 1, 3, 5, 7-Tetranitro-1, 3, 5, 7-tetraazacyclooctane(HMX)crystals[J].
J Phys Chem B, 1998, 1: 6692-6695. |
| [5] |
Bedrov D, Smith G D, Sewell T D. Thermal conductivity of liquid octahydro-1, 3, 5, 7-tetranitro-1, 3, 5, 7-tetrazocine (HMX) from molecular dynamics simulations[J].
Chem Phys Lett, 2000, 324(1): 64-68. |
| [6] |
Bedrov D, Ayyagari C, Smith G D, et al. Molecular dynamics simulations of HMX crystal polymorphs using a flexible molecule force field[J].
Comput Aided Mater Des, 2001, 8: 77-85. DOI:10.1023/A:1020046817543 |
| [7] |
Manaa M R, Fried L E, Melius C F, et al. Decomposition of HMX at extreme conditions: a molecular dynamics simulation[J].
J Chem Phys, 2002, 106: 9024-9029. DOI:10.1021/jp025668+ |
| [8] |
Sewell T D, Menikoff R, Bedrov D, et al. A molecular dynamics simulation study of elastic properties of HMX[J].
J Chem Phys, 2003, 119: 7417-7426. DOI:10.1063/1.1599273 |
| [9] |
肖继军, 马秀芳, 黄玉成, 等. TATB/氟聚物PBX力学性能的分子动力学模拟[J].
含能材料, 2004(增刊): 488-492. XIAO Ji-jun, MA Xiu-fang, HUANG Yu-cheng, et al. Molecular dynamics simulation of mechanical properties of TATB/fluorine-polymers PBX[J]. Chinese Journal of Energetic Materials(Hanneng Cailiao), 2004(Suppl.): 488-492. |
| [10] |
Xiao J J, Fang G Y, Ji G F, et al. Simulation investigations in the binding energy and mechanical properties of HMX-based polymer-bonded explosives[J].
Chin Sci Bull, 2005, 50(1): 21-26. DOI:10.1360/982004-147 |
| [11] |
Xu X J, Xiao H M, Huang H, et al. Molecular dynamics simulations on the structures and properties of ε-CL-20-based PBXs—Primary theoretical studies on HEDM formulation design[J].
Sci China Ser B, 2007, 50: 737-745. |
| [12] |
朱伟, 肖继军, 郑剑, 等. 高能混合物的感度理论判别——不同配比和不同温度AP/HMX的MD研究[J].
化学学报, 2008, 66(23): 2592-2596. ZHU Wei, XIAO Ji-jun, ZHENG Jian, et al. A theoretical criterion for sensitivity of composites——molecular dynamics studies on AP/HMX systems at various concentrations and temperatures[J]. Acta Chim Sinica, 2008, 66(23): 2592-2596. DOI:10.3321/j.issn:0567-7351.2008.23.004 |
| [13] |
Zhu W, Wang X J, Xiao J J, et al. Molecular dynamics simulations of AP/HMX composite with a modified force field[J].
J Hazard Mater, 2009, 167(1): 810-816. |
| [14] |
马秀芳, 赵峰, 肖继军, 等. HMX基多组分PBX结构和性能的模拟研究[J].
爆炸与冲击, 2007, 27(2): 109-115. MA Xiu-fang, ZHAO Feng, XIAO Ji-jun, et al. Simulation study on structure and property of HMX-based multi-components PBX[J]. Explo Shock Waves, 2007, 27(2): 109-115. DOI:10.11883/1001-1455(2007)02-0109-07 |
| [15] |
Ma X F, Xiao J J, Zhu W H, et al. Computational study of structure and performance of four constituents HMX-based composite material[J].
J Mol Struct (Theochem), 2008, 851: 22-29. DOI:10.1016/j.theochem.2007.10.044 |
| [16] |
于艳春, 朱伟, 肖继军, 等. 四组份高能体系结合能和力学性能的分子动力学模拟[J].
化学学报, 2010, 68(12): 1181-1187. YU Yan-cun, ZHU Wei, XIAO Ji-jun, et al. Molecular dynamics simulation of binding energies and mechanical properties systems with four components[J]. Acta Chim Sinica, 2010, 68(12): 1181-1187. |
| [17] |
Xiao J J, Wang W R, Chen J, et al. Study on structure, sensitivity and mechanical properties of HMX and HMX-based PBXs with molecular dynamics simulation[J].
Comput Theor Chem, 2012: 999, 21-27. |
| [18] |
Xiao J J, Zhao L, Zhu W, et al. Molecular dynamics study on the relationships of modeling, structural and energy properties with sensitivity for RDX-based PBXs[J].
Sci China Ser B, 2012, 55(12): 2587-2594. DOI:10.1007/s11426-012-4797-1 |
| [19] |
Xiao J J, Li S Y, Chen J, et al. Molecular dynamics study on correlation between structure and sensitivity of defective RDX crystals and their PBXs[J].
J Mol Model, 2013, 19(2): 803-809. DOI:10.1007/s00894-012-1607-9 |
| [20] |
Andersen H C. Molecular dynamics simulations at constant pressure and/or temperature[J].
J Chem Phys, 1980, 72(4): 2384-2393. DOI:10.1063/1.439486 |
| [21] |
Ewald P P. Evaluation of optical and electrostatic lattice potentials[J].
Ann Phys, 1921, 64: 253-287. |
| [22] |
Allen M P, Tildesley D J.
Computer Simulation of Liquids[M]. Oxford: Oxford University Press, 1989 |
| [23] |
肖鹤鸣, 许晓娟, 邱玲.
高能量密度材料的理论设计[M]. 北京: 科学出版社, 2008.
XIAO He-ming, XU Xiao-juan, QIU Lin. Theoretical Design of High Energy Density Materials[M]. Beijing: Science Press, 2008 |
| [24] |
肖鹤鸣, 朱卫华, 肖继军, 等. 含能材料感度判别理论研究——从分子、晶体到复合材料[J].
含能材料, 2012, 20(5): 514-527. XIAO He-ming, ZHU Wei-hua, XIAO Ji-jun, et al. Theoretical studies on sensitivity criterion of energetic materials——from molecules, crystals, to composite materials[J]. Chinese Journal of Energetic Materials(Hanneng Cailiao), 2012, 20(5): 514-527. |
| [25] |
费业泰, 赵静. 材料热膨胀系数的发展与未来分析[J].
中国计量学院学报, 2002, 13(4): 259-263. FEI Ye-tai, ZHAO Jing. Analysis of the development and future tendency of material thermal expansion coefficient[J]. J China Univ Metrol, 2002, 13(4): 259-263. |
| [26] |
Qiu L, Xiao H M, Zhu W H, et al. Ab initio and molecular dynamics studies of crystalline TNAD (trans-1, 4, 5, 8-Tetranitro-1, 4, 5, 8-tetraazadecalin)[J].
J Phys Chem B, 2006, 110: 10651-10661. |
| [27] |
肖继军, 黄辉, 李金山, 等. HMX热膨胀系数的分子动力学模拟研究[J].
含能材料, 2007, 15(6): 622-625. XIAO Ji-jun, HUANG Hui, LI Jin-shan, et al. A MD simulation study of the coefficients of thermal expansion for β-HMX crystal[J]. Chinese Journal of Energetic Materials(Hanneng Cailiao), 2007, 15(6): 622-625. |
| [28] |
Guillermet A F. The pressure dependence of the expansivity and of the Anderson-Grüneisen parameter in the murnaghan approximation[J].
J Phys Chem Solids, 1986, 47: 605-607. DOI:10.1016/0022-3697(86)90163-0 |
| [29] |
严祖同, 孙振华. Anderson-Grüneisen参数、热膨胀系数与压强的普遍关系[J].
物理学报, 1989, 38(10): 1634-1641. YAN Zu-tong, SUN Zhen-hua. The pressure dependence of the expansivity and of the Anderson-Grüneisen parameter in the general condition[J]. Acta Phys Sinca, 1989, 38(10): 1634-1641. DOI:10.7498/aps.38.1634 |
| [30] |
Swenson R J. Comments on viral theorems for bounded systems[J].
Am J Phys, 1983, 51: 940-942. DOI:10.1119/1.13390 |
| [31] |
Watt J P, Davies G F, O′Connell R J. The elastic properties of composite materials[J].
Reviews of Geophysics, 1976, 14: 541-563. DOI:10.1029/RG014i004p00541 |
| [32] |
Pugh S F. Relation between the elastic moduli and the plastic properties of polycrystalline pure metals[J].
Phil Mag, 1954, 45: 823-843. DOI:10.1080/14786440808520496 |
| [33] |
Pettifor D G. Theoretical predictions of structure and related properties of intermetallics[J].
Mater Sci Technol, 1992, 8: 345-349. DOI:10.1179/mst.1992.8.4.345 |
The models of a five-component system ((PEG/NG/BTTN)/AP/HMX) with four sets of proportion and a six-component system ((PEG/NG/BTTN)/AP/HMX/Al) were designed and established. The safety performance, thermal expansion and mechanical properties of these complicated systems were explored by molecular dynamics simulation.