邵 鈳，劉遠(yuǎn)超*，鐘建斌，鄒 玉，黃志東

(1.北京石油化工學(xué)院機(jī)械工程學(xué)院，北京 102617； 2.遼寧科技學(xué)院機(jī)械工程學(xué)院， 遼寧 本溪 117004)

隨著高新技術(shù)的迅速發(fā)展，電子器件的特征尺寸已縮小到微納米數(shù)量級(jí)。同時(shí)隨著5G時(shí)代的到來(lái)，智能設(shè)備的功能將更加復(fù)雜，其核心部件——芯片的尺寸也將急劇減小，導(dǎo)致元器件的熱流密度迅速增加[1]。碳納米管(Carbon Nanotube，CNT)自從1991年被日本學(xué)者Iijima[2]發(fā)現(xiàn)以來(lái)，因其具有超常的導(dǎo)熱性能和極好的應(yīng)用前景，在微電子領(lǐng)域具有巨大的應(yīng)用潛力，碳納米管的導(dǎo)熱性能研究也成為了研究熱點(diǎn)。

碳納米管導(dǎo)熱性能研究方法有實(shí)驗(yàn)測(cè)量和分子動(dòng)力學(xué)模擬2種方法[3-10]。其中，分子動(dòng)力學(xué)模擬方法是一種成本低且能在計(jì)算機(jī)上實(shí)現(xiàn)模擬的“實(shí)驗(yàn)”方法，已成為另一種研究碳納米管的熱物性的有效手段。

碳納米管熱導(dǎo)率的影響因素主要有溫度、管徑、管長(zhǎng)和各種缺陷(空位、摻雜和Stone-Wales缺陷)等。迄今為止，人們對(duì)于碳納米管熱導(dǎo)率隨溫度及管徑變化規(guī)律的研究相對(duì)較多，而對(duì)于管長(zhǎng)及各種缺陷的影響研究并不多。關(guān)于管長(zhǎng)影響規(guī)律的代表性研究有：Maruyama[11]模擬發(fā)現(xiàn)CNT熱導(dǎo)率隨長(zhǎng)度呈指數(shù)變化的關(guān)系；王照亮等[12]采用四焊盤(pán)-3ω實(shí)驗(yàn)法測(cè)量發(fā)現(xiàn)熱導(dǎo)率隨長(zhǎng)度增加表現(xiàn)出微尺度效應(yīng)；曹炳陽(yáng)等[13]通過(guò)模擬發(fā)現(xiàn)碳納米管熱導(dǎo)率與管長(zhǎng)之間呈冪指數(shù)關(guān)系，即隨著管長(zhǎng)的增加，熱導(dǎo)率呈冪指數(shù)的規(guī)律衰減。關(guān)于缺陷影響的代表性研究有：Kondo等[14]模擬計(jì)算發(fā)現(xiàn)，在室溫下1%的空位缺陷會(huì)使碳納米管熱導(dǎo)率下降60%；馮妍卉等[15]模擬研究發(fā)現(xiàn)，含SW缺陷碳納米管的熱導(dǎo)率下降顯著。迄今為止，針對(duì)含空位缺陷和SW缺陷的碳納米管導(dǎo)熱性能的研究相對(duì)較多，而關(guān)于摻雜缺陷對(duì)碳納米管熱導(dǎo)率影響的研究較少。因此，有必要對(duì)管長(zhǎng)和摻雜缺陷對(duì)碳納米管熱導(dǎo)率的影響規(guī)律進(jìn)行進(jìn)一步研究。

筆者采用反向非平衡分子動(dòng)力學(xué)(Reversed Non-Equilibrium Molecular Dynamics，RNEMD)方法，針對(duì)手性為(10，10)的單壁碳納米管(Single-walled carbon nanotube， SWNT)，分別研究其熱導(dǎo)率隨著管長(zhǎng)和氮摻雜缺陷濃度的變化規(guī)律，探討碳納米管內(nèi)的熱輸運(yùn)機(jī)理，為相關(guān)碳納米材料熱物性研究提供理論基礎(chǔ)。

1 計(jì)算模型

