天然橡膠的拉伸行為模擬

(青島科技大學機電學院 山東 青島 266000)

引言

橡膠材料擁有很好的彈性且容易變形，自19世紀中葉就成為了一種十分重要的工程材料，廣泛用于彈簧、減震襯墊、輪胎、密封件、聯軸器和載重結構的座駕。對于橡膠材料的研究已有很長的歷史，在過去50多年里人們試圖對橡膠材料做真實可靠的描述，對其進行了大量實驗，根據橡膠材料在不同環境條件、應變歷程、應變率和加載速率等條件下的力學性能變化規律，建立了許多基于不同理論的本構模型[1-4]。但是橡膠材料的宏觀性質是由單個高分子鏈的性質決定的(如，單鏈的熵和焓彈性或共價鍵的力學穩定性等)。利用傳統的研究方法，我們只能得到橡膠材料在宏觀尺度的力學性質信息。這些信息中既包含著不同高分子鏈之間的分子間相互作用又包含著單個聚合物鏈對材料宏觀性質的貢獻。而要想從這復雜的信息中提取出某一個因素(單個高分子鏈的彈性或分子間相互作用)對材料力學性質的影響幾乎是不可能的。但是隨著計算機的發展和分子動力學技術的發展，使對橡膠材料力學性能微觀機理的研究成為可能。近年來已有越來越多的科研小組開始應用分子模擬方法對橡膠等高分子材料進行了研究并取得了很大的進展[5-6]。

分子動力學方法(MD)是基于牛頓力學的一種模擬計算方法，有很高的準確性和有效性,已被廣泛應用于化工、生物、醫學、材料和能源等各個領域。對于材料研究，它省掉了很多繁瑣的實驗過程，可以較低成本對材料進行篩選、設計、配比等。本研究用MD方法中的cvff力場對天然橡膠的單向拉伸過程進行模擬，以研究實驗中難以控制的分子鏈聚合度對NR力學性能的影響，分析在拉伸過程中的各項能量變化，計算分子鏈的均方位移、均方回轉半徑，揭示聚合度對天然橡膠力學性能影響的微觀機理。

一、模擬方法

(1) 本文采用的是cvff(consistent-valence forcefield)力場，它是Discover程序提供的原始力場，是一個廣義價鍵力場(Dauber-Osguthorpe等，1988)。為氨基酸，水和各種其他官能團提供參數。cvff適用于有機小分子(酰胺，羧酸等)晶體和氣相結構，可以處理多肽，蛋白質和各種有機體系。作為Discover的默認力場，已被廣泛使用了很多年。它主要用于研究結構和結合能，可以很好地預測振動頻率和構象能。

(2) 模型構建：NR的主要成分是順式-1，4-聚異戊二烯，所以我們用順式-1，4-異戊二烯為基本重復單元來構建NR分子模型。構建NR分子模型的軟件我們選用Material Studio 7.0軟件。為了保證模擬過程的可靠性，首先要確定能夠真實反映聚合物的聚合度。通過改變聚合物的聚合度，計算聚合物不同聚合度時的溶度參數(δ)，當溶度參數不再隨著聚合度的改變而變化時，就說明該聚合度能保證計算結果的準確性。劉倩[7]等人已經通過模擬將得到的結果與文獻值進行對比得到了能夠真實反映NR的最小聚合度在20左右。

采用Material Studio 7.0軟件中的Visualizer模塊分別構建聚合度為20、50、100和200的分子模型，依據在298K、1atm下NR的密度0.92g/cm3，搭建其在三維周期性邊界條件下的尺寸為80?×80?×80?無定形分子模型。然后進行退火處理，選用NPT系綜(指控制體系的粒子數N、壓強P、溫度T以及體系的總動量不變)保持體系壓強為0.101MPa，使體系從200K溫度起以50K間隔逐步直至達到500K，每一個溫度保持200ps，然后在意相同的溫度間隔將溫度再逐步降回到200K，如此反復3次，體系經過了多次高溫的松弛膨脹，基本可以消除建立的模型中的局部不合理結構，從而為后續的模擬提供比較合理的幾何結構。最后利用Discover Minimizer模塊的Smart Miniminzer方法進行能量最小化處理。

