聚丙烯／聚苯乙烯共混體系相態(tài)結構的耗散粒子動力學模擬

陸未谷，張玉梅，陳彥模

（東華大學a.材料科學與工程學院；b.纖維材料改性國家重點實驗室，上海 201620）

耗散粒子動力學（dissipative particle dynamics，DPD）模擬是對于具有動態(tài)和流變性質的簡單及復雜的流體的一種計算模擬方法.在DPD模擬中構建一系列粒子（即所謂粗粒化），其相互間的碰撞規(guī)則和運動規(guī)律遵循牛頓運動定律，這些粒子不是代表單一原子，而是代表整個分子或流體的區(qū)域，且計算過程與原子的細節(jié)無關.對于小分子而言，粒子構建較為容易，一般為整個分子；對于高聚物，情況變得復雜，DPD模擬采用“珠－簧”模型[1].近年來，DPD模擬已被應用于物理和化學研究的眾多領域.

在聚丙烯（PP）纖維實際生產與研究中，其致密的超分子結構成為一個棘手的問題，如纖維可染性和細旦化的集合共混成型加工.因此研究共混改性PP，特別是共混體系相態(tài)結構，具有重要的意義.

本文采用DPD模擬對PP／聚苯乙烯（PS）二元共混體系進行系統(tǒng)研究，將溫度、PS體積分數、PS相對分子質量以及剪切速率作為考察因素進行分析，為共混體系的相態(tài)結構控制提供較為全面的依據和參考.

1 DPD模擬基本原理

DPD模擬中，粒子間有3種成對出現(xiàn)的作用力：保守力（conservative force）、耗散力（dissipative force）、隨機力（random force）[2－3].耗散力是指當前模擬尺度下體系消耗掉進入更微觀尺度下的能量，隨機力則指從更微觀尺度進入當前尺度下的能量.因此，用牛頓運動定律計算粒子狀態(tài)，必須通過引入耗散理論，設置耗散力參數和隨機力經驗值以保證模擬體系的動量守恒，即溫度恒定.文獻[4]將DPD模擬與Flory－Huggins平均場理論相結合，建立高聚物體系的原子和分子信息與DPD粒子之間的映射關系.

就相對分子質量確定的高聚物分子而言，分子尺寸會隨著構象的變化而變化[5－6].文獻[7]認為，應采用高分子鏈的持久長度所涉及的范圍來構建高聚物分子的粒子.因為長度長于持久長度的高聚物，其性質可以用統(tǒng)計學方式描述，就像三維隨機游動，可以不考慮構象的影響.由于DPD模擬并不涉及分子細節(jié)，因此在這種尺度上構建粒子是合理的.

在DPD模擬中，為了簡化單位，通常采用模擬構建的單個粒子的質量和截斷半徑rc（粒子間相互作用范圍）來作為單位量綱，本文也采取這一做法.

2 DPD模擬參數計算

2.1 特征比

對真實的高聚物分子PP和PS進行粗?；?，關鍵步驟是建立真實分子結構和介觀參數之間的聯(lián)系.本文中用高分子鏈持久長度中所包含的重復單元數即特征比來建立這一聯(lián)系，式（1）對真實高聚物分子鏈進行粗?；玫酱至；Ｐ蚚8].

其中：NDPD為高分子單鏈上DPD粒子的個數；Mp為高聚物的相對分子質量；Mm為單元結構的相對分子質量；Cn為特征比.

考慮到計算機的運算能力，選定相對分子質量比較適中的PP和PS樣品，其相對分子質量分別為9 000和10 000，特征比分別為6.7和9.9[9].為了運算方便，選擇的特征比為10.按照式（1），具體的粗粒化方法[10]以PP為例：相對分子質量為9 000，重復單元的相對分子質量為42，其數均聚合度N≈214，因此NDPD＝214／10≈21.通過此方法得到每條PP分子鏈上含有的DPD粗?；Ｗ訑的浚笤儆脧椈蓪⑾噜廌PD粒子鏈接起來，得到線性PP分子的DPD粗粒化模型.采用同樣的方法對PS分子進行粗?；?

2.2 保守力排斥強度選擇