模擬的對(duì)象是(10，10)SWNT，環(huán)境溫度為300 K。為了研究管長(zhǎng)對(duì)SWNT熱導(dǎo)率的影響規(guī)律，分別建立了10～100 nm(每個(gè)模型間隔10 nm)的完整無(wú)缺陷的SWNT物理模型。此外，為了研究氮摻雜缺陷濃度對(duì)SWNT熱導(dǎo)率的影響規(guī)律，利用將SWNT點(diǎn)陣結(jié)構(gòu)中的碳原子替換為氮原子來(lái)實(shí)現(xiàn)氮摻雜的方法，分別向管長(zhǎng)為10 nm的SWNT內(nèi)隨機(jī)摻雜了不同數(shù)量(10～100個(gè))的氮原子(以10個(gè)氮原子為增加值逐漸遞增)，對(duì)應(yīng)摻雜質(zhì)量分?jǐn)?shù)分別為0.31%、0.61%、1.22%、1.83%、2.44%、3.1%、3.66%、4.27%、4.88%、5.49%和6.10%。

計(jì)算模擬的物理模型如圖1所示。管內(nèi)碳-碳原子間的鍵長(zhǎng)為1.42 ?；碳納米管的橫截面積由式S=πDδ計(jì)算得到，其中：D為碳納米管直徑，δ為碳納米管壁厚(取范德瓦爾斯厚度δ=3.4 ?[16])。

勢(shì)函數(shù)的選取是分子動(dòng)力學(xué)模擬的關(guān)鍵，將直接決定熱導(dǎo)率預(yù)測(cè)的準(zhǔn)確性[17]。采用AIREBO勢(shì)函數(shù)描述碳—碳原子間的相互作用，形式如下：

(1)

2 模擬方法

傳統(tǒng)的非平衡分子動(dòng)力學(xué)方法是通過(guò)先在體系兩端施加溫差再產(chǎn)生熱流，而反向非平衡分子動(dòng)力學(xué)方法則是通過(guò)先施加恒定熱流進(jìn)而生成溫度梯度[19]。實(shí)踐表明，反向非平衡分子動(dòng)力學(xué)方法在非平衡狀態(tài)下的收斂速度更快，因此擬采用該方法來(lái)計(jì)算碳納米管的導(dǎo)熱系數(shù)。

采用Velocity-Verlet算法求解模擬中的運(yùn)動(dòng)方程，積分步長(zhǎng)為0.5 fs。首先，采用正則系綜(NVT)模擬50萬(wàn)步，調(diào)整系統(tǒng)的溫度并使其達(dá)到預(yù)定值；然后，采用微正則系綜(NVE)馳豫50萬(wàn)步，使系統(tǒng)達(dá)到平衡狀態(tài)；最后，繼續(xù)采用微正則系綜(NVE)模擬50萬(wàn)步，待系統(tǒng)穩(wěn)定后，統(tǒng)計(jì)平均求得系統(tǒng)的溫度分布與熱流密度。

溫度T由Boltzmann能量均分原理給出：

(2)

式中：N為模擬系統(tǒng)內(nèi)的原子數(shù)目；kB為玻爾茲曼常數(shù)，kB=1.38×10-23J/K；mi和vi分別是原子的質(zhì)量(kg)和速度(m/s)；尖括號(hào)表示系綜平均。

系統(tǒng)中熱流密度q(W/m2)可表示為：

(3)

式中：A為碳納米管軸向橫截面積(m2)；t為總模擬時(shí)間(s)；NE為交換原子的對(duì)數(shù)；vh和vc分別為熱浴中運(yùn)動(dòng)最慢的原子和冷浴中運(yùn)動(dòng)最快的原子的速率(m/s)；m為原子的質(zhì)量(kg)。

最后，根據(jù)傅里葉導(dǎo)熱定律計(jì)算熱導(dǎo)率：

(4)

式中：λ為熱導(dǎo)率，W/(m·K)；T為碳納米管的溫度梯度。