(3) 將構建好的模型導入LAMMPS軟件進行模擬。首先在298K、101.325kPa條件下利用恒溫恒壓(NPT)系綜對體系進行平衡，使其結構達到穩定。然后分別對建立的四種分子模型進行在相同拉伸速度條件下的單向拉伸模擬。

二、結果與討論

(一)聚合度對拉伸強度的影響

模型沿盒子的x軸方向進行單軸拉伸，步長ΔL=Lx,i-Lx,i-1=0.008nm,Lx,i和Lx,i-1分別是在迭代次數為i-1和i時盒子沿X方向的長度。每10步獲取一次體系的性質信息。圖1所示為幾種模型的應力-應變曲線。

圖1 不同聚合度的天然橡膠模型的應力-應變曲線

一般來說，彈性體在單軸彈性變形范圍內得到的彈性應力-應變響應可以表示為

σ(ε)=E(ε)ε

(2-1)

式中，拉伸應變ε=(Lx,i-Lx,0)/Lx,0，σ(ε)為拉伸應力，E(ε)為彈性模量[8]。

圖1中我們可以看到，高聚合度體系的應力-應變行為表現出與在無定型聚合物的實驗測試中觀察到的具有相同趨勢，即彈性狀態，然后是軟化和應變硬化狀態。隨著聚合度的增加，體系的彈性模量和大應變時的應力上升有較為明顯的增加，聚合度為200的體系的彈性模量大概是聚合度為20時的兩倍，同時高聚合度體系的強度極限要大于低聚合度體系，表明了天然橡膠的力學性能對聚合度有較強的依賴性。聚合度為20的體系彈性階段和屈服階段較為明顯，但無明顯的強化階段，推測是由于分子鏈長度過短使分子鏈剛性過強，分子鏈不易發生纏結導致拉伸過程后期結構松散、自由體積過大而使分子鏈間的相互作用過小。

(二)分子中各項能量對拉伸過程的影響

圖2 拉伸過程中個能量的變化曲線：鍵能E_bond，二面角勢能E-dihed，非鍵合作用E_non-bond，鍵角能E-angle

經過模擬我們發現各聚合度模型的能量變化基本一致，為了研究在拉伸過程中各種能量的變化規律，我們取聚合度為50的體系進行分析，如圖2所示。彈性狀態下，能量變化曲線表明，鏈間非鍵合(范德瓦爾斯)相互作用發揮了重要作用。在體系發生彈性變化過程中，體系中的鏈間非鍵合(范德瓦爾斯)相互作用能在持續增長，并在屈服階段應力頂峰處附近達到最大，但隨著拉伸過程的延續，體系中的自由體積逐漸變大，使分子鏈間的粒子距離逐漸增加，分子鏈間的鏈間非鍵合相互作用也就相應減弱。在整個拉伸過程中鍵合能與二面角勢能發生了較大變化，到模擬過程結束時鍵合能增長了將近5000kcal/mol的能量，而二面角勢能減少了大約2000kcal/mol的能量。鍵合能的增長主要是因為分子鏈間的相互纏結導致在拉伸過程中分子鏈中的鍵長被拉長而導致，而二面角能量的減少與體系中大規模構象轉變為較低能量的旁氏構象有關。鍵角能在整個拉伸過程中并沒有發生較為明顯的變化。

(三)拉伸過程中分子結構的變化

模型在拉伸過程中分子鏈的活動狀況可以用均方位移(MSD)和均方回轉半徑(Rg)來表征[9]。前者是指在一段時間內系統中的粒子相對于初始時刻的位移的平方，其定義式為式(2-2)。后者是衡量分子鏈在空間伸展程度的一種尺度，能表征分子鏈的柔順性和鏈的活動性，其定義是為式(2-3)。

