液氧相容改性环氧树脂体系低温力学性能分子动力学模拟*

高 畅，董旭峰，严 佳，徐 浩, 武湛君

(1. 大连理工大学 材料科学与工程学院，辽宁 大连 116024 2. 大连理工大学 航空航天学院，辽宁 大连 116024)

0 引 言

近年来，轻质纤维增强复合材料在航空航天工程中吸引了广泛关注[1]。然而，航空航天复合材料的安全性与完整性收到日益恶劣的服役环境挑战。例如，在液氧复合材料贮箱中，材料需同时面临低温(-183 ℃)与液氧不相容的挑战(即材料与液氧之间产生的燃烧或爆炸反应)[2-5]。由于复合材料中的纤维具有优异的物理性能和化学稳定性，因此聚合物基体的性能对于评估复合材料整体性能至关重要。环氧热固性树脂在机械性能、易加工性和低成本等方面显示出独特优势。最近研究结果还表明，向环氧树脂中添加各种类型的阻燃剂(FR)分子可以大大改善材料的液氧相容性[6-9]。其中，添加含磷阻燃剂在阻燃效率和无毒性方面具有明显优势。Satio等[10]合成了9，10-二氢-9-氧代-10-磷杂菲-10-氧化物(DOPO)的含磷阻燃剂，具有良好的热稳定性与树脂相容性。使用DOPO作为中间体，Cho等[11]合成了10-(2,5-二羟基苯基)-10-氢-9-氧杂-10-磷杂菲-10-氧化物(ODOPB)，表现出更高的热稳定性。李家亮等[12-14]通过在环氧分子上分别接枝ODOPB官能团，对环氧树脂进行了改性，成功通过标准液氧相容性试验。

然而，接枝阻燃剂将改变环氧热固性树脂中交联网络的分子结构，引起力学性能变化。实验结果表明，磷含量高于某个阈值，材料的模量、韧性等性能都会出现退化[14]。因此，同时提高环氧热固性树脂的力学性能和液氧相容性是一项重要而具有挑战性的任务。

深入了解改性环氧树脂在承载过程中的微观机理，不仅可以揭示宏观性能的潜在机制，而且还可以通过采取实验措施来改善材料性能。然而由于现有的分子尺度实验表征技术准确性有限，相关的研究仍然很少。分子动力学(MD)模拟为揭示聚合物结构-性能关系提供了有利工具[15-17]。本文利用分子动力学(MD)模拟方法，研究了10-(2，5-二羟基苯基)-10-氢-9-氧杂-10-磷杂菲-10-氧化物(ODOPB)改性环氧树脂的低温力学性能。首先采用两步反应方案生成改性的环氧分子，然后生成交联网络。在不同的磷含量下，计算了材料常温与低温的弹性模量，并进行了拉伸变形模拟。联合使用通用非反应性DREIDING力场[18]和ReaxFF反应力场[19]，进行了截断键长校准，使键断裂模拟具有高精度和计算效率。研究发现，材料体系中链段长度分布特征在控制力学性能方面起着主要作用。更重要的是，在环氧树脂改性阶段采用的链长调节能够显著提高材料的承载与变形能力，为工程实践中设计和制备高强、高韧、液氧相容环氧树脂基复合材料提供了有价值的参考。

1 改性交联环氧树脂体系分子动力学建模方法

使用Materials Studio进行MD建模，ODOPB，双酚A二缩水甘油醚(DGEBA)与4,4-二氨基二苯基甲烷(DDM)固化剂单体的化学结构如图1所示。为了与实验中的化学改性过程保持一致，采用两步反应方案，首先生成改性环氧树脂，然后生成交联热固性树脂。具体步骤为：(1)环氧树脂分子ODOPB枝接改性。将ODOPB和活化DGEBA环氧树脂按照5种磷含量，即0、0.5%、1%、1.5%和2%(质量分数)，放在周期性盒子中进行混合，每个模型中的原子数控制在50 000左右。MD模拟在NVT系综下进行，为加速反应速率，温度控制在600 K。使用开源软件包LAMMPS和DREIDING力场进行模拟，ODOPB和DGEBA之间的反应使用Abbott等开发的Polymatic算法进行模拟[20]。由于ODOPB和DGEBA单体各自具有两个反应点位(即图1中ODOPB单体两端的-OH基团与DGEBA单体两端的-CH2基团)，因此改性后的环氧树脂分子链长(用n表示)将显现出较大的分散性。为了模拟链长调控(chain length regulation(CLR))，在2%(质量分数)磷含量下另外构建MD模型，使其仅包括n=0和1两种分子链长。图2(a)和(b)分别给出了未经调控(unregulated)与CLR体系中不同分子链长的数量比例分布。可以观察到，相比于未调控体系，CLR体系中包含更少的超短链(即n=0的未改性分子链)，以及更多n=1的分子链。这是因为较短分子链中的ODOPB分子会消耗更多的DGEBA分子。(2)交联固化反应。在不同磷含量下，将ODOPB改性环氧树脂分子按照化学计量在周期盒子中与DDM固化剂分子混合。原子数控制在60 000左右。使用Polymatic算法在600 K温度下模拟交联过程，交联度控制在90%，与实验情况相当。