模擬系統(tǒng)的結(jié)構(gòu)如圖2所示，模擬計(jì)算中將碳納米管的物理模型沿著軸向分成N層，其中第1層為冷浴，第(N+1)/2層為熱浴；在模擬過(guò)程中，每隔一定的步數(shù)，就對(duì)冷浴中速度最大的原子和熱浴中速度最小的原子的動(dòng)能進(jìn)行一次交換，即能產(chǎn)生熱流。

計(jì)算中取10或40步，主要取決于碳納米管內(nèi)溫度是否會(huì)產(chǎn)生非線性效應(yīng)。由于在模擬中，冷浴和熱浴中的聲子會(huì)出現(xiàn)失配的現(xiàn)象，如圖3所示，從而使得溫度發(fā)生跳躍。因此，通過(guò)對(duì)兩邊線性溫度分布段進(jìn)行擬合來(lái)計(jì)算溫度梯度。

3 結(jié)果及討論

3.1 管長(zhǎng)對(duì)單壁碳納米管熱導(dǎo)率的影響

3.1.1 模擬結(jié)果

為了研究管長(zhǎng)對(duì)于SWNT熱導(dǎo)率的影響，分別計(jì)算了室溫下10～100 nm長(zhǎng)單壁碳納米管的熱導(dǎo)率，得到熱導(dǎo)率隨管長(zhǎng)的變化關(guān)系如圖4所示。由圖4中可以看出，當(dāng)外界溫度為300 K時(shí)，10 nm長(zhǎng)(10，10)SWNT的熱導(dǎo)率為89.7 W/(m·K)，碳納米管越長(zhǎng)，其熱導(dǎo)率越高。當(dāng)管長(zhǎng)到達(dá)100 nm時(shí)，SWNT的熱導(dǎo)率升高至537.4 W/(m·K)，升高了近5倍，熱導(dǎo)率有顯著提高。當(dāng)管長(zhǎng)低于40 nm時(shí)，SWNT的熱導(dǎo)率與管長(zhǎng)成正比；當(dāng)管長(zhǎng)高于40 nm后，隨著管長(zhǎng)增加，熱導(dǎo)率上升的速率逐漸減小。

同時(shí)，圖4中還給出了同樣以(10，10)SWNT為研究對(duì)象的Lukes[20]的模擬結(jié)果進(jìn)行對(duì)比分析。從圖4中可以看出，與Lukes的模擬結(jié)果相比，筆者模擬出的SWNT熱導(dǎo)率隨管長(zhǎng)的變化規(guī)律更明顯：不僅熱導(dǎo)率的數(shù)值更高，且熱導(dǎo)率的增加速率也更快。這是由于分子動(dòng)力學(xué)的模擬方法和勢(shì)函數(shù)模型的選取不同所致，筆者采用的是NEMD方法和AIREBO勢(shì)，而Lukes選擇的則是平衡態(tài)分子動(dòng)力學(xué)(EMD)方法和REBO勢(shì)。NEMD方法比EMD方法所需要的計(jì)算時(shí)間更短，精確度更高；同時(shí)，相比REBO勢(shì)，AIREBO勢(shì)的優(yōu)勢(shì)在于存在能夠準(zhǔn)確描述原子間長(zhǎng)程相互作用力的E_LJ項(xiàng)，可以全面地描述C—C原子間的作用力，提高了計(jì)算精度。

3.1.2 機(jī)理探討

計(jì)算得到的室溫下管長(zhǎng)對(duì)SWNT熱導(dǎo)率的影響規(guī)律為：熱導(dǎo)率隨著管長(zhǎng)增加而升高，但升高速率卻逐漸減緩；當(dāng)管長(zhǎng)低于40 nm時(shí)，隨著管長(zhǎng)增加，熱導(dǎo)率呈近似的線性升高(該階段斜率約為7.85)，可以認(rèn)為導(dǎo)熱處于彈道輸運(yùn)階段；當(dāng)管長(zhǎng)從40 nm增至100 nm時(shí)，熱導(dǎo)率依然繼續(xù)升高，但升高速率卻逐漸減小(該階段斜率從原來(lái)的7.85降至3.54)，說(shuō)明導(dǎo)熱處于彈道—擴(kuò)散輸運(yùn)階段。

