三維石墨烯-碳納米管復合結構熱導率的分子動力學模擬

于澤沛，馮妍卉，2，馮黛麗，2，張欣欣，2

（1 北京科技大學能源與環境工程學院，北京100083； 2 北京科技大學冶金工業節能減排北京市重點實驗室，北京100083）

引 言

隨著電子技術的不斷發展，電子元器件體積不斷縮小，而其工作頻率及集成度卻在不斷提高，這就對電子元器件的散熱提出了更高的要求。電子元器件的散熱方式主要是通過導熱介質將熱傳遞給散熱器，并通過進一步的風冷或者水冷的形式對散熱器進行散熱。而如何提高電子元件與散熱器之間的散熱效率，尋找高效的導熱介質，使得電子元件能夠穩定高效長時間的工作，成為當今研宄的熱點[1-10]。

碳材料是一類具有優良傳熱等特性的材料[11]，而碳納米管及石墨烯是其中更具代表性的材料，自問世以來，這兩種材料因具有極其優良的傳熱特性引起了學術界的廣泛關注[12-15]。但是由于這兩種材料的熱導率具有明顯的各向異性，所以其應用得到了很大限制。2008 年Froudakis 等[16]首次提出三維石墨烯-碳納米管復合結構結構(3D graphene-CNT hybrid structure, GCHs)，并探究其在氫氣儲存領域的應用。這種結構不僅彌補了石墨烯及碳納米管在傳熱方向機械特性的弱勢，還兼顧了二者的優良傳熱特性。隨后大量學者針對這種結構開展不同領域的探索和研究[17-23]。

2013 年，Park 等[20]研究了GCHs 的傳熱特性，研究結果表明此三維納米結構具有良好的面內熱導率和法向熱導率。Loh 等[24]研究表明長波長法向聲子模式對GCHs 熱傳導貢獻顯著。此外，Loh 等[25]使用分子動力學方法研究了GCHs 機械應變對界面處聲子動力學的影響。2014 年，Kamaliya 等[26]研究出了一種可以大規模制備PGNs 的方法，并對其光電特性進行了探討。2017 年，Zhang 等[27]針對GCHs 不同形狀的錐型界面開展研究，表明錐型界面的熱導率高于垂直界面。雖然近年來GCHs 的傳熱特性引起很多研究者的關注[28]，但是其微觀的傳熱機理還不明確，碳納米管的添加對整體傳熱的影響，其相關理論還不夠完善，需要進一步的研究。

本文采用非平衡分子動力學方法模擬了GCHs的熱導率，通過聲子態密度和重疊能等機理分析方法探究了GCHs 內部的傳熱機制以及碳納米管的添加對傳熱的影響。本研究將利于GCHs 的應用及高熱導材料的開發，為三維納米材料傳熱機理提供新的啟示，以拓展高導熱三維納米材料的設計與應用。

1 模型與模擬方法

1.1 模型的構建

為了將多層石墨烯與碳納米管(carbon nanotube，CNT)以化學鍵的方式相結合，先在石墨烯片層上挖掉略大于CNT直徑的孔洞，然后將CNT放到對應的位置，通過重新成鍵的形式將二者結合起來[19]。

本文采用Materials Studio 軟件建模。首先，使用軟件自有的石墨模型，優化計算得到30.1 ?×48.8 ?(1 ?=0.1 nm)、層間距24 ? 的多層石墨烯單晶結構，然后在石墨烯每一層上挖去相應孔洞得到多層石墨烯納米網的單晶結構，之后將長度為23.5 ? 的(6,6)扶手椅型CNT 添加移動到孔洞位置，重新成鍵得到三維GCHs 的單晶結構。沿著z 方向建立3 晶胞的超晶胞結構，如圖1所示。

1.2 模擬方法

本文采用非平衡分子動力學的方法模擬計算GCHs 的熱導率。該方法的熱傳遞是通過對體系施加擾動的方式來實現的，通過計算獲得熱流與溫度梯度，然后依據傅里葉導熱定律得到模擬對象的熱導率

式中，J為熱通量；λ為熱導率；?T為溫度梯度。

