江 浩，岳 紅，劉 倩

（西北工業(yè)大學理學院，陜西 西安710129）

分子動力學模擬研究NR/BR力學性能和界面相互作用

江 浩，岳 紅，劉 倩

（西北工業(yè)大學理學院，陜西 西安710129）

運用分子動力學（MD）方法，模擬研究了天然橡膠（NR）/順丁橡膠（BR）混合膠體的力學性能，從微觀上解釋了NR和BR共混后性能得到改善的原因，并通過界面結合能和徑向分布函數(shù)分析揭示了混合物組分之間的相互作用的本質。結果表明，NR/BR的平均結合能為21.35kJ/mol，除去程序自行產(chǎn)生的能量校正值（－3.64kJ/mol），結合能在數(shù)值上近似等于非鍵能（17.91kJ/mol），表明NR/BR組分界面間的相互作用主要由范德華力構成。

分子動力學；力學性能；界面結合能；徑向分布函數(shù)

0 前言

橡膠自問世以來，隨著工業(yè)生產(chǎn)的不斷發(fā)展和對橡膠需求的日益增加，人們對橡膠的改性研究一直沒有間斷過。目前，對橡膠改性已經(jīng)有了諸多改性手段，其中橡膠共混并用就是常用改性方法之一。長期以來，橡膠并用的配方設計主要依靠實驗，缺乏基礎理論指導，需要進行大量的實驗研究，耗費大量的人力、物力和財力。與實驗研究相比，分子模擬是一種更直接的研究方法，可以通過分子模擬技術先于實驗預測和了解材料的性能，縮短新材料的研制周期，大大降低實驗的成本和危險性。在橡膠工業(yè)中，NR和BR是應用最為廣泛的并用體系，以其制備的輪胎胎面膠、胎側膠以及橡膠輸送帶的覆蓋膠等都具有良好的技術效果，其加工、形態(tài)和性能都已被廣泛研究。通過分析NR/BR的共混原理有利于推進橡膠并用的改性手段。近年來，對共混物力學性能與界面間相互作用進行分子模擬已有相關研究［1－4］，但利用分子模擬研究 NR和BR共混物的微觀機理，至今尚未有這方面的報導。

本文運用MD方法，使用Material Studio 4.0程序和 Compass力場［5－6］分別模擬計算BR和 NR/BR的原子運動軌跡，并比較了其力學性能（彈性系數(shù)、柯西壓、模量和泊松比等）。求解出混合物組分間的界面相互作用（結合能），并以徑向分布函數(shù)揭示了相互作用的本質，從微觀尺度上解釋了NR/BR共混后力學性能得到改善的原因。我們知道，通過分子模擬來解釋NR和BR共混物的一切實驗現(xiàn)象是不可能的，尤其是在如此復雜的界面相互作用的情形下有很多因素（如混合膠填料、添加劑的影響等）都不能考慮在內。所以，我們只希望得到能解釋或預測特定性能的結果，就已經(jīng)達到模擬的目的了。這些研究為本課題組后續(xù)研究“杜仲膠/NR/BR共混膠的結構與性能”和多種橡膠并用的分子模擬建立了一定的基礎數(shù)據(jù)，所得結果可用于定性分析，具有一定的指導意義。

1 模型構建與方法

NR主要成分是順式聚異戊二烯，BR主要成分是聚丁二烯，這兩種聚合物分子式如圖1所示。

圖1 順式聚異戊二烯和順式聚丁二烯分子結構圖Fig.1 Molecular structure of cis－polyisoprene and cis－polybutadiene

運用Materials Studio 4.0軟件包中Visualizer模塊分別構建聚合度為20的NR分子模型和聚合度為15的BR的分子模型，然后依據(jù)在298K和1個標準大氣壓下NR、BR相應的密度（NR相對密度為0.96g/cm3，BR相對密度為0.91g/cm3），利用 Amorphous cell模塊分別搭建NR、BR純物質和NR/BR共混物（質量比為2∶1［7］）三維周期性邊界條件下的無定形高分子模型，如圖2、3所示。

圖2 BR與NR三維無定形結構模型Fig.2 Three－dimensional structural model of amorphous for BR and NR

圖3 NR/BR三維無定形結構模型Fig.3 Three－dimensional structural model of amorphous of NR/BR blends

