ε-CL-20/含能黏结剂复合体系结合能及力学性能的模拟

陶 俊， 王晓峰， 赵省向， 王彩玲， 刁小强， 韩仲熙

(西安近代化学研究所， 陕西 西安 710065)

1 引 言

实现武器系统的高效毁伤，可以采用先进的战斗部设计和新型特种材料，但是更为直接有效的技术途径是采用更高能量的装药，因此高能量密度化合物的合成与应用研究受到世界各国的高度重视[1-3]。六硝基六氮杂异伍兹烷(HNIW，CL-20)的能量输出比奥克托今(HMX)高10%～15%，以其为基(主体炸药)添加少量高聚物黏结剂，形成的高聚物黏结炸药是最有应用前景的高能量密度材料之一[4-5]。含能黏结剂已经广泛应用于炸药、推进剂及发射药中，在提高能量水平的同时可以降低危险等级，提高使用的安全性。

含能黏结剂粘结CL-20混合炸药的研究具有极大的吸引力。Wagstaff[6]制备了一种含能增塑剂包覆的ε-CL-20 混合炸药，并称含能增塑剂的添加可使ε-CL-20 的感度低于高氯酸铵，这些增塑剂包括丁三醇三硝酸酯(BTTN)、三羟甲基乙烷三硝酸酯(TMETN)、聚叠氮缩水甘油醚(GAP)和双(2, 2-二硝基丙醇)缩乙醛/双(2, 2-二硝基丙醇)缩甲醛(BDNPA/F) 等； 美国采用乙酸丁酸纤维素(CAB)/ BDNPA/F或聚邻苯二甲酸二丁酯(HyTemp)/ GAP作为黏结剂/增塑剂体系并与混合炸药LX-19进行了性能比较，其气流冲击载荷(ABL)撞击感度提高了22.5 cm，摩擦感度也有一定程度的降低[7]。目前，炸药配方的设计及确定还主要依赖于实验，主要有以下缺点: 实验周期长、实验费用高、难以预测配方性能、存在安全问题。运用分子动力学模拟高聚物粘结炸药及其组分的结构和性能有助于指导其配方设计[8-11]。