上述管長(zhǎng)對(duì)熱導(dǎo)率的影響規(guī)律從機(jī)理上可解釋為：碳納米管內(nèi)熱能的承載和傳遞主要是依靠聲子的運(yùn)動(dòng)，聲子的平均自由程一般直接決定碳納米管熱導(dǎo)率數(shù)值，當(dāng)管長(zhǎng)小于40 nm時(shí)，管內(nèi)聲子平均自由程大于整個(gè)系統(tǒng)的尺寸，故聲子傳熱出現(xiàn)彈道輸運(yùn)的特點(diǎn)；當(dāng)管長(zhǎng)大于40 nm后，隨著管長(zhǎng)的增加，系統(tǒng)尺寸逐漸接近聲子平均自由程大小，聲子發(fā)生邊界散射，尺寸效應(yīng)明顯，熱導(dǎo)率受到了限制，因此熱導(dǎo)率上升的速率逐漸減緩，導(dǎo)熱過(guò)程由彈道輸運(yùn)向擴(kuò)散輸運(yùn)轉(zhuǎn)變。

3.2 氮摻雜對(duì)單壁碳納米管熱導(dǎo)率的影響

3.2.1 模擬結(jié)果

分別對(duì)300 K下的不同質(zhì)量分?jǐn)?shù)氮摻雜的單壁碳納米管熱導(dǎo)率進(jìn)行了分子動(dòng)力學(xué)模擬，如圖5所示。由圖5中可以看出，隨著氮摻雜質(zhì)量分?jǐn)?shù)的增加，單壁碳納米管的熱導(dǎo)率先急劇下降，后逐漸趨于穩(wěn)定。當(dāng)?shù)獡诫s質(zhì)量分?jǐn)?shù)低于3.1%時(shí)，SWNT熱導(dǎo)率下降趨勢(shì)明顯。尤其是氮摻雜質(zhì)量分?jǐn)?shù)在0～0.61%時(shí)，熱導(dǎo)率急劇下降，摻入的氮原子質(zhì)量分?jǐn)?shù)為0.61%時(shí)，熱導(dǎo)率為74.0 W/(m·K)，相比于相同條件下完整無(wú)缺陷時(shí)的熱導(dǎo)率(89.7 W/(m·K))降低了約20%。當(dāng)?shù)獡诫s質(zhì)量分?jǐn)?shù)達(dá)到3.1%時(shí)，SWNT的熱導(dǎo)率低至37.06 W/(m·K)，相比于相同條件下完整無(wú)缺陷時(shí)的熱導(dǎo)率(89.7 W/(m·K)降低了約60%。值得注意的是，當(dāng)?shù)獡诫s質(zhì)量分?jǐn)?shù)超過(guò)3.1%以后，單壁碳納米管的熱導(dǎo)率隨摻氮質(zhì)量分?jǐn)?shù)的增加而變化越來(lái)越平緩，其數(shù)值在37 W/(m·K)附近小幅浮動(dòng)。

3.2.2 機(jī)理探討

聲子平均自由程是研究單壁碳納米管熱輸運(yùn)機(jī)理的重要參考數(shù)據(jù)。根據(jù)經(jīng)典晶格熱輸運(yùn)理論，固體材料的熱導(dǎo)率κ可表示為：

κ=Cυl

(5)

式中：C為聲子的單位體積比熱；υ為聲子的速度；l為聲子的平均自由程。

模擬單壁碳納米管的長(zhǎng)度取為10 nm，遠(yuǎn)小于碳納米管聲子的平均自由程，而此時(shí)SWNT內(nèi)的熱傳導(dǎo)過(guò)程處近似于彈道輸運(yùn)階段。由于在傳統(tǒng)彈道輸運(yùn)理論中，聲子的比熱C和速度υ對(duì)固體內(nèi)部線性聲子的熱傳輸特性影響可以忽略不計(jì)，因此聲子平均自由程l將是單壁碳納米管的熱導(dǎo)率的決定性因素。同時(shí)由式(5)可知，單壁碳納米管的熱導(dǎo)率與聲子平均自由程成正比，而聲子平均自由程的大小又由2個(gè)過(guò)程決定：一是聲子間的相互碰撞；二是固體中缺陷對(duì)聲子的散射。