為了使體系更好的越過勢能面上局部極小值之間的位壘，我們首先利用Smart Minimization方法對所構建的無定形分子模型進行結構優(yōu)化，之后做退火模擬，使體系從200K開始逐漸升溫，間隔為50K，選擇NPT系綜（恒定分子數(shù)、恒定壓力以及恒定溫度的系綜），在體系壓強為0.101MPa時對每個溫度進行50ps的MD模擬，直至溫度達到600K時為止；然后再以50K為間隔，將溫度逐漸降至200K；進行2次循環(huán)，體系在經(jīng)歷了高溫下的多次松弛，空間結構被逐漸優(yōu)化，這一過程使建造的模型中產(chǎn)生的局部不合理結構基本消除。之后進行MD模擬計算分析力學性能等，具體模擬細節(jié)見表1.

表1 MD模擬細節(jié)Tab.1 Details of MD simulation

力學性能通過靜態(tài)模型分析獲得。關于應力－應變關系的廣義虎克定律如式（1）所示。

式中 Cij——6×6彈性系數(shù)矩陣

σi——應力

εj——應變

由于應變能的存在，且是對稱的，因此描述任意材料的應力－應變行為最多只需21個系數(shù)。在原子水平計算中，靜態(tài)模型用維里方程計算內應力張量σ如式（2）所示。

式中 mi——原子質量

vi——原子速度

V0——無形變時體系的體積

施加于體系的應力使體系內粒子相對位置發(fā)生變化。對平行六面體（如模擬所取各邊長分別為a、b、c的周期箱），若以列向量a0、b0、c0表征參比狀態(tài)，以向量a、b、c表征形變狀態(tài)，則應變張量為：

式中 h0——列向量a0，b0，c0組成的矩陣

A——度量張量

h——a，b，c組成的矩陣

因此通過計算拉伸和剪切形變的斜率，即可求得彈性系數(shù)矩陣，然后依據(jù)彈性力學原理，即可求得模量和泊松比等各向同性力學參數(shù)。各向同性材料的拉伸模量（E）、泊松比（v）、體模量（K）和剪切模量（G）可以用拉梅系數(shù)λ和μ 表示為式（5）～式（8）［10］。

2 結果與討論

2.1 力學性能分析

拉伸模量、體模量和剪切模量均可以用來衡量材料的剛性強弱，即在外力作用下抵抗變形的能力。泊松比是材料力學和彈性力學中的名詞，定義為材料受拉伸或壓縮力時，材料會發(fā)生變形，而其橫向變形量與縱向變形量的比值，是一無量綱的物理量。柯西壓（C12－C44）通?？捎脕砗饬矿w系的延展性，即材料發(fā)生形變而不產(chǎn)生裂縫的能力。一般來講，一個具有較好延展性的物質，其（C12－C44）為正值，反之則為負值［11］。根據(jù)靜態(tài)力學分析，對 MD模擬所得體系的平衡運動軌跡分析得到BR和NR/BR的彈性系數(shù)、模量和泊松比等，如表2所示。

表2 BR與NR/BR的彈性系數(shù)和有效各向同性力學性能Tab.2 Elastic constants and effective isotropic moduli for BR and NR/BR

對比表2中BR和NR/BR的彈性系數(shù)，發(fā)現(xiàn)后者的C11、C22和C33的數(shù)值相差不大，且C12、C13和C23組以及C44、C55和C66組內數(shù)值彼此更加接近，表明NR/BR更加趨近于各向同性彈性體［12］。另外，BR的拉伸模量、體模量和剪切模量值均大于混合物NR/BR的模量值，表明BR剛性較強；NR/BR共混后，體系模量下降，表明混合物剛性減弱，柔性增強，較易于發(fā)生彈性形變。BR的柯西壓值大于NR/BR，即BR的延展性強于NR/BR，表明共混后膠料強度增大。這些都表明共混后體系的力學性能有較為顯著的改善。

2.2 徑向分布函數(shù)分析