MSD=<|r(t)-r(0)|2>

(2-2)

式中，r(t)與r(0)分別表示在t時刻與初始時刻粒子的坐標。

(2-3)

式中，mi為分子鏈中第i個鏈單元的質量，ri為分子鏈質心到第i個鏈單元的質量。

圖3 (a)均方位移；(b)均方回轉半徑

如圖3(a)所示。圖中在彈性狀態范圍內，幾種模型的均方位移在彈性階段都出現了平臺區域。由圖2能量變化曲線可知，由于在彈性階段分子之間的非鍵合相互作用力較大，使橡膠結構較為穩定而不易發生變形。圖3(b)中可以清晰的看到隨著鏈長的增加，分子鏈的均方回轉半徑有很明顯的變化。聚合度為20的模型在拉伸過程中均方回轉半徑基本未發生變化，說明其分子鏈剛性較大，不易產生變形，驗證了我們之前的推測。在模擬過程中我們觀測到隨著聚合度的增加，拉伸過程中分子鏈的具有更好的柔順性和活動性，分子鏈之間更易產生纏結和牽拉而使分子結構中的自由體積不會產生過大的變化，與圖中MSD和Rg2的變化趨勢相一致。

三、結語

通過本研究得到一下結論：

(1) 在一定范圍內隨著聚合度的增加，NR的結構強度更高，具有更加優異的力學性能；

(2) 由能量變化分析可知，拉伸過程中得分子鏈非鍵合作用能、鍵能和二面角勢能變化較大，對NR的力學性能起決定性作用。其中，非鍵合作用能對NR的彈性狀態起重要作用。鍵角能在整個拉伸過程中并未發生較大變化，對NR的力學性能沒有太大影響；

(3) 通過Rg2和MSD分析表明，分子鏈的聚合度的增加使其柔順性和活動性增加，分子鏈間更加容易產生纏結和牽拉，這是使NR具有更加優異的力學性能的重要原因。

【參考文獻】

[1]朱艷峰, 劉鋒, 黃小清,等. 橡膠材料的本構模型[J]. 橡膠工業, 2006, 53(2):119-125.

[2]陳家照, 黃閩翔, 王學仁,等. 幾種典型的橡膠材料本構模型及其適用性[J]. 材料導報, 2015(s1):118-120.

[3]Amnuaypornsri S, Toki S, Hsiao B S, et al. The effects of endlinking network and entanglement to stress-strain relation and strain-induced crystallization of un-vulcanized and vulcanized natural rubber[J]. Polymer, 2012, 53(15):3325-3330.

[4]Sliozberg Y R, Hoy R S, Mrozek R A, et al. Role of entanglements and bond scission in high strain-rate deformation of polymer gels[J]. Polymer, 2014, 55(10):2543-2551.

[5]Xie F, Lu Z, Yang Z, et al. Mechanical behaviors and molecular deformation mechanisms of polymers under high speed shock compression: A molecular dynamics simulation study[J]. Polymer, 2016, 98:294-304.

[6]Hooper J B, Bedrov D, Smith G D, et al. A molecular dynamics simulation study of the pressure-volume-temperature behavior of polymers under high pressure[J]. Journal of Chemical Physics, 2009, 130(14):144903.

[7]劉倩, 岳紅, 江浩,等. 杜仲膠/天然橡膠共混物的分子動力學模擬和耗散粒子動力學模擬[J]. 材料導報, 2012, 26(6):141-145.

[8]Simes R, Cunha A M, Brostow W. Molecular deformation mechanisms and mechanical properties of polymers simulated by molecular dynamics[J]. e-Polymers, 2004, 4(1):761-783.

[9]Gooch J W. Radius of Gyration[J]. 2011.