交聯對硅橡膠熱導率影響的分子動力學模擬*

徐文雪梁新剛? 徐向華祝淵

1)(清華大學航天航空學院,熱科學與動力工程教育部重點實驗室,北京100084)

2)(南方科技大學深港微電子學院,深圳518055)

(2020年5月17日收到;2020年6月15日收到修改稿)

1 引 言

電子器件功率密度的增加使得如何有效散熱成為一個關鍵問題[1].同時,電子元器件若不能做到有效散熱,必然會導致溫度升高,降低可靠性.熱失效是電子設備最主要的失效模式,據相關數據統計[2],超過半數的電子設備失效是由不良的冷卻導致的.因此電子元器件的降溫散熱十分重要.

電子器件的散熱過程中,利用界面材料降低界面熱阻是一個重要環節.硅橡膠具有絕緣、耐熱等優勢,是一種常用的熱界面材料基體.硅橡膠基體通常和高熱導率材料粒子結合形成復合材料,其擁有許多產品形態[3].影響硅橡膠復合材料導熱性能的因素主要有硅橡膠基體、導熱填料以及加工工藝等.影響硅橡膠基體導熱性能的主要因素有取向度、交聯度、分子內部結構的致密性和溫度等.導熱填料及加工工藝對硅橡膠熱導率的影響有大量的實驗研究[4,5],但缺少硅橡膠基體結構對硅橡膠導熱影響研究.硅橡膠是由高黏滯性的塑性態硅橡膠生膠經過硫化反應得到的,硫化的目的是通過加成、縮合或者氧化等化學反應過程形成交聯結構.因此本文研究交聯對硅橡膠基體熱導率的影響.

有大量關于無定形材料內部交聯形成空間網狀結構來提升聚合物熱導率的研究.Kanari和Ozawa[6]報道了不同交聯密度下環氧樹脂的熱導率,研究表明熱導率隨著交聯密度的增大線性增大.Kikugawa等[7]用分子動力學方法研究了交聯對聚乙烯和聚苯乙烯熱導率的影響,結果表明熱導率會隨著交聯密度增大而增大,但是兩種材料提升的效果不同;此外他們還發現在低分子量形成的聚乙烯中,交聯能對熱導率有更大的提升.Li等[8]報道了交聯官能團密度對有機硅基材料熱導率的影響,結果表明交聯密度對熱導率有提升的效果,文章認為高交聯度的復合材料能提供更多傳播途徑且自由體積更小,所以熱導率會更高.目前尚未見硅橡膠交聯密度對熱導率影響的報道,因此開展硅橡膠交聯密度對熱導率的研究是有必要的.

2 計算方法

本文采用基于傅里葉導熱定律的非平衡分子動力學方法計算熱導率—雙溫度方法[9](dualthermostat method):該方法需要建立一個溫度分布,統計維持穩定的溫度梯度需要的熱流.系統沿熱流方向被分為若干層,設置兩端作為熱端區域,中間層作為冷端區域,分別控制冷/熱區域的溫度,并統計控溫過程中的能量輸入與輸出,能量的統計是通過粒子的動能變化來計算的,熱流方向為圖1中箭頭的方向.熱導率的計算根據以下公式進行[9]:

其中DEtotal表示時間間隔Dt內總能量變化,DEH表示熱端區域的動能變化量,DEC表示冷端區域的動能變化量,DEerr表示計算誤差導致的能量變化;(2)式中,l是熱導率,下標i取值為H或C,分別表示熱端和冷端,DEi分別按區域求解即可得到熱端/冷端區域的動能變化;A是橫截面積,dT(z)/dz表示沿熱流z方向的溫度梯度.計算中硅橡膠采用COMPASS力場[10],控溫方法采用Nosé-Hoover方法[11],非鍵截斷半徑與庫侖截斷半徑均為12.5?,時間步長為0.05 fs[12].

圖1雙溫度方法示意圖[13]Fig.1.Schematic diagram of dual-thermostat method[13].