图1 ODOPB、DGEBA、DDM单体与改性环氧树脂分子结构图(虚线标明了参与反应的基团)Fig 1 Monomers of ODOPB, DGEBA, DDM and the molecular structure of modified epoxy resin (the reaction functional groups are marked by dashed circles)

表1 磷含量2%(质量分数)下未经链长调控与链长调控ODOPB改性环氧树脂分子的数量比例分布

2 结果与讨论

基于ODOPB改性交联环氧树脂分子动力学模型，开展材料体系的常低温力学性能计算。

2.1 常低温弹性模量

首先将MD体系温度从600 K降温至300 K(室温)与90 K(接近液氧温度),降温速率控制在0.05 K/ps。采用NPT系综，在标准大气压条件下对材料体系分别进行300与90 K温度条件下的拉伸模拟，应变率为5×108，最大拉伸应变为5%。图2(a)与(b)分别给出了0%与2%(质量分数)磷含量下交联环氧树脂材料体系的拉伸曲线。可以看出，5%应变范围内材料应力与应变线性关系显著。而低温下的应力水平明显高于常温。通过比较图2(a)与(b)，可以看出在相同温度下，两种磷含量材料体系的应力水平并无显著差异。

图2 (a)0%与(b)2%(质量分数)磷含量下材料体系常低温拉伸曲线Fig 2 Room-temperature and low-temperature tensile curves of material systems containing (a) 0 wt% and (b) 2 wt% P content, respectively

通过应力-应变关系曲线的斜率计算了常低温下的材料弹性模量。图3给出了5种不同磷含量体系的常低温模量分布。从图中可以看出，材料弹性模量数值随磷含量增加的变化趋势并不明显，即在同一温度下，含磷量的不同并未显著改变材料弹性模量。另一方面，温度对材料模量的影响十分突出，在90 K低温环境下(接近液氧温度)，材料的弹性模量对比常温环境显著增加，增幅达30%以上。以上现象与实验观测结果基本一致[5]。

图3 不同磷含量材料体系的常温与低温弹性模量Fig 3 Room-temperature and low-temperature elasticity modulus at different P contents

2.2 低温拉伸断裂行为模拟

为了准确复现材料破坏过程中的化学断裂过程，联合使用反应与非反应力场进行截断键长校准。首先，分别使用非反应DREIDING力场和反应型ReaxFF力场计算小型材料体系(约8 000原子)的应力应变响应。由于ReaxFF力场有效融合了元素的化学反应参数，获得的响应被视为准确结果，并用作参考数据。通过不断调整DREIDING体系中的截断键长，直到生成的应力应力响应与参考数据高度一致。然后，基于DREIDING的截断键长模拟大体系(约60 000原子)的拉伸和断裂行为。校准过程可以极大地提高DREIDING力场断键复现的准确性，同时保证计算效率，解决ReaxFF力场计算耗时的难题。8 000原子体系下，基于两种力场的未改性环氧树脂体系拉伸曲线如图4所示，可以看到应力应变曲线的吻合很好。

图4 分别采用Dreiding力场与ReaxFF力场的拉伸曲线Fig 4 Tensile curves calculated using Dreiding and ReaxFF force field