在DPD模擬中，最為重要的參數是DPD粒子之間的保守力排斥強度（同類型粒子間的保守力排斥強度為aii、不同粒子間的保守力排斥強度為aij，其中aii與模擬體系的數密度和Flory-Huggins相互作用參數有關）.為了有效關聯(lián)DPD粒子體系的壓縮系數與高分子熔體的壓縮系數，且能合理地描述流體的密度變化，aii的選取應遵循式（2）[2]：

其中：kBT反映出溫度信息（kB為波爾茲曼常數，T為絕對溫度）；ρ為數密度.在模擬中，選取kBT＝1，ρ＝3，計算可得同類型粒子間的保守力排斥強度aii＝25.不同粒子間的保守力排斥強度aij可使用Flory－Huggins相互作用參數χ計算.

對于PP／PS的二元共混體系，通過 Material Studio軟件的Blend模塊進行計算，擬合直線如圖1所示.

圖1 PP／PS共混體系的Flory-Huggins參數χ與 溫度倒數1／T之間的關系Fig.1 Interaction parameterχas a function of inverse temperature for a series of PP／PS blending system

由圖1可得，擬合直線為

依據式（3）可求得450K時，χ為0.541.但需要注意，通過式（3）求得的是丙烯單體和苯乙烯單體之間的Flory-Huggins相互作用參數，為了得到粗?；蟮腄PD粒子與粒子之間的相互作用參數，還需要在單體間相互作用參數的基礎上乘以特征比Cn.如前所述，選取的特征比為10，因此，DPD粒子與粒子之間的Flory-Huggins相互作用參數應為5.41.根據式（4）可以求得不同種類DPD粒子之間的保守力排斥強度[4]：

不同溫度下的PP和PS的保守力排斥強度如表1所示.

表1 PP和PS分子參數Table 1 PP and PS molecular parameters

2.3 其他模擬參數的選擇

本 文 采 用 正 則 系 綜 （canonical ensemble，NVT，表示具有確定的粒子數N、體積V、溫度T）作為宏觀約束條件，并應用周期性邊界條件進行模擬.計算盒子的大小為20×20×20r3c，體系數密度為3，因此體系包含DPD粒子24 000個.為了便于模擬，截斷半徑rc、質量m以及溫度kBT的單位都設置為l，彈簧常數C為4.0，耗散因子、噪聲強度均選擇DPD方法的經典值（4.5和3）.由于所采用的積分方法為velocity-verlet算法，為了使在該積分方法中能較好地控制溫度，所選用的時間步長為0.05.每一次模擬均從均勻共混狀態(tài)開始，運行至少1×104步直至相態(tài)結構達到平衡.

3 PP／PS相態(tài)結構的形成

圖2所示為PS相在共混體系中的聚并過程，通過與文獻[11]中的掃描電鏡圖樣（如圖3所示）的對比，發(fā)現(xiàn)隨著動力學運行時間的增加分散相顆粒的體積都會因聚并而逐漸變大，但是，由于模擬中看到的分散相情況與掃描電鏡照片[11]所反映出的現(xiàn)實中共混體系的結構在空間尺度上的差異，模擬的有效性和真實性需要驗證.

圖2 PS體積分數為80%，kBT＝1，每1 000步模擬截取的圖形Fig.2 PS volume fraction is 80%，kBT＝1，the figure of each 1 000steps in simulation

圖3 PP／PS（20／80）不同聚并時刻的電鏡圖片[11]Fig.3 PP／PS（20／80）electron microscopic images of different coalescence moments

在DPD模擬中看到的分散相情況與掃描電鏡照片有共同點并非偶然.由于計算機性能的局限性，模擬的體系不能很大，而現(xiàn)實中的共混體系無論是平衡時間或是微相尺寸都要比模擬的體系大好幾個數量級.另外計算機模擬時的初始狀態(tài)兩相共混得非常均勻，分子鏈之間的距離很短，隨后開始的聚并過程就很迅速；而現(xiàn)實中的共混不可能完全均勻，PP和PS在相對較大的尺度中依然是兩相分離的；PS顆粒之間的距離相對較遠，熱擾動的效果不明顯，碰到一起就需要很長的時間，如果沒有剪切誘導，完全靜態(tài)地進行聚并則時間更長.

分子模擬過程提供了一個能在小尺度下直觀考察相行為的機會.因為，無論是宏觀聚并還是微觀聚并，本質上是一樣的，都是分子在一定的初始熱擾動的基礎上開始運動而產生的結果.