徑向分布函數(shù)g（r）是反映材料微觀結構的特征物理量，表示在一個分子周圍距離為r的地方出現(xiàn)另一個分子的概率密度相對于隨機分布概率密度的比值。分子內徑向分布函數(shù)可以提供模擬體系有序度的信息，一般出現(xiàn)大于0.3nm的峰表明分子鏈長程有序，屬于結晶體系；而出現(xiàn)小于0.3nm的峰則表明分子鏈短程無序，屬于無定形結構。分子間徑向分布函數(shù)可以揭示非鍵原子間相互作用方式和本質，氫鍵作用范圍為0.26～0.31nm，范德華作用范圍為0.31～0.5nm。我們對整個平衡軌跡進行徑向分布函數(shù)g（r）分析，做出NR，BR與NR/BR分子內與分子間徑向分布函數(shù)，如圖4和圖5。

圖4 BR和NR分子內及分子間徑向分布函數(shù)Fig.4 Plots for pair correlation function computed for atoms of carbon of BR and NR

圖5 NR/BR分子內及分子間徑向分布函數(shù)Fig.5 Plots for pair correlation function Computedfor atoms of carbon of NR/BR blends

從圖4可見，分子內徑向分布函數(shù)中0.3～0.7nm范圍內兩種純物質均沒有出現(xiàn)強的峰，說明BR、NR的分子結構屬于無定形結構；分子間徑向分布函數(shù)則揭示，BR、NR分子間主要作用方式為范德華力作用。二者的徑向分布函數(shù)曲線非常接近，根據(jù)相似相容的原理，可以說明其共混體系的相容性較好。從圖5可見，分子內徑向分布函數(shù)中0.3～0.7nm范圍內兩種純物質均沒有出現(xiàn)強的峰，說明BR/NR混合物屬于無定形聚合物；分子間徑向分布曲線顯示NR/BR混合物內氫鍵作用不大，而形成了較強的范德華力。表明混合物內各組分間的相互作用主要由范德華力構成。

2.3 界面相互作用和結合能

界面相互作用的本質是各組分分子間的相互作用，而各組分間的相互作用直接影響其力學強度等性能。分子間相互作用能或分子間結合能可定量地表征分子間相互作用的強度［13－15］。為細致考察混合物組分間的相互作用，我們以穩(wěn)定構型下各體系的總能量進行分子間相互作用能計算，NR與BR的平均相互作用能（Einter）可表示為：Einter＝Etotal－（EBR＋ENR），其中Etotal為NR/BR的平均總能量，EBR和ENR分別為BR和NR的平均單點能，定義結合能（Ebind）為相互作用能的負值，即Ebind＝－Einter。

經(jīng)MD模擬計算，所得結果如表3。

表3 NR/BR界面的相互作用 kJ/molTab.1 Interface interaction of NR/BR blends kJ/mol

由表3可見，NR/BR 的平均結合能為21.35kJ/mol，若不考慮MS程序自行產(chǎn)生的能量校正值（－3.64kJ/mol），則結合能在數(shù)值上近似等于非鍵能（17.91kJ/mol），即為范德華力（17.71kJ/mol）和靜電力（0.20kJ/mol）之和。其中只有微弱的靜電力，范德華力起主導作用，表明NR/BR組分界面間的相互作用主要由范德華力構成，這與徑向分布函數(shù)分析結果相一致。

3 結論

（1）與BR相比較，BR與NR共混可顯著提高生膠強度，其動態(tài)拉伸疲勞性能極佳，其結論與實驗結果一致。表明可以通過分子模擬來計算橡膠共混后的特定性能；

（2）比較NR與BR純物質的分子內、分子間徑向分布函數(shù)發(fā)現(xiàn)，二者分子鏈短程無序，屬于無定形結構，同時揭示了NR、BR非鍵原子間主要作用方式為范德華力作用。分析NR/BR共混膠徑向分布函數(shù)發(fā)現(xiàn)混合物內各組分間的相互作用主要由范德華力構成；

（3）在NR/BR共混膠中，界面結合能主要由范德華力構成。

［1］ Abou－Rachid H，Lussier L S，Ringuette S.On the Correlation Between Miscibility and Solubility Properties of Energetic Plasticizers/Polymer Blends： Modeling and Simulation Studies［J］.Propellants Explosives Pyrotechnics，2008，33（4）：301－310.