圖1 GCHs單晶胞及其3晶胞結構Fig.1 Single lattice and 3 lattices structure of GCHs

本文所采用的非平衡分子動力學方法通過外加熱流的形式來對模擬體系施加擾動，即將模擬體系劃分為熱端、冷端和熱傳遞區域三部分，每間隔一定的模擬時間就從冷端取走部分能量，并同時從熱端注入等量的能量，從而使模擬體系產生恒定的熱流，而熱流則由注入和取出的能量來決定。系統通過改變冷端和熱端原子的動量來實現能量轉移。一定的模擬時間后，模擬體系達到平衡狀態，系統通過統計所有原子的能量和溫度獲得模擬體系的溫度梯度，進而通過傅里葉定律求出熱導率。與平衡分子動力學方法相比，非平衡分子動力學方法的優勢在于計算速度快，體系不易溫漂，但是其模擬的難點在于選擇合適的統計層數，只有當統計層內部原子的分熱流均處于平衡狀態時，才能統計所有原子的能量和溫度。

由于本文的研究對象由碳元素組成，因此采用Tersoff勢函數[29]來約束C—C鍵。

1.3 模擬過程

本文所有的模擬計算均在開源軟件LAMMPS[30-32]中進行，為了防止體系過分漂移，將z軸兩端的兩層原子進行固定，同時將兩端除固定原子外的部分原子分別作為冷端和熱端。首先，在NVT(體積恒定，溫度恒定)系綜下對整個模擬體系進行能量最小化，得到優化后的穩定結構。然后在NVE(體積恒定，能量恒定)系綜下，模擬溫度300 K，時間步長0.1 fs，模擬步數1500 萬步，共計1.5 ns 的模擬時間。如果模擬體系最終的平衡溫度接近300 K，則證明體系處于平衡狀態。整個模擬過程采用Nose-Hoover[32-33]恒溫恒壓方法來控制體系的壓強與溫度。為了保證計算的準確性，兩個模擬盒子大小的體系均進行了3次計算并取平均值。

2 結果與分析

2.1 熱導率的計算

系統平衡以后，得到兩個模擬體系平衡狀態下的溫度分布，如圖2所示。利用傅里葉導熱定律的公式來計算熱導和界面熱阻[34]

式中，κ 為體系的熱導率；J 為熱通量；A 為沿熱流方向的橫截面積；?T 為冷端與熱端之間的溫度梯度。

式中，R為界面熱阻；ΔT為界面處溫度階躍。

體系中，沿Z軸，也就是垂直于石墨烯片層方向所施加的熱流為1 eV/ps，晶胞沿z 軸方向的橫截面積即為熱流方向的橫截面積，根據式(2)計算得到GCHs 的熱導率為(1.42±0.05)W/(m·K)。相比較多層石墨烯的層間熱導率，提升了一個量級。界面熱阻為2.15×10-9(m2·K)/W，相比碳管的接觸熱阻6.46×10-8(m2·K)/W 降低了一個量級[34]，證明碳納米管的存在為GCHs 法向方向的傳熱構建了傳熱通道，提高了其法向熱導率。此外，從圖中可以看出，GCHs 溫度的階躍主要集中在石墨稀與碳納米管的界面處，說明法向傳熱的限制在于石墨烯片層和碳納米管界面形變的存在。這是由于聲子由碳管傳遞到石墨烯再傳遞到碳管這一過程中，聲子能量會產生重新分配，振動的不匹配影響了聲子能量的傳遞；界面處出現較強的聲子散射，尤其是低頻聲子散射，而低頻聲子的平均自由程比高頻聲子大，故其對熱導率的貢獻更大，其散射越強，熱導率就越低。

2.2 振動態密度分析

納米材料的熱導率與晶格的振動相關，振動態密度(vibrational density of states，VDOS)可以得到原子的振動狀態及晶體結構中能量的分布，從而可以得到限制能量傳遞的原子或結構。本文通過對比計算GCHs 和多層石墨烯的VDOS 探究了熱導率與結構之間的關系。β 類原子的VDOS(DP,β(ω))是通過其速度自相關函數(velocity autocorrelation function，VACF)的傅里葉變換和權重因子cβ得到。VDOS 的計算公式如下[13,35]