4 PP／PS相態(tài)結構的影響因素分析

4.1 模擬參數選擇

靜態(tài)聚并過程模擬（剪切速率為0）：kBT＝1.0，1.1，1.2，1.3；PP體積分數＝0.1，0.2，0.3，0.5；PS分子鏈粒子數NDPD＝5，10，20，50；其他模擬參數參照2.3節(jié)選取.

剪切聚并過程模擬（剪切速率為x軸方向上）：kBT＝1.2；剪切速率r′＝0.2，0.4，0.6；PS體積分數＝0.3；PS分子鏈粒子數NDPD＝10；其他模擬參數參照2.3節(jié)選取.

4.2 結果討論

4.2.1 溫度對共混體系相形態(tài)的影響

不同溫度條件下的模擬結果如圖4所示，圖中的顆粒為PS分散相.由圖4可以看出，溫度越高，則PS相的分散性越好，顆粒尺寸越小.因為溫度的升高加快了高分子及其鏈段運動速度，同時降低了兩相組分的黏度，降低了體系內界面的表面張力.

圖4 不同溫度下PS分散相的最終結構Fig.4 The PS final structure at different temperatures

4.2.2 剪切速率對共混體系相行為的影響

圖5所示為不同剪切速率下共混體系動力學模擬的情況.從圖5中可以看出，隨著剪切速率的增大，PS分散相的顆粒明顯變小，同時沿剪切方向，分散相被拉長，說明剪切破壞了分散相的聚并過程.

圖5 不同剪切速率下PS分散相的最終結構Fig.5 The PS final structure at different shear rates

另外，分析分散相的形成過程發(fā)現(xiàn)，隨著剪切速率的增加，穩(wěn)定的分散相顆粒形成的時間提前了.圖6列出了圖5所示模擬過程中前3幀過程的PS相態(tài)結構.從圖6可以看出，初始結構相同，剪切速率為0.6的體系在3幀后基本形成最終結構，而沒有剪切的體系，則很難看出最終結構情況，說明剪切誘導聚并的發(fā)生.

圖6 不同剪切速率下PS分散相的演變過程Fig.6 The PS evolution process of the dispersed phase at different shear rates

在外界力場作用下，分相顆粒的分散過程和聚并過程共同決定了共混體系的相態(tài)結構.條件不同則這兩個過程對最后形成的相態(tài)結構的影響程度也不同，對此可用Elemendorp相圖（圖7所示）[12]來研究.從圖7可知，相圖包含了顆粒破裂曲線和顆粒聚并曲線，在全圖中僅在區(qū)域Ⅱ內，相態(tài)結構的演變僅由聚并過程控制.

在前述模擬過程中，不同剪切速率下顆粒的平均尺寸和粒間距離都變大，這一現(xiàn)象與靜態(tài)的聚并過程相似；不同的是，最終穩(wěn)定的分散相顆粒形成的速率明顯大于靜態(tài)聚并的情況，說明剪切作用的存在能夠有效地加快顆粒間的聚并.

4.2.3 PS體積分數對共混體系相行為的影響

圖7 Elemendorp相態(tài)圖[12]Fig.7 Elemendorp phase schemes

通過觀察PS相的聚并過程，分析PS體積分數對共混體系相態(tài)結構的影響，不同PS體積分數下共混體系的相態(tài)結構如圖8所示.試驗開始默認PS均勻分布在PP相中，隨著時間的推移，粒子在運動過程中發(fā)生黏結和聚并，達到平衡時，形成相態(tài)結構.隨著PS體積分數的增加，PS相的尺寸變大，且分布變寬，聚并的傾向明顯，PS分散相的數目減少.由于體積分數增加，PS相之間接觸的機會增加，由于PS與PP之間的相容性不好，所以都朝著各自組分的方向前進，最終均形成連續(xù)相.

圖8 不同PS體積分數下（kBT＝1.3）共混體系進行模擬后最終的相態(tài)結構Fig.8 The PS final phase morphologies structure at different PS volume fraction（kBT＝1.3）

4.2.4 PS相對分子質量對共混體系相行為的影響

相對分子質量對共混體系的影響較為特殊.由于研究體系發(fā)生變化，特別是PS的分子結構發(fā)生變化，分子間的Flory-Huggins相互作用參數也發(fā)生了變化，這會導致保守力排斥強度aij發(fā)生變化，模擬的結果可能不是很精確，因此對圖9所示的不同相對分子質量PS的共混體系最終相態(tài)結構的圖像變化趨勢描述可能存在著誤差.