因此，關(guān)于氮摻雜缺陷造成SWNT熱導(dǎo)率整體下降的機(jī)理可以解釋為：第一，氮原子質(zhì)量比碳原子的質(zhì)量大，向單壁碳納米管內(nèi)摻氮后會(huì)導(dǎo)致晶格振動(dòng)的非簡(jiǎn)諧效應(yīng)，從而降低了聲子的移動(dòng)速度，同時(shí)減小了聲子平均自由程，故熱導(dǎo)率下降；第二，氮摻雜會(huì)加強(qiáng)聲子與氮原子間的相互作用，導(dǎo)致聲子與氮原子間的散射增多，使氮原子周圍的局部熱阻增大，且Umklapp散射會(huì)隨著摻氮缺陷質(zhì)量分?jǐn)?shù)的增大而增強(qiáng)，進(jìn)而抑制熱量的傳輸，導(dǎo)致熱導(dǎo)率下降；第三，聲子在傳輸熱能的過(guò)程中是連續(xù)的，當(dāng)有2種不同種類的原子存在時(shí)(氮原子和碳原子)，2種不同聲子模式之間會(huì)發(fā)生相互轉(zhuǎn)換，在轉(zhuǎn)換過(guò)程中會(huì)出現(xiàn)聲子的倒逆、反射和散射現(xiàn)象，聲子的傳熱能力也因此下降，故熱導(dǎo)率也隨之下降。

總之，氮摻雜缺陷破壞了單壁碳納米管的理想晶格結(jié)構(gòu)，不可避免地造成聲子平均自由程減小，最終都致使碳管熱導(dǎo)率整體下降；當(dāng)摻入的氮原子質(zhì)量分?jǐn)?shù)達(dá)到了一定程度后，氮原子的聲子傳熱模式將逐步取代碳原子的聲子模式，待氮原子的聲子模式在碳納米管內(nèi)的聲子傳熱中占主導(dǎo)地位后，SWNT熱導(dǎo)率也最終隨之趨于穩(wěn)定。

3 結(jié)論

采用反向非平衡分子動(dòng)力學(xué)方法，在室溫300 K下計(jì)算了(10，10)單壁碳納米管的熱導(dǎo)率，研究了管長(zhǎng)和氮摻雜對(duì)單壁碳納米管熱導(dǎo)率的影響規(guī)律，主要結(jié)論如下：

(1)當(dāng)管長(zhǎng)從10 nm增大到100 nm時(shí)，單壁碳納米管的熱導(dǎo)率從89.7 W/(m·K)增至537.4 W/(m·K)，升高了近5倍，熱導(dǎo)率有顯著提高；當(dāng)管長(zhǎng)小于40 nm時(shí)，SWNT的熱導(dǎo)率與管長(zhǎng)成正比；當(dāng)管長(zhǎng)大于40 nm后，隨著管長(zhǎng)增加，熱導(dǎo)率上升的速率逐漸減小，管內(nèi)導(dǎo)熱從彈道輸運(yùn)逐漸向彈道—擴(kuò)散輸運(yùn)轉(zhuǎn)變。

(2)隨著氮摻雜缺陷質(zhì)量分?jǐn)?shù)的增加，單壁碳納米管的熱導(dǎo)率先急劇下降后逐漸趨于穩(wěn)定；氮摻雜缺陷明顯降低了SWNT的熱導(dǎo)率，氮摻雜質(zhì)量分?jǐn)?shù)較低時(shí)熱導(dǎo)率的下降趨勢(shì)尤為顯著。當(dāng)?shù)獡诫s質(zhì)量分?jǐn)?shù)為3.1%時(shí)，單壁碳納米管的熱導(dǎo)率相比相同條件下完整無(wú)缺陷時(shí)顯著下降，達(dá)到了60%；而隨著氮摻雜缺陷質(zhì)量分?jǐn)?shù)的繼續(xù)增加，氮原子的聲子模式在管內(nèi)聲子傳熱中占主導(dǎo)地位，熱導(dǎo)率逐漸趨于穩(wěn)定。