3 交聯硅橡膠的建模

甲基乙烯基硅橡膠和交聯劑DBPMH的反應[14]是一種常見的交聯反應,兩者的分子結構如圖2所示,圖中R表示端基,可以是甲基或者乙烯基,DBPMH含有兩個酯鍵,通過斷開酯鍵與硅橡膠高分子鏈形成交聯鍵,從而得到交聯硅橡膠,在圖3中用ROOR表示.反應過程如圖3所示,在交聯劑的作用下,硅橡膠產生兩種交聯反應:一種是硅橡膠兩個甲基中的碳原子形成交聯鍵,如圖3(a)所示;另一種是雙鍵中的碳原子與甲基中的碳原子形成交聯鍵,如圖3(b)所示.

圖2甲基乙烯基硅橡膠和交聯劑DBPMH的分子式Fig.2.Molecular formula of methyl vinyl silicone rubber and crosslinker DBPMH.

圖3甲基乙烯基硅橡膠和DBPMH的反應(a)硅橡膠兩個甲基中的碳原子形成交聯鍵;(b)雙鍵中的碳原子與甲基中的碳原子形成交聯鍵Fig.3.Reaction of methyl vinyl silicone rubber with DBPMH.

COMPASS力場[15]是收費力場,公開的COMPASS力場缺少C=C雙鍵的力場參數,且甲基乙烯基硅橡膠含有極少量的C=C雙鍵,通常低于1%;再者,本文主要研究內容是交聯形成的結構對硅橡膠的影響.綜上所述,采用圖3(a)所示交聯反應過程,忽略乙烯基形成的交聯結構和甲基形成的交聯結構的區別,針對該反應本文做出如下簡化假設:

1)與硅氧鏈高分子相比,交聯劑的分子相對較小,交聯劑活動性遠遠大于硅氧鏈高分子,即交聯反應不會由于交聯劑的缺少而進行不下去;

2)交聯分子非常容易形成自由基,建模過程中直接設置活性原子形成;

3)各個自由基反應活性相同,反應是擴散控制的,即距離決定是否反應.

圖4交聯反應的流程Fig.4.Process of crosslinking reaction.

圖5 設定的反應活性原子R1和R2 Fig.5.Set reactive active atoms R1 and R2.

為了讓交聯過程較為合理,采用交聯-松弛交替運行的策略模擬交聯過程,該流程要求兩個反應活性原子的間隔不能相距太遠,交聯反應的流程如圖4所示.交聯第一步需要構建含有反應活性原子的無定形硅橡膠模型.不同鏈上的反應活性原子分別被命名為R1或R2,為了區別非活性原子,R1和R2被加大顯示,如圖5所示.在Materials Studio軟件的Amorphous模塊中封裝含有反應活性原子的高分子鏈得到計算模型,接下來在Forcite Plus模塊中采用幾何優化和退火對模型進行初步的能量最小化,充分弛豫后得到初始結構進入交聯反應.設置計算的基本參數,統計最大可能交聯鍵數目和已經形成的交聯鍵數目,并判斷需要交聯的鍵數目是否小于最大可能交聯鍵數目,不滿足直接退出,滿足則進入下一步.對模型進行初步優化,在NVT,NPT系綜下松弛.判斷當前搜索半徑是否小于設定的最大搜索半徑,不滿足則直接退出交聯反應,滿足則進入下一步優化平衡,并執行建立新鍵命令.建立新鍵命令是以每一個活性原子為原點,判斷是否有未反應的活性原子在搜索半徑之內,距離小于搜索半徑則建立兩原子的近鄰接觸,判斷兩個原子是否均為未成鍵的原子,如果是的話則施加一個逐漸加大的簡諧約束力,減少近鄰接觸原子間的距離,使原子對的距離縮短至成鍵距離,建立新鍵;否則跳過該近鄰接觸判斷下一個近鄰接觸,若未找到滿足條件的近鄰接觸,則跳出該循環,增大搜索半徑進入下一個循環.建立新鍵以后對模型進行幾何優化并施加微擾,統計模型的熱力學性質,判斷當前交聯密度enow和目標交聯密度etarget的大小:如果enow