［2］ Karl G Wagner，Martin Maus，Andreas Kornherr.Glass Transition Temperature of a Cationic Polymethacrylate Dependent on the Plasticizer Content－simulation vs Experiment［J］.Chemical Physics Letters，2005，406（1/2/3）：90－94.

［3］ Jaidann M，Abou－Rachid H，Lafleur－Lambert X.Mode－ling and Measurement of Glass Transition Temperatures of Energetic And Inert Systems［J］.Polymer Engineering ＆Science，2008，48（6）：1141－1150.

［4］ Morita H，Tanaka K，Kajiyama T.Study of the Glass Transition Temperature of Polymer Surface by Coarse－Grained Molecular Dynamics Simulation［J］.Macromolecules，2006，39（18）：6233－6237.

［5］ H Sun.COMPASS：an Abinitio Force Field Optimized for Condensed Phase Applications Overview with Details on Alkane and Benzene Compounds［J］.Phys Chem，1998，102（38）：7338－7364.

［6］ S W Bunte，H Sun.Molecular Modeling of Energetic Materials：The Parameterization and Validation of Nitrate Esters in the COMPASS Force Field ［J］.Phys Chem B，2000，104（11）：2477－2489.

［7］ 姚鐘堯，寇利濤.NR/BR并用膠的共混方法與微觀結構和強度之間的相互關系［J］.世界橡膠工業(yè)，2007，34（7/8）：16－24.Yao Zhongyao，Kou Litao.NR/BR Rubber′s Blends Method and the Relationship Between Microstructure and Strength［J］.World Rubber Industry，2007，34（7/8）：16－24.

［8］ Andersen H C.Molecular Dynamics Simulations at Constant Pressure and/or Temperature［J］，Chem Phys，1980，72（4）：2384－2393.

［9］ Parrinello M，Rahman A.Polymorphic Transitions in Sin－gle Crystals：A New Molecular Dynamics Method［J］.J App Phys，1981，52（12）：7182－7190.

［10］ 張慧軍，岳 紅，劉 倩.形狀記憶高分子材料性能評價的分子模擬研究［J］.材料導報，2010，24（9）：116－119.Zhang Huijun，Yue Hong，Liu Qian.Performance Evaluation of Shape Memory Polymer by Molecular Simulation［J］.Materials Review，2010，24（9）：116－119.

［11］ 夏 露.高能材料結構和性能的分子動力學模擬［D］.蘇州：蘇州大學化學化工學院，2008.

［12］ Xu Z L.Elasticity Mechanics［M］.Beijing：Higher Education Press，1990：328－376.

［13］ Radovan T，Alessandro C，Marco F.Computer Simulation of Polypropylene/Organoclay Nanocomposites：Characterization of Atomic Scale Structure and Prediction of Binding Energy［J］.Polymer，2004，45（23）：8075－8083.

［14］ Paolo C，Giulio S，Sabrina P.Many－scale Molecular Simulation for ABS－MMT Nano Composites：Up Grading of Industrial Scraps ［J］.Microporous and Mesoporous Materials，2008，107（1/2）：169－179.

［15］ Subhra M，Golok B，Nando K.Mechanical and Dynamic Mechanical Properties of Miscible Blends of Epoxidized Natural Rubber and Poly（ethylene－co－acrylic acid）［J］.Polymer，1996，37（24）：5387－5394.

Molecular Dynamics Simulation of Mechanical Properties and Surface Interaction for NR/BR Blends

JIANG Hao，YUE Hong，LIU Qian
（School of Science，Northwestern Polytechnical University，Xi’an 710129，China）

Molecular dynamics simulation was applied to investigate the mechanical properties of NR/BR blends，and the reasons of their superior performance were interpreted in microstructure.Finally，by calculating binding energy and pair correlation function，the interface interaction and its essence of the blends were elucidated.The results showed that the average binding energy of NR/BR blends was 21.35kJ/mol，the number of average binding energy was approximately the nonbond energy（17.91kJ/mol）except energy calibration （－3.64kJ/mol）produced by program，which meant that interaction between NR and BR was formed by Vander Waals force.

molecular dynamics；mechanical property；interface binding energy；pair correlation function

TQ332

B

1001－9278（2012）05－0064－05

2011－12－12

聯(lián)系人，russelljiang＠163.com

（本文編輯：劉本剛）