圖9 不同相對分子質量PS的（kBT＝1.2）共混體系進行模擬后最終的相態(tài)結構Fig.9 The PS final phase morphologies structure at different PS molecular weight（kBT＝1.2）

從圖9可以看到，相對分子質量較小的PS相形成的分散相尺寸較均勻，相對分子質量較大的PS相形成的分散相尺寸差異較大.這是由于在分散相聚并過程中，相對分子質量越大的PS分子整體上受到的力就越多，而且各個方向都有，使得其所受合力較小，所以運動速度較慢，因此相對分子質量大的PS相發(fā)生聚并的概率就小，但是由于越長的高分子鏈之間的纏結點就越多，相互結合越緊密，分子間的作用力越大，表現(xiàn)出黏度很大，這導致其與PP相的黏度差異變大，因此PS相一旦發(fā)生聚并就很難分開，因而相對分子質量大的PS相會形成尺寸差異較大的分散相顆粒.

5 結 語

本文依據DPD模擬的原理，通過Material Studio 5.0軟件的Blend模塊計算了不同溫度下Flory－Huggins相互作用參數χ及DPD粒子之間的排斥參數.通過對PP／PS二元共混體系的相態(tài)模擬，研究了溫度、剪切速率、PS體積分數以及PS相對分子質量對PP／PS共混體系的影響，得到下述結論.

（1）共混體系溫度越高，分散相PS相的分散性越好、顆粒尺寸越小.

（2）剪切速率增大，分散相的顆粒明顯變小，分散性變好，同時通過對聚并過程的觀察，了解到剪切又能誘導聚并發(fā)生.

（3）隨著PS體積分數的增加，分散相尺寸逐漸增加，最終成為連續(xù)相.隨著PS相對分子質量的增加，PS的黏度增加，造成PP與PS的黏度比變大，形成尺寸差異較大的分散相顆粒.

參 考 文 獻

[1]SCHLIJPER A G，HOOGERBRUGGE P J，MANKE C W.Computer-simulation of dilute polymer-solutions with the dissipative particle dynamics method[J].Journal of Rheology，1995，39（3）：567-579.

[2]ESPANOL P，WARREN P.Statistical-mechanics of dissipative particle dynamics [J]. Europhysics Letters，1995，30（4）：191-196.

[3]朱孑系.液體自擴散系數與溫度的關系[J].化工學報，1959（2）：172-175.

[4]GROOT R D，WARREN P B.Dissipative particle dynamics：Bridging the gap between atomistic and mesoscopic simulation[J].Journal of Chemical Physics，1997，107（11）：4423-4435.

[5]MA P W，LIU W C，WOO C H，et al.Large-scale molecular dynamics simulation of magnetic properties of amorphous iron under pressure[J].Journal of Applied Physics，2007，101（7）：073908-9.

[6]RUBINSTEIN M.Polymer physics of the perciliary layer[J].Pediatric Pulmonology，2007，30（sl）：112-113.

[7]GROOT R D，MADDEN T J.Dynamic simulation of diblock copolymer microphase separation [J].Journal of Chemical Physics，1998，108（20）：8713-8724.

[8]GAI J G，LI H L，SCHRAUWEN C，et al.Dissipative particle dynamics study on the phase morphologies of the ultrahigh molecular weight polyethylene／polypropylene／poly（ethylene glycol）blends[J].Polymer，2009，50（1）：336-346.

[9]BRANDUP J，IMMERGUT E H，GRUIKE E A.Polymer handbook[M].4th ed.New York：John Wiley and Sons，1999.

[10]苑世領，劉成卜，徐桂英，等.聚合物PVP與表面活性劑AOT相互作用的介觀模擬[J].高等學?；瘜W學報，2003，24（6）：1048-1051.

[11]李云巖.聚丙烯／聚苯乙烯共混物相結構形成演變及其動力學研究[D].天津：天津大學材料科學與工程學院，2007.

[12]ELMENDORP J J，VANDERVEGT A K.A study on polymer blending microrheology 4：The influence of coalescence on blend morphology origination [J].Polymer Engineering and Science，1986，26（19）：1332-1338.