圖6建立新鍵的過程示意圖Fig.6.Process diagram of creating new bond.

關于交聯密度的定義目前不是很統一.實驗研究方面由于無法得到內部分子鏈的具體結構,所以只能通過靜態模量、平均分子質量等宏觀參數來表示交聯密度.基于不同的模型或者簡化理論,存在不同的交聯密度表達式:Popovic和王韶暉[16]通過靜態模量方程計算交聯密度g,表達式為g =r/Ma,r是聚合物的密度,Ma是兩交聯鍵之間的數均分子質量;胡冬紅等[17]用Mooney-Riving方程求解高分子鏈上兩個交聯點間的有效分子量Mc,用單位質量高分子鏈中交聯鍵的數目表示交聯密度D,D=1/(2Mc).此外不同的材料、不同摻雜時,顆粒大小、顆粒表面性質等都會影響交聯的程度,所以不同材料間的交聯密度難以比較.數值模擬中交聯密度是用來表征分子交聯的程度,表征相對意義,所以本文交聯密度的定義是已經反應的交聯鍵數目與最大可能交聯鍵數目的比值,與Kikugawa等[7]的定義本質相同.

4 計算與討論

4.1 交聯鍵位置對硅橡膠熱導率的影響

圖7是按照上述交聯流程,起始搜索半徑設置為3.5?,終止搜索半徑設置為8.0?,目標交聯密度是100%的時候得到的若干個交聯結構的示意圖,綠色加粗的原子對即為交聯鍵上的原子.每個計算盒子內保證105個反應活性點R1和105個反應活性點R2,所以每個計算模型中添加了圖5所示的長鏈共30條, 模型的尺寸在32? ×32 ? ×70?左右.由于建模中會用到NPT和NVT系綜進行松弛,所以每個模型的密度和尺寸會略有不同.

圖7 不同交聯密度下的結構示意圖Fig.7.Structure diagram under different crosslinking density.

圖8冷熱端能量隨模擬時間的變化Fig.8.Energy variation at hot and cold ends.

根據上述交聯方法得到了不同交聯密度的硅橡膠結構.按照雙溫度模擬方法,計算得到了不同交聯密度下硅橡膠的熱導率,圖8是交聯密度為0.95%時能量隨時間的變化趨勢圖,可以看出散熱已進入穩定階段,其他交聯密度下結果類似.圖9是交聯密度對熱導率的影響,可以看出,隨著交聯密度的變大,熱導率逐漸升高,交聯密度達到80%的時候,熱導率提高了40%.由于硅橡膠密度變大會導致單位體積內原子數量變多,所以激發的聲子模式也會變多,因此單位體積比熱也會變大,而這對熱導率是有利的.圖10是不同交聯密度下結構的質量密度,可以看出,交聯密度為80%的時候,質量密度只提高了5%.如果把硅橡膠密度對熱導率的影響近似看成是線性的,把密度的倍數乘以交聯密度為1%時的熱導率得到一個由質量密度增大導致的熱導率的大小,將該值與模擬值對比,結果如圖11所示.可以看出,由質量密度增大導致的熱導率增大遠遠不及計算得到的熱導率.所以可以推測,交聯形成的空間網狀結構導致傳熱模式的改變,進而導致熱導率有一個較大的提升.這進一步驗證了Li等[8]認為提高有機硅基復合材料的交聯度能提供更多傳播途徑且自由體積更小,所以對熱導率有提升效果的結論.

圖9交聯密度對熱導率的影響Fig.9.Influence of crosslinking density on thermal conductivity.

圖11由密度導致的熱導率增大與模擬對比Fig.11.Comparison between the increase of thermal conductivity caused by density and simulation results.

4.2 交聯鍵位置對硅橡膠熱導率的影響