圖2 單晶胞GCHs(a)和3晶胞GCHs(b)的溫度分布圖(300 K)Fig.2 Temperature distribution of single lattice GCHs(a)and three lattice GCHs(b)(300 K)

式中，ω為頻率。

式中，Γβ為β 類原子的速度自相關函數；uiβ(t)為β類原子的速度；〈〉為時間平均。

圖3 為GCHs 的VDOS 計算結果。從圖中可以看出，相比較多層石墨烯，GCHs 的VDOS 的低頻區原子有所降低，而在中高頻區出現多個峰值，證明有更多的聲子被激活參與振動，進而使得其層間熱導率得到大幅度提升。對圖3 數據進行進一步處理，得到不同頻率段對GCHs 熱導率的貢獻，其結果如圖4 所示。從圖中可以看出，雖然中高頻聲子得到不同程度激發而參與傳熱，但是激發的0～30 THz的低頻聲子依然占據主導，對導熱貢獻最大。

2.3 重疊能分析

在VDOS 定性分析的基礎上，為了定量分析GCHs 熱導率與結構之間的關系，本文引入了重疊能的概念，其計算公式如下[31,36]

式中，Eoverlap為重疊區域go(v)的重疊能；h 為Planck 常 數；v 為 頻 率；exp(hv / kBT)- 1 為Boltzmann 分布；T 為熱力學溫度；kB為Boltzmann 常數。重疊能越高表明兩種原子類型的聲子振動越協調，越有利于能量的傳遞。其計算結果如圖5 所示，紅色為碳管和界面區域的重疊能，藍色為石墨烯片層與界面區域的重疊能。

圖3 GCHs及多層石墨烯的歸一化振動態密度(300 K)Fig.3 Normalized VDOS of GCHs and multilayers graphene(300 K)

圖4 GCHs不同頻率段聲子對熱導率的貢獻(300 K)Fig.4 Contribution of phonons at different frequencies to thermal conductivity of GCHs(300 K)

圖5 GCHs不同區域不同頻率段的重疊能Fig.5 Overlapping energy of different frequency in different regions of GCHs

可以看出，在熱導率貢獻占主導的低頻聲子中，石墨烯片層和界面區域的重疊能大于碳管與界面區域的重疊能。這與上述計算中溫度分布在界面處出現階躍，熱阻陡增的結論相一致，并進一步明確，是碳管-界面區域的聲子能量局域化最為嚴重，聲子振動的協調性最差，影響了能量的傳遞，是限制GCHs熱導率的關鍵因素[37]。

3 結 論

本文基于非平衡分子動力學模擬了三維石墨烯-碳納米管復合結構的法向熱導率。

(1)其法向熱導率為1.42 W/(m·K)，相比多層石墨烯的法向熱導率提高了一個量級，其界面熱阻為2.15×10-9(m2·K)/W，相比碳納米管的接觸熱阻降低了一個量級，證明碳納米管的存在為石墨烯層間構建了傳熱通道，提高了其層間熱導率。

(2)通過溫度階躍位置對應結構，發現界面是熱阻產生的主要原因，界面形變的存在卻影響了熱導率的進一步提高。

(3)通過計算多層石墨烯和GCHs 的振動態密度，發現碳納米管的添加使得更多頻率范圍的聲子得到了激發，進而提高了法向熱導率，但是對熱傳導起主導作用的依然是低頻聲子。

(4)通過計算界面結構的重疊能，探究了三維石墨烯-碳納米管復合結構結構能量的傳遞及聲子的局域化情況。在低頻聲子區域，碳納米管與界面原子的振動不匹配度，是影響熱傳導的最主要的限制因素，進一步驗證了同一元素的三維結構中，界面處的形變是限制熱量傳遞的主要因素，進而為后續高熱導材料的改進及開發提供了一定的方向性指導：同種元素的三維結構中，界面突變越少，界面越平滑，結構原子類型越少，原子間振動的協調性越好，聲子散射越少，能量局域化的程度越低，熱導率就越高。