本研究构建了不同软硬段比例((3-叠氮甲基-3-甲基氧丁环(AMMO): 3，3-二叠氮甲基氧丁环(BAMO))的含能高聚物及其与ε-CL-20/含能黏结剂的PBX体系，对含能黏结剂以及ε-CL-20/含能黏结剂复合体系进行分子动力学模拟，对其性能进行分析和比较，相关研究未见报道，研究结果为含能黏结剂粘结ε-CL-20混合炸药的配方设计提供理论指导和设计思路。

2 计算方法

2.1 力场选择

选择Compass力场对ε-CL-20/含能黏结剂复合体系进行模拟。选择原因如下: 首先，Compass力场中多数力场参数的调试确定都基于从头算数据，它是通过从头计算方法,获得分子内的键参数，同时又采用以凝聚态分子动力学为基础的经验方法优化范德华非键合参数，Compass力场是当前较为完善的分子力场； 其次，Compass 对ε-CL-20分子的优化构型与量子力学优化构型符合得很好，且基于Compass力场预测的ε-CL-20晶胞结构也与实验相符[8,12-13]。

2.2 含能高聚物的构建

构建了不同软硬段比例(BAMO∶AMMO)的含能高聚物，控制软段的长度和硬段的长度，使得PBAMO(12)、PBAMO(9)-g-PAMMO(2) 、PBAMO(5)-g-PAMMO(7)、PAMMO(17)占含能黏结剂/CL-20复合体系的质量百分数基本保持一致，分别为: 8.07%、8.28%、8.03%及7.97%。PBAMO(9)-g-PAMMO(2) 、PBAMO(5)-g-PAMMO(7)软段硬段之间的氨酯键齐聚醇扩链剂[14]结构如图1所示，均聚合物及共聚物的分子式如图2所示。

图1 扩链剂的结构

Fig.1 Structure of the chain extender

图2 PBAMO、PAMMO及PBAMO(m)-g-PAMMO(n)结构式

(m代表PBAMO中BAMO个数，n代表PAMMO中AMMO个数，g代表嵌段共聚物)

Fig.2 The structure formulas of PBAMO, PAMMO and PBAMO (m)-g-PAMMO (n)

(mrepresents the number of BAMO in PBAMO,nrepresents the number of AMMO in PAMMO,grepresents the block copolymer)

用Materials studio软件中的Discover模块，以Compass力场进行MD模拟，选取NVT系综，在真空条件下，采用Anderson恒温器，温度为298 K，时间步长为1 fs，总模拟时间为300 ps，获得的最终构象为高聚物链的平衡构象。

2.3 ε-CL-20基PBXs的构建

ε-CL-20晶体结构取自X-ray衍射结果[15]，构建的原胞如图3所示。由MS软件构建其(3×3×2)晶胞模型，晶胞原子总个数为2592个。因ε-CL-20的(001)晶面分子堆积较为紧密，与各黏结剂作用较强[16]，故本研究将超晶胞沿(001)晶面方向切割，真空层设为20 Å。将含能高聚物的平衡构象分别加入该真空层中，并尽可能使接近ε-CL-20分子，从而PBX初始构型，经过能量优化和动力学模拟在(Compass力场下进行NVT系综MD模拟，温度设为298 K，选择Andersen控温方法，步长1 fs，模拟步数为20万步，前10万步用于平衡，后10万步用于统计分析，每50步保存一次轨迹文件，用于分析结合能)得到其平衡(体系的平衡可由温度和能量的同时平衡来确定，当温度和能量在5%～10%范围内波动即可认为体系已达到平衡)构型如图4所示。

2.4 ε-CL-20基PBXs的力学性能模拟及计算

将2.3中得到的平衡结构，边压缩边进行动力学模拟，经过多次压缩，得到密度接近ε-CL-20基PBXs理论密度的平衡结构[8,16]，其中，PBAMO(12)、PBAMO(9)-g-PAMMO(2)、PBAMO(5)-g-PAMMO(7)、PAMMO(17)与ε-CL-20复合体系平衡结构的密度分别为: 1.897，1.891，1.895，1.892 g·cm-3，平衡结构中总原子个数分别为2835、2833、2843及2901个。通过MS对密度接近含2.2中的质量分数含能黏结剂的ε-CL-20基PBXs理论密度的平衡结构进行形变量为0.05%的拉伸与纯切形变操作，得到弹性系数应变各方向应力分量—弹性系数Cij(i，j=1～6)的矩阵然后分析模型并计算可得拉伸模量E、剪切模量G和体积模量K，以及泊松比γ等力学性能参数。

图3ε-CL-20单晶胞

Fig.3 Single cell ofε-CL-20 crystal

3 结果与讨论

3.1 含能黏结剂与ε-CL-20的相互作用及相容性

PAMMO及PBAMO是两种较为常见的含能黏结剂，PAMMO-g-PBAMO的含能黏结剂由PBAMO的硬段及PAMMO的软段通过与扩链剂反应得到，由于其良好的能量性质及力学性能，在发射药和推进剂中得到了广泛的应用。

图4为ε-CL-20/含能黏结剂复合体系的平衡结构，从图4中可以看出，相比于其他3种黏结剂体系，PAMMO由于分子柔性好，能较为平坦地铺展在ε-CL-20晶体表面(图4d)。

结合能是表征共混体系组分间相互作用强度的特征参数对共混体系的力学性能有着重要影响。聚合物黏结剂在ε-CL-20晶体表面的结合能可表示为[17-19]:

Ebind=-Einter=-[Etotal-(Ecell+Epolymer)]

Ebind是聚合物与晶面的结合能，通常Ebind越大，聚合物与晶面的相互作用越强；Etotal是由平衡结构求得的单点能；Ecell为去掉高聚物后计算所得ε-CL-20的单点能；Epolymer为去掉ε-CL-20分子后计算所得高聚物的单点能。

表1列出了ε-CL-20与不同软硬段比例黏结剂分子之间的结合能，可以发现四个体系的结合能排序为:ε-CL-20/PBAMO(9)-g-PAMMO(2)>ε-CL-20/PBAMO(5)-g-PAMMO(7)>ε-CL-20/PAMMO(17) >ε-CL-20/PBAMO(12)，除了ε-CL-20/PBAMO(12)外，其余三个体系的结合能均为正值。结合能的大小决定了高聚物与ε-CL-20晶体作用的强度，结合能越大，作用强度越大。吸附能为正值表明吸附稳定； 吸附能为负表明聚合物与晶体的吸附不稳定、可逆，存在着吸附与解吸附现象。PBAMO与ε-CL-20的结合能之所以为负值，可能是因为PBAMO的分子链中的BAMO结构单元具有两个强极性、大体积的叠氮甲基，空间位阻较大，使分子链的柔顺性受到影响，不易在ε-CL-20晶体较好的铺展以发生分子间相互作用。四个体系中PBAMO(9)-g-PAMMO(2)与ε-CL-20晶体的结合能最大，在PBAMO分子链上引入少量的PAMMO分子，可改善分子链的柔顺性，使得整个分子链可以变形和运动，易于与ε-CL-20晶体接触并发生相互作用。

对于PBXs而言，高聚物黏结剂是惰性成分，一般不会使主体炸药发生化学变化，故它们之间只有物理相容性，物理相容性的好坏，与主体炸药和黏结剂之间的相互作用有关。结合能(Ebind)愈大表明形成的PBX越稳定，反映了主体炸药和黏结剂之间的彼此相容性越好。因此，从模拟结果可以看出四种含能黏结剂与ε-CL-20的相容性排序为:ε-CL-20/PBAMO(9)-g-PAMMO(2)>ε-CL-20/PBAMO(5)-g-PAMMO(7)>ε-CL-20/PAMMO(17)>ε-CL-20/PBAMO(12)。

图4ε-CL-20/含能黏结剂体系的平衡结构

Fig.4 The equilibrium structures ofε-CL-20 / energetic binder

表1中列出了范德华力及静电相互作用对各部分能量的贡献，图5为范德华力对体系结合能的贡献与黏结剂链中BAMO摩尔分数的关系图。从图5中可以看出随着黏结剂链中BAMO摩尔分数的增加，范德华力对ε-CL-20与黏结剂之间的结合能贡献越大。PAMMO分子的范德华力对结合能的贡献为64.6%，PBAMO分子的范德华力对结合能的贡献接近100%，这表明: 四个体系的结合能主要是由黏结剂与ε-CL-20晶体之间的范德华力作用引起的，BAMO结构单元基本全靠范德华力与ε-CL-20晶体发生相互作用，AMMO结构单元与ε-CL-20晶体之间的相互作用则以范德华力为主，还有少量的静电力等其他相互作用力。

表1 CL-20与黏结剂分子之间的结合能

Table 1 The binding energy between CL-20 and binder kJ·mol-1

Note:Etotalis the single-point energy of the equlibrium strcture，ECL-20is the single-point energy ofε-CL-20，Epolymeris the single-point energy of polymer，Eis the total energy of each structure，vdW is the energy of each structure obtained by vdW interaction, electrostatic is the energy of each structure obtained by electrostatic interaction.

图5 范德华力对结合能的贡献与黏结剂链中BAMO摩尔分数的关系

Fig.5 The relationship between the contribution of van der waals force on the binding energy and the molar fraction of BAMO on the binder chain

3.2 ε-CL-20/含能黏结剂复合体系的力学性能

图6为经动力学模拟得到的密度接近体系理论密度的ε-CL-20/含能黏结剂复合体系的平衡结构，表2为图6中的平衡结构经过力学性能模拟得到的弹性系数，表中未列入的弹性系数的数值等于或者接近于0。材料的应力是应变的线性组合，其线性组合系数为应变各方向应力分量——弹性系数Cij(i,j=1～6)的矩阵，原则上材料的所有力学性能均可从它的弹性系数矩阵得到[20-21]。

通过对弹性系数Cij(i，j=1～6)进行计算，得到拉梅系数λ和μ，计算方法为:

材料在拉伸、剪切和压缩情况下应力与应变的比值，即拉伸模量E、剪切模量G和体积模量K，以及泊松比γ均可用λ和μ表示为:

式中，E为拉伸模量，GPa；G为剪切模量，GPa；K为体积模量，GPa；γ为泊松比。由于应变能的存在，有Cij=Cji，即矩阵对称，因此对于极端各向异性体，只有21个独立的弹性系数。随着物体对称性的提高，独立的弹性系数减少，对各向同性体，只有两个独立的弹性系数C11和C12。从表2可以看出，ε-CL-20的三组弹性系数C11、C22、C33，C44、C66，C12、C13；ε-CL-20/含能黏结剂体系的三组弹性系数C22、C33，C44、C55、C66，C12、C13的组内弹性系数都很接近。

从表2可以看出，与ε-CL-20相比，三组ε-CL-20/含能黏结剂复合体系的弹性系数及模量都发生了不同程度的下降，表明含能黏结剂的加入均能够减弱ε-CL-20的刚性，而使其塑性增强，在不同的受力情况下更易发生弹性形变。五个体系拉伸模量及剪切模量的排序为:ε-CL-20>ε-CL-20/PBAMO(5)-g-PAMMO(7)>ε-CL-20/PBAMO(12)>ε-CL-20/PBAMO(9)-g-PAMMO(2)>ε-CL-20/PAMMO(17)。相对于其他三个ε-CL-20/含能黏结剂复合体系，ε-CL-20/PAMMO体系的弹性系数、模量值最小，这是由于相对于PBAMO，PAMMO分子链较为柔顺、刚性小，加入ε-CL-20中可以使得整个体系的塑性增强。

图6 接近理论密度的ε-CL-20/含能黏结剂复合体系平衡结构

Fig.6 Equilibrium structures ofε-CL-20/energetic binder systems with the densities close to the theoretical values

表2ε-CL-20/含能黏结剂的弹性系数和模量

Table 2 The elasticity coefficient and modulus ofε-CL-20/energetic binder

parametersε-CL-20ε-CL-20/PBAMOε-CL-20/PBAMO(9)-g-PAMMO(2)ε-CL-20/PBAMO(5)-g-PAMMO(7)ε-CL-20/PAMMOC1171.88 1.03-1.62 1.07 0.22C2254.93 46.2649.85 47.40 43.85C3370.30 42.6436.15 43.43 46.20C4412.62 5.864.20 3.97 5.19C55-4.85 5.064.59 5.84 0.056C6611.9 5.847.12 8.16 6.93C1223.66 15.0915.98 16.89 15.64C1330.07 17.1916.52 16.50 16.90C2343.69 31.5631.62 31.33 33.65E/GPa18.95 15.5114.59 18.65 11.55K/GPa56.96 29.2021.05 16.51 24.68G/GPa6.56 5.605.30 5.99 4.06γ0.44 0.390.38 0.38 0.42(C12-C44)/GPa11.039.2311.78 12.92 10.44

Note:Eis tensile modulus,Kis bulk modulus,Gis shear modulus,γis Poissom′s ratio,C12-C44is Cauchy.

在弹性系数中，C12只与张应力有关，C44只与切应力有关。当C12大于C44时，材料易发生剪切形变，此时材料的延展性较好； 当C12小于C44时，材料易发生弹性形变，材料显脆性。因此可用C12-C44(柯西压，Cauchy)评估材料的延展性和脆性[20-21]。从表2中可知，ε-CL-20及ε-CL-20/含能黏结剂复合体系的C12-C4值都大于零，表现为韧性，其中PBAMO及PAMMO两种均聚物的加入并不能有效的改善ε-CL-20的延展性，相比较而言，ε-CL-20/PBAMO(9)-g-PAMMO(2)、ε-CL-20/PBAMO(5)-g-PAMMO(7)两个体系的C12-C44值比起其他几个体系相对较大。因此，含能热塑性弹性体的加入改善了ε-CL-20的力学性能，使其延展性增强。

从力学性能的分析可以看出，相比BAMO、AMMO的两种均聚物，具有一定软硬段比例的含能热塑性弹性体并不等同于相同分子量PBAMO和PAMMO的物理混合，含能热塑性弹性体具有均聚物所不具有的一些力学性能。

3.3 含能黏结剂的能量性能

基团加和法计算ETPE的生成焓与实际符合较好，能较好的预测热塑性弹性体的生成焓[23]。因此，决定采用Van Krevelen和Chemin基团估算法对含能黏结剂的生成焓进行评估[24]，不同基团对聚合物分子生成焓的贡献见表3。

表3 不同基团对聚合物分子的生成焓的贡献[22]

Table 3 Contributions of groups to enthalpy of formation of polymer

groupΔHθf/kJ·mol-1groupΔHθf/kJ·mol-1—CH3-46-3—CH2—-22—N3+356C+20C-132COO-337100O+58ONO2-88CH2-132CH+58

因为不同的实验方法合成的聚合物分子量可能有差别，因此这里用单位质量物质的生成焓来描述和比较物质生成焓的大小。图7为不同软硬段比例含能黏结剂的生成焓，从图7中可以看出生成焓的排序为: PBAMO(12)>PBAMO(9)-g-PAMMO(2)>PBAMO(5)-g-PAMMO(7)>PAMMO(17)。随着高分子链中BAMO含量的增加，单位质量含能黏结剂的生成焓逐渐增加。这是因为BAMO中含有两个叠氮基团，AMMO只有一个叠氮基团，叠氮基团对生成焓的贡献很大(为+356 kJ·mol-1)，多出的一个叠氮基团能显著提升体系的能量。

硝基类化合物由于硝基的生成焓为负值，最终导致以其为软段的ETPE的生成焓为负值。单从能量角度考虑，在设计含能热塑性弹性体时，最好选用叠氮类化合物为软段。BAMMO、AMMO的嵌段共聚物在能量方面还是相当有优势的，软硬段都含有对生成焓的贡献很大的叠氮基团。

图7 不同软硬段比例含能黏结剂的生成焓

Fig.7 The enthalpy of formation of energetic binders with different hard and soft segments proportion

另外，可以发现BAMO的分子链中的BAMO结构单元具有两个强极性、大体积的叠氮甲基与ε-CL-20的物理相容性不好；力学性能的模拟中PAMMO分子链较为柔顺、刚性小，加入ε-CL-20中可以使得整个体系的塑性增强。但是，BAMO含有两个叠氮基团能显著提升体系的能量。因此，含能黏结剂存在着能量性能和力学性能及其他性能无法综合平衡的困难，具体应用时可根据应用对象的要求进行选择。

4 结 论

对ε-CL-20与四种不同AMMO/BAMO比例的含能黏结剂构成的高聚物粘结炸药的结构和性能进行了计算模拟，得到以下结论:

(1)ε-CL-20与软硬段比例不同的四种黏结剂之间的结合能的排序为:ε-CL-20/PBAMO(9)-g-PAMMO(2)>ε-CL-20/PBAMO(5)-g-PAMMO(7)>ε-CL-20/PAMMO(17) >ε-CL-20/PBAMO(12)。PBAMO的分子链中的BAMO结构单元的叠氮甲基，空间位阻较大，不易在ε-CL-20晶体较好的铺展以发生分子间相互作用； 分子链上引入少量的PAMMO分子，可改善分子链的柔顺性，易于与ε-CL-20晶体接触并发生相互作用。

(2) BAMO结构单元基本全靠范德华力与ε-CL-20晶体发生相互作用，AMMO结构单元与ε-CL-20晶体之间的相互作用则以范德华力为主，还有少量的静电等其他相互作用力。

(3)ε-CL-20及ε-CL-20/含能黏结剂复合体系平衡构型力学性能较为接近各向同性； 含能黏结剂的加入均能够减弱ε-CL-20的刚性，而使其塑性增强；ε-CL-20/PAMMO体系的弹性系数、模量值最小，这是由于相对于PBAMO，PAMMO分子链较为柔顺、刚性小。

(4) 不同软硬段比例含能黏结剂的生成焓的排序为: PBAMO(12)>PBAMO(9)-g-PAMMO(2)>PBAMO(5)-g-PAMMO(7)>PAMMO(17)，随着高分子链中BAMO含量的增加，单位质量含能黏结剂的生成焓逐渐增加。

参考文献:

[1] Agrawal P J. Recent trends in high-energetic materials[J].ProgEnergyCombustSci, 1998, 24: 1-30.

[2] Nedelko V V, Chukano N V, Raevski A V, et al. Comparative investigation of thermal decompsition of various modification of hexnitrohexaazaisowurtzitane (CL-20)[J].PropellExplosPyrotech, 2000, 25: 255-259.

[3] Zhang M X, Eaton P E, Gilardi R. Hepta and octanitrocubanes[J].AngewChemIntEd, 2000, 39(2): 401-404.

[4] Lobbecke S, Bohn M A, Pfeil A, et al. Thermal behavior and stability of HNIW (CL-20)[C]∥Proceedings of 29th Internationsl Conference of ICT, Karsruhe, 1998.

[5] Paval Vavra. Procedure for selection of molecular structures of explosives having high performance[C]∥Proceedings of 30th International Conference of ICT, Karsruhe, 1999.

[6] Wagstaff D C. Desensitization of energetic materials: GB, 2374867A[P]. 2002-10-30.

[7] Donald A，Geiss J．Additional characterization of high performance CL-20 formulation[C]∥Insensitive Munitions and Energetic Materials Technology Symposium. San Diego: NDIA，1999: 129-140.

[8] 许晓娟，肖继军，黄辉，等.ε-CL-20基PBX结构与性能的分子动力学模拟—HEDM理论配方设计初探[J].中国科学 B辑: 化学，2007，37(6): 556-563.

[9] Sewell T D, Menikoff R, Bedrov D, et al. A molecular dybamics simulation study of elstic properties of HMX[J].JChemPhys, 2003, 119: 7417-7426.

[10] Qiu L,Xiao H M, Zhu W H, et al. Ab initio and molecular dynamics study of crystalline TNAD (trans-1,4,5,8-tetrabitro-1,4,5,8-tetraazadecalin)[J].JPhysChemB, 2006, 110: 10651-10661.

[11] 杨小震. 分子模拟与高分子材料[M].北京: 科学出版社，2002: 1-10.

[12] Sun H, Ren P, Fried J R. The compass force field: parameterization and validation for phosphazenes[J].ComputTheorPolymSci, 1998, 8(1/2): 229-246.

[13] Sun H. Compass: An ab initio force-field optimized for condense-phase applications-overview with detials on alkanes and benzene compounds[J].JPhysChemB, 1998, 102: 7338-7364.

[14] 甘孝贤，李娜，卢先明，等.BAMO/AMMO基ETPE的合成与性能[J].火炸药学报，2008，31(2): 81-85.

GAN Xiao-xian, LI Na, LU Xian-ming，et al. Synthesis and properties of ETPE based on BAMO/AMMO[J].ChineseJournalofExplosives&Propellants, 2008, 31(2): 81-85.

[15] 赵信歧，施倪承.ε-六硝基六氮杂异戊兹烷的晶体结构[J]. 科学通报，1995，40: 2158-2160.

[16] Xu X J， Xiao J J, Zhu W, et al. Molecular dynamics simulations for pureε-CL-20 andε-CL-20-based PBXs[J].JPhysChemB, 2006, 110: 7203-7207.

[17] 徐茜.RDX/氟聚合物体系与PEG/硝酸酯体系分子动力学模拟[D].南京: 南京理工大学，2013.

XUN Qian. Molecular Dynamics Simulation Study on RBX/Fluoropolymer and PEG/Nitrate Composite[D]Nanjing: Nanjing University of Science & Technology, 2013.

[18] 翟頔.HNIW基高聚物粘结炸药结构和性能的MD模拟研究[D].太原: 中北大学，2012.

ZAI Ye. HNIW Based Polymer Bonded Explosive Structure and Properties of MD Simulation Study[D]. Taiyuan: North University of China, 2012.

[19] 满梅玲，陆春海，陈文凯，等.一氧化碳分子在Pt/t-ZrO2(101)表面的吸附性质[J].物理化学学报，2012，28(1): 51-57.

MAN Mei-ling, LU Chun-hai, CHEN Wen-kai, et al. Adsorption properties of CO molecules on Pt/t-ZrO2(101) surfaces[J].ActaPhysChimSin, 2012, 28 (1), 51-57.

[20] 郭昕，南海，齐晓飞，等. RDX和HMX晶体力学性能的分子动力学模拟及其撞击加载响应[J].含能材料，2013，21(4): 485-489.

GUO Xin, NAN Hai, QI Xiao-fei, et al. Molecular dynamics simulation on mechanical properties of RDX and HMX crystals and their impacting load response[J].ChineseJournalofEnergeticMaterials(HannengCailiao), 2013, 21 (4): 485-489.

[21] 朱伟. HMX、TATB基高聚物粘结炸药(PBX)力学性能的MD模拟研究[D].南京: 南京理工大学，2005.

ZHU Wei. Molecular dynamics simulations of HMX and TATB based polymer bonded explosives (PBXs) [D]. Nanjing: Nanjing University of Science & Technology, 2005.

[22] 夏露. 高能材料结构和性能的分子动力学模拟[D].苏州: 苏州大学，2008.

XIA Lu. Molecular dynamics simulations of the structures and properties of highly energetic materials[D].Soochow: Soochow University, 2008.

[23] 吕勇，罗运军，葛震.基团加和法估算含能热塑性聚氨酯弹性体的生成焓[J].含能材料，2009，17(2): 131-136.

Lü Yong, LUO Yun-jun, GE Zhen. Estimation of enthalpy of formation for energetic thermoplastic polyurethane elastomers by group additivity method[J].ChineseJournalofEnergeticMaterials(HannengCailiao), 2009, 17 (2): 131-136.

[24] 范克雷维纶.聚合物的性质: 性质的估算及其与化学结构的关系[M].许元泽，赵得赫，吴大诚,译. 北京: 科学出版社，1981: 336-340.