接下来基于DREIDING力场，采用NVT系综计算了0，1%，2%(质量分数)3种不同磷含量材料90 K温度下的大应变范围拉伸曲线，应变率为5×108。在NVT系综下，材料体系沿拉伸方向伸长，但其横截面在拉伸过程中保持不变。这种应力状态类似于I型裂纹张开试件实验中，裂纹尖端的三向应力状态。图5给出2%(质量分数)改性树脂在100%，200%与300%应变下原子体系的拉伸形貌。可以看到，随应变增加，体系中的分子链被不断拉伸，并逐渐出现了内部孔洞以及分子链的断裂现象。3种不同磷含量的应力-应变拉伸曲线结果如图6 所示。材料的拉伸响应随磷含量的增加发生明显变化。在后屈服阶段，3种材料的平台应力水平(即S)相差不明显，但最大应力Smax存在很大差异。其中，未改性体系和1%(质量分数)磷含量的改性体系相比，Smax差别不大。而随磷含量进一步增加到2%(质量分数)，Smax值明显变小，说明材料的承载能力随磷含量增高呈现下降趋势。

图5 90 K温度下2%(质量分数)交联环氧树脂原子体系100%，200%与300%应变水平下的拉伸形貌Fig 5 Snapshots of the crosslinked epoxy resin systems at 100%, 200% and 300% strains under 2 wt% P content and 90 K

图6 0、1%与2%(质量分数)磷含量下材料的低温拉伸曲线Fig 6 Low-temperature tensile curves at 0, 1 wt% and 2 wt% P contents

需要指出的是，应力-应变曲线的MD仿真结果与常规热固性树脂试件实验曲线存在较大差异。例如，MD体系拉伸的应力水平和变形水平都远远高于实验数据，其主要原因是仿真与实验的尺度差异。MD模型的空间尺度在纳米量级，因此模型中无任何材料缺陷，其仿真结果可近似理解为理想情况下无缺陷材料的极限力学响应。虽然存在仿真与实验的差异，但MD的仿真在定性分析角度仍具有重要意义，例如可有效实现不同材料组分下力学性能的比较和筛选，有效揭示材料变形过程中的分子链行为等。

2.3 链长调控对宏观力学性能的影响及其微观机理分析

图7(a)给出了2%(质量分数)磷含量下经过链长调控(CLR)与未经链长调控(unregulated)材料的低温应力-应变曲线。可以明显发现，CLR体系具有更大的应力峰值(Smax)以及更强的变形能力(应力峰值对应的应变水平)。材料拉伸性能的差异归因于交联网络结构的不同。受拉伸载荷的影响，交联网络中较短的分子链段先于较长链段被拉伸(直至断裂)。与CLR体系相比，未经链长调控的体系在链长分布分散性更大。在材料的同一应变水平下，不同链长无法有效协同承载。相比之下，CLR体系链长分布更为集中，能够更有效地协同传递载荷，进而提高材料的应力和变形水平。为了进一步量化拉伸过程中分子链段的伸长行为，计算了环氧树脂链段的平均端对端距离(end-to-end distance，Le)，并将调控与未调控体系进行了比较，如图7(b)所示。如图所示，Le值随应变增大而增加，代表材料中的环氧树脂链段从最初的卷曲与松弛状态逐渐被拉伸至断裂。在拉伸过程中，链长调控体系的Le值比未调控体系增加的幅度更大，说明调控体系中的链段整体伸长行为更为充分，不同长度链段伸长的协同效应也更为显著。

图7 (a) 2%(质量分数)磷含量下改性与未改性材料的低温拉伸曲线与(b)随应变增大，两种体系环氧树脂链段的平均端对端距离(Le)变化Fig 7 (a) Low-temperature tensile curves of CLR and unregulated systems at 2 wt%, and (b) variations of the Le values of the two material systems along tensile deformation

3 结 论

通过两步反应法，建立了ODOPB改性交联环氧树脂分子动力学模型，计算了不同磷含量材料体系的常低温力学性能，为探究改性环氧树脂力学性能的微观机理，以及高性能材料的结构设计与性能优化提供了有价值的参考。

磷含量的增加对材料的拉伸模量并无明显影响，而温度的降低使材料的弹性模型发生明显增大。在材料的拉伸与破坏研究中，通过ReaxFF力场得到DREIDING力场的断键截断半径，从而在保证计算精度的同时，大幅提高了分子动力学计算效率。从拉伸曲线可以看出，在未改性体系中，磷含量的增加显著降低了材料的承载能力。而在同一磷含量(2%，质量分数)下，通过链长调控能够较显著改善材料的力学性能，其原因是调控后材料体系中的链长能够在拉伸过程中表现出更好的协同承载效应。研究成果为液氧相容、高强高韧树脂基体及其复合材料的设计与制备提供了有益参考。