上一節中設置的反應活性點的位置是均勻地散布在鏈長上的,如圖12(a)所示.為了研究反應活性點的位置對熱導率是否有影響,本節對只在端部設置反應活性點和只在中間設置反應活性點的兩種情況進行比較.圖12中每幅圖下面的紅色點代表設置了反應活性點的位置,藍色點代表未設置反應活性點.每個計算盒子內保證有100個反應活性點R1和100個反應活性點R2,所以每個計算模型中共添加了50條硅氧鏈.端部設置反應活性點的分子鏈形(圖12(b))在交聯反應中形成的交聯密度分別有:12%,19%,27%,38%,51%,61%和71%;中間位置設置反應活性點的分子鏈形(圖12(c))在交聯反應中形成的交聯密度分別有:1%,9%,20%,31%,41%,52%和63%.反應活性點在不同位置時交聯率對熱導率的影響結果如圖13所示,可以看出,交聯點的位置不論在端部還是中部,熱導率都會隨著交聯密度的增大而增大,端部和中部設置交聯點對熱導率的影響差異較小.由此可推得反應活性鍵的位置對熱導率影響較小,選擇不同位置上的交聯劑對熱導率影響較小.

圖12反應活性原子不同位置分布(a)均勻;(b)端部;(c)中間Fig.12.Positions distribution of reactive atoms.

圖13反應活性點在不同位置時交聯密度對熱導率的影響Fig.13.Effect of crosslinking density on thermal conductivity at different positions of reactive atoms.

4.3 反應活性點間隔對硅橡膠熱導率的影響

前面兩節的建模中,每隔2個硅原子設置一個活性反應點.為了驗證反應活性點的間隔對熱導率是否有影響,本節設計了反應活性點間隔稍大的兩種分子鏈結構,即每隔4個和6個Si原子設置一個反應活性點,如圖14所示.間隔為4個Si原子的分子鏈上有4個活性反應點,間隔為6個Si原子的分子鏈上有3個活性反應點.為了保證反應活性點的數目基本一致,間隔為4的計算盒子內封裝了50條高分子鏈,保證有100(4×25=100)個反應活性點R1和100個反應活性點R2,模型的尺寸在33?×33?×118?左右;間隔為6的計算盒子內封裝了66條高分子鏈, 保證有99(3 ×33 =99)個反應活性點R1和99個反應活性點R2,模型的尺寸在37?×37?×117?左右.間隔為4時,模擬交聯反應形成的交聯密度分別有:16%,31%,49%,65%;間隔為6時,模擬交聯反應形成的交聯密度分別有:3%,20%,39%,55%,72%.

圖14不同反應活性點間隔的分子鏈(a)間隔為4;(b)間隔為6Fig.14.Molecular chains with different reactive atoms intervals.

由于不同間隔下系統的體積不一樣,單位體積內的交聯點數不一樣.間隔為2的時候,活性點的設置最密集,間隔為4和6時,活性點的設置則稀疏一些.為了保證交聯密度的可比性,采用表1所列的方式對交聯密度進行歸一化處理.

表1歸一化交聯密度的計算Table 1.Calculation of normalized crosslinking density.

不同反應活性點間隔對熱導率的影響如圖15所示,間隔為2和4的情況下,熱導率隨著歸一化交聯密度的提高而提高.間隔為6的情況下,當歸一化交聯密度小于20%時,熱導率不隨歸一化交聯密度變化;歸一化交聯密度大于20%時,熱導率增加.對比間隔2,4,6時的熱導率,發現相同歸一化交聯密度下間隔為6的熱導率最大.交聯間隔越大,熱導率越高,為了提高交聯硅橡膠的熱導率,應當適當增加活性點的間距.

圖15不同反應活性點間隔對熱導率的影響Fig.15.The effect of reaction active atoms intervals on thermal conductivity.

不同的交聯間隔下的密度有所不同,密度隨歸一化交聯密度的變化趨勢如圖16所示,對比熱導率隨著歸一化交聯密度的變化趨勢可以看出,密度的變化趨勢大體與熱導率的變化趨勢一致,通常是較高的密度對應較高的熱導率.把質量密度和熱導率放在一張圖中,如圖17所示,可以看出熱導率和質量密度正相關,相關系數為0.838.實際應用中應將硅橡膠材料緊密加載,這是因為緊密加載會使其密度增高進而使其熱導率增加,此外緊密加載時接觸熱阻也會較小.

圖16不同歸一化交聯密度下的質量密度Fig.16.Density under different normalized crosslinking density..

圖17 質量密度和熱導率的相關性Fig.17.Correlation between mass density and thermal conductivity.

4.4 交聯密度對聲子態密度的影響

晶格振動的量子稱為聲子,通常用分子動力學模擬的信息得到聲子的相關信息[18],最常見的是由速度自相關函數的傅里葉轉換得到聲子態密度[19].還有一種方法是基于振動散射理論得到聲子的信息[20],它是利用分子動力學計算中的位移進行傅里葉變換得到倒空間的格林函數[21],再進一步得到力常數矩陣和動力學矩陣,通過計算動力學矩陣的特征值得到聲子態密度[22],該方法能考慮到聲子的非諧波效應.本文聲子態密度的計算方法采用后者, 這一部分工作在LAMMPS軟件的PHONON命令中完成,結合后處理工具PHANA可以得到聲子態密度的信息.

圖18間隔為2時不同交聯密度下的聲子態密度分布(a)Si原子;(b)O原子Fig.18.Distribution of phonon density of states under different crosslinking densities at atoms interval 2:(a)Si;(b)O..

圖19間隔為4時不同交聯密度下的聲子態密度分布(a)Si原子;(b)O原子Fig.19.Distribution of phonon density of states under different crosslinking densities at atoms interval 4(a)Si;(b)O.

圖18所示為間隔為2的情況下交聯密度分別為1%,10%,75%,80%時Si和O原子振動的聲子態密度;圖19所示為間隔為4的情況下交聯密度為16%和65%時的聲子態密度.可以發現,隨著交聯密度增加低頻區聲子數變少,在高頻區聲子數增加.交聯密度增加,使得原子間的相互作用力加強,高頻聲子增加;但是另外一方面,交聯密度的增加,使得從低溫區向高溫區熱量傳遞的路徑更短,因而使得熱導率隨交聯密度的增加而增加.

5 結 論

本文根據甲基乙烯基硅橡膠的交聯反應原型作出三條簡化假設,并根據交聯反應的流程,建立了不同交聯密度下的模型,計算了熱導率.研究發現,熱導率隨交聯密度的提高而提高,交聯密度達到80%的時候,熱導率提高了40%;通過對比分析發現,由結構的質量密度導致的熱導率增大遠遠不及計算得到的熱導率,因此可以認為交聯形成的空間結構網絡對熱導率的促進作用更大.對不同交聯鍵位置的硅橡膠熱導率的比較表明,其熱導率都隨交聯密度變大而變大.當設置的總的最大交聯數相同時,端部設置反應活性點和中間位置設置反應活性點對熱導率基本沒有影響,選擇不同位置上的交聯劑對熱導率影響較小.不同交聯點間隔對熱導率影響的研究表明,間隔較小的情況下,熱導率隨著交聯密度的提高較明顯;間隔較大的情況下,熱導率隨交聯密度先不變后增加.交聯間隔越大,熱導率越高,為了提高交聯硅橡膠的熱導率,應當適當增加活性點的間距. 分析不同間隔下的質量密度與熱導率發現兩者正相關,相關系數為0.838.對聲子態密度進行了計算,發現交聯密度變大,低頻聲子占比變少,高頻聲子占比變多,這是因為交聯使原子間相互作用增強,所以頻率變高.熱導率隨交聯密度變大而變大的原因是交聯鍵的加入縮短了熱傳導的路徑.