氮化硼納米管熱輸運性能的分子動力學模擬

高宇飛 孟慶元 張 璐 劉甲秋 荊宇航,*

(1哈爾濱工業大學航天學院,哈爾濱150001;2中國醫科大學,沈陽110001;3哈爾濱玻璃鋼研究院,哈爾濱150001)

氮化硼納米管熱輸運性能的分子動力學模擬

高宇飛1孟慶元1張 璐2劉甲秋3荊宇航1,*

(1哈爾濱工業大學航天學院,哈爾濱150001;2中國醫科大學,沈陽110001;3哈爾濱玻璃鋼研究院,哈爾濱150001)

采用基于聲子散射理論的Boltzmann-Peierls聲子傳輸方程(BTE)和非平衡態分子動力學模擬(NEMD)方法研究了氮化硼納米管(BNNT)的熱輸運性能.分析了BNNT的熱力耦合效應,通過BTE與NEMD兩種方法相結合,分析了溫度和長度對BNNT熱輸運性能的影響,并應用量子修正擴大了NEMD的研究范圍.結果表明:隨著拉伸或壓縮應變的增加,BNNT熱輸運性能均呈降低的趨勢.通過計算聲子態密度(PDOS)在理論上分析了以上結果,發現在拉伸狀態下,聲子模式的變化是決定BNNT熱輸運性能變化的主要因素;在壓縮狀態下,熱導率變化是由于模型發生明顯的屈曲變形引起的.在低溫段,BNNT的熱輸運性能受量子效應影響最初有一個線性增加的過程,當溫度超過一定值時,其開始顯著地降低;當BNNT長度小于120 nm時,隨著長度的增加,其彈道性能逐漸減弱,但仍主要體現為彈道特征,其熱導率(κ)與長度(L)基本滿足κ∝Lα這一關系.

氮化硼納米管;熱輸運性能;Boltzmann-Peierls聲子傳輸方程;非平衡態分子動力學模擬;聲子態密度;量子修正

1 引言

納米材料因其獨特的結構和優良的力學、電學、熱學、光學等性能,表現出廣闊的應用前景.至1991年發現碳納米管(CNT)以來,1各種不同種類的納米管相繼被發現,對它們的研究也隨之展開.

1994年Rubio等2預測到一定存在氮化硼納米管(BNNT).之后人們在實驗中得到了多壁BNNT. 1996年Loiseau等3采用電弧放電法得到了單壁BNNT.有關BNNT及氮硼化合物的研究,如電子結構和力學性能的研究已經展開,4但有關其熱學性能的研究卻較少,在力學和熱學交叉領域,如熱力耦合效應方面還未見報道.BNNT作為一種典型的納米材料,它的很多結構參數都對其熱導率,尤其是軸向熱導率有著顯著的影響,5因此,可以利用BNNT的這種優良特性制成可控制的導熱設備,這與傳統的通過選擇不同的材料來實現對器件導熱性能控制的方法相比有了革命性的突破.因此有必要對BNNT的熱傳導機理和導熱性能進行系統的研究.

目前針對BNNT導熱性能的研究成果指出:對BNNT的熱輸運起主要貢獻的是聲子,這與CNT的熱輸運機理不同(對CNT熱導率起主要作用的是聲子和電子).BNNT的熱導率低于CNT,但BNNT比CNT有著更強的熱穩定性和化學穩定性;Liew和Yuan6通過分析CNT和BNNT混合物發現,BNNT具有更好的穩定性,較之CNT更適宜做管狀結構的外壁;然而也有部分研究7指出,在低溫段,BNNT的熱導率要高于CNT.此外,同位素對BNNT熱導率有著顯著的影響,Zettl等8,9的研究指出,當BNNT中含有高濃度的B同位素時,其熱導率可與CNT相比;純的11B同位素可將BNNT熱導率提升50%;10而Barman11的研究表明,同位素10B會削弱氮化硼晶體的熱導率.BNNT具有較好的導熱性能,在室溫下,外直徑為30-40 nm的BNNT的熱導率可達到350 W·m-1·K-1;Stewart等12通過計算表明(8,0) BNNT的熱導率為200 W·m-1·K-1.Nakayama等13發現氮化硼納米帶可以增強部分金屬氧化物的導熱性能.Mingo和Broido14通過研究輸運方程和彈道函數發現納米管的熱傳導性能表現出彈道特征.以上結果表明:BNNT在掃描探針顯微鏡、原子力顯微鏡的針尖,高溫陶瓷以及納米尺度熱量計和微電子處理器等方面有著很大的潛在應用.但由于實驗測量中的實際困難和材料在納米尺度下的誤差,目前理論和實驗針對納米管導熱性能的定量研究結果具有較大的差別,如有部分研究結果15,16表明BNNT具有極低的熱導率(熱導率約1.5-14 W·m-1·K-1).總體來說,目前對于BNNT熱輸運機理和性能的定量研究成果仍很缺乏,而且隨著聲子散射理論的逐步完善,一些對BNNT導熱性質產生重要影響的較復雜的因素,如納米管的彈道傳輸與擴散傳輸的分界性能、應變作用下動態納米管結構的熱輸運性能等逐漸被提出,并逐漸成為研究的熱點.

本文首先采用非平衡態分子動力學模擬(NEMD)方法研究了BNNT在軸向應變作用下的熱輸運性能,并通過計算聲子態密度(PDOS)的方法對模擬結果進行了理論上的補充和驗證.進一步通過求解Boltzmann-Peierls聲子傳輸方程(BTE)和NEMD相結合的方法分析了溫度和長度對BNNT熱輸運性能的影響,通過量子修正擴大了NEMD的研究范圍.為研究BNNT彈道傳輸與擴散傳輸的分界以及動態納米管結構的熱輸運性能提供了一些理論與計算模擬的參考數據.

2 計算模型與方法

2.1 非平衡態分子動力學模擬方法

本文在模擬過程中應用Tersoff多體勢17來描述BNNT中原子間的相互作用.求解原子的運動方程時采用預測—校正算法,積分步長為0.5 fs.邊界條件選擇為自由邊界條件.構建BNNT模型時,首先根據BNNT的空間結構特點得到其初始構型的空間排列,然后將BNNT的初始構型在0 K下弛豫50000步,隨后逐漸升溫,每次升溫100 K,并弛豫相同的步數,直到溫度升到1600 K,最后進行退火處理來得到其穩定結構.

計算模型的熱導率采用NEMD,18其基本思想為:將BNNT沿軸向等分為N層,設0層和N/2層為“冷區”和“熱區”(見圖1).每隔一定步數提取0層中動能最大的原子和N/2層中動能最小的原子,并交換兩者速度.一般情況下冷區中動能最大的原子的動能大于熱區中動能最小的原子的動能,這樣便產生了一個從0層到N/2層的熱流.又因為熱區的溫度高于冷區的溫度,因此便存在一個與熱流方向相反的溫差.基于上述模型,經過統計可獲得軸向各層的溫度,進一步計算可得溫度梯度(▽T),再結合能量交換過程中統計得到的熱流密度(J),由Fourier導熱定律J=-κ?T,可計算得到BNNT的熱導率(κ).整個非平衡態分子動力學模擬過程是通過Lammps開源軟件包19實現的.

圖1 非平衡態分子動力學模擬過程示意圖Fig.1 Schematic diagram of NEMD process

圖2為BNNT在300 K溫度下軸向方向的溫度分布圖,其中各層溫度是通過統計各層中所有原子的動能,通過公式(1)獲得

其中N為層中的原子數,T代表溫度,kB為Boltzmann常量,mi和vi分別為模型中i原子的質量和速度.從圖中可發現中間各層的溫度呈線性分布,而在冷層和熱層附近區域,則有明顯的非線性溫度響應,這是由于熱源和熱槽的強相互作用產生的.20因此在計算BNNT熱導率過程中,需要統計中間各層的溫度,進而獲得BNNT軸向方向的溫度梯度以及熱導率.

圖2 BNNT軸向(L)的溫度分布Fig.2 Temperature distribution of BNNT along axial direction(L)

模擬的整個過程分為五步:首先用Conjugate Gradient(CG)方法弛豫整個構型,以獲得BNNT的穩定結構;然后將納米管兩端固定,在NPT系綜下運行10000步(5 ps)標定壓力使其處于最小狀態;再在最小壓力狀態下運行100000步(50 ps),使系統達到穩定狀態;然后在NVT系綜下運行200000步(100 ps)來弛豫整個系統;最后在NVT系綜下運行1000000步(500 ps)統計熱流密度和溫度梯度,計算得出BNNT的熱導率(在最后一步過程中每100000步統計一次熱導率,通過計算10個熱導率數值的偏差分布可獲得模型熱導率的不確定度).

2.2 Boltzmann-Peierls聲子傳輸方程方法

微觀固體中的熱量主要是靠晶格振動產生的聲子傳輸的,聲子是一些被抽象出的離散粒子,單一聲子的能量為?ω,?=h/2π,h為普朗克常量,ω為頻率.Boltzmann-Peierls聲子傳輸方程14,21(BTE)可以完整地描述固體中聲子間的相互作用.但獲得BTE的精確解是困難的,目前可以通過一些近似的方法來獲得不同模型對應的BTE的近似解.求解BTE的基本思想是通過不同聲子的散射時間和處于非平衡態聲子的分布,并結合模型的邊界條件來獲得模型熱流與溫度梯度,最后求得模型的熱導率.

2.3 低溫段量子修正

一般情況下,經典分子動力學模擬(MD)方法在溫度接近和高于德拜溫度時被認為是有效的,因為在此溫度區域內所有的聲子模式都被完全激發出來.與此相對的,當溫度顯著地低于德拜溫度時,大量的聲子沒有被激發,模型的量子效應不能被忽略,此時經典的MD方法無效.為突破這一限制,量子修正方法(QCs)22-25被提出,目前它被廣泛應用于熱力學和力學的研究中.

QCs的基本思想是由經典力學系統中的模擬溫度來計算量子系統中對應的修正溫度,根據經典力學系統與量子系統能量相等,具體為系統的動能等于聲子能量的一半這一原則.模擬溫度(TMD)與修正溫度的關系可表示為

其中,N0為平衡態聲子分布函數,Tq為量子修正溫度,g(ω)為聲子態密度.結合模擬中獲得的TMD可計算Tq.通過統計不同時間的軸向熱流密度和溫度梯度,可獲得修正后的熱導率.

除了基于MD的QCs修正外,還存在兩種基于勢函數的量子修正方法,Feynman-Hibbs勢函數(FH勢)和Wigner-Kirkwood勢函數(WK勢)方法.FH勢是由量子均分方程的Feynman-Hibbs估算變量獲得的,是一種依賴于溫度的二體勢,可以方便地植入分子動力學或蒙特卡羅算法的計算程序中.應用這種方法,Sese26計算了Lennard-Jones系統的熱力學和結構特性,Tchouar等27應用由FH方法演化的Lennard-Jones勢研究了氖的熱力學、結構和傳輸特性,Goharshadi等28應用Hartree-Fock dispersion近似勢函數和FH勢研究了氦的熱力學性質.WK勢是由對均分方程中德布羅意波長的擴展獲得的,可以被用來研究不同物理特性系統的量子效應.

3 計算結果與討論

3.1 熱力耦合效應

本文選取20 nm的20_20扶手椅型BNNT(原子數為6320個)為研究對象,分別對模型施加在0-0.1范圍的軸向拉伸和壓縮應變.施加位移載荷過程中,應變是通過重新標定模型軸向原子坐標實現的,具體為首先選擇BNNT軸向方向最外部兩層原子,分別對其施加反方向的軸向位移載荷,每施加完一次位移載荷讓系統弛豫一定步數,使系統達到平衡狀態.加載速率為每5000個步長施加0.0025 nm的位移.在此期間統計不同時間BNNT軸向方向的熱流密度和溫度梯度,最終計算得到熱導率隨應變變化趨勢如圖3所示.由圖3可發現,在拉伸和壓縮應變作用下,BNNT的軸向熱導率隨著應變的增加而逐漸減小.

圖3 BNNT軸向熱輸運性能隨軸向應變的變化曲線Fig.3 Variation curve of axial thermal conductivity of BNNT under axial strain

圖4 BNNT在軸向拉伸應變作用下的聲子態密度圖Fig.4 PDOS of BNNT under axial tensile strain

聲子態密度(PDOS)可以完整地描述各種模式聲子的分布,可用于研究動態結構材料的熱輸運性能.PDOS是通過對速度自相關函數進行傅里葉變換得到的.29-32速度自相關函數可以表示為

其中

對(4)式進行傅里葉變換可得到PDOS的表達形式為

其中t代表時間,k=x,y,z代表模型的三個方向.通過上式,結合MD中統計的原子速度,可以獲得不同應變作用下BNNT的聲子態密度分布圖.BNNT在拉伸和壓縮應變作用下的聲子態密度分布如圖4、 5所示.

圖5 BNNT在軸向壓縮應變作用下的聲子態密度高頻波峰圖Fig.5 High frequency peak in PDOS of BNNT under axial compressive strain

由圖4可見,PDOS中有三個較高的波峰,其中前兩個位于低頻段,可被視為代表低頻聲子模式的波峰,后一個位于高頻段,被認為是代表高頻聲子模式的波峰.隨著應變的變化,相比于高頻波峰,兩個低頻波峰的變化較小.這意味著高頻聲子對應變更加敏感.由于波峰的變化直接反映了聲子數量的增加或減少,以及振動模式的增強或減弱,對比圖4和圖3我們可以發現,高頻波峰的左移趨勢伴隨著模型熱導率的下降,因此可以認為高頻聲子對模型熱輸運性能有著主要影響.

在0-0.1的拉伸應變作用下,隨著拉伸應變的增加,從圖4可以看出高頻聲子波峰逐漸左移.這意味著高頻聲子數量逐漸減少,部分高頻聲子逐漸轉化為低頻聲子,它將導致聲子的平均頻率降低,因此此過程可以視為一個紅移過程.Xu和Buehler33研究CNT在軸向應變作用下的熱輸運性能,發現紅移過程將會減小聲子的速度.在聲子的彈道輸運過程中,模型的熱導率(κ)與聲子速度有以下關系

其中c為比熱,l為平均自由程,v為聲子速度.因此圖4中高頻聲子的變化趨勢將導致BNNT的熱導率隨著拉伸應變的增加而逐漸減小.觀察圖5可以發現,在壓縮應變作用下,不同應變的PDOS變化較小,幾乎沒有明顯的變化.結合圖6中0.1的拉伸和壓縮應變作用下BNNT的構型圖可以發現,在0.1的拉伸應變作用下BNNT被逐漸拉長,但其沒有顯著的塑性變形,而在0.1的壓縮應變作用下,BNNT發生了嚴重的塑性變形.隨著明顯的塑性變形的產生,系統的能量將顯著降低,進而導致模型中聲子平均速度的降低,根據公式(6),在壓縮應變下, BNNT的熱導率隨著應變的增加而逐漸減小.這也說明了在壓縮應變下BNNT的熱導率發生變化的機理不同于拉伸應變下的情況:在拉伸應變下,是聲子模式的變化導致了BNNT熱導率的變化;在壓縮應變下,是BNNT結構的顯著變化導致了其熱導率的變化.

圖6 BNNT在(a)0.1的拉伸應變作用和(b)0.1的壓縮應變下的構型圖Fig.6 Structures of BNNT under(a)0.1 tensile strain and (b)0.1 compressive strain

Buehler33和Zhang32等分別研究了CNT在應變作用下導熱性能,本文在拉伸應變作用下得到的PDOS與Buehler和Zhang的結果有著相同的變化趨勢,高頻聲子波峰均隨著拉伸應變的增加而逐漸左移.在壓縮應變的作用下,Zhang計算了-0.02-0壓縮應變下CNT的PDOS,發現其高頻聲子波峰有部分右移的趨勢.本文認為此處的不同是由于CNT的抗壓縮性能強于BNNT導致的.在-0.02的小壓縮應變作用下CNT較之BNNT更不容易產生較大的塑性變形,因此CNT的能量減少將體現在聲子模式的變化上,而BNNT的能量減少體現在結構的變化上.

3.2 溫度影響

本文應用NEMD,結合低溫段溫度量子修正方法研究了10-1600 K溫度范圍內長度為20 nm的10_10扶手椅型BNNT(原子數為3160個)的導熱性能.首先統計模擬過程中BNNT的溫度TMD,結合公式(2),可獲得量子修正溫度Tq.量子修正溫度Tq與模擬溫度TMD的關系見圖7.其中實線為Tq隨TMD的變化曲線,虛線用來對比Tq隨TMD的變化率,斜率為1.從圖中可見在低溫段,?Tq/?TMD遠大于1,隨著溫度的升高,其值逐漸減小,當溫度達到1000 K時,?Tq/?TMD≈1,這說明此時Tq與TMD的變化率一致.結合Tq,可獲得BNNT的熱導率隨溫度的變化趨勢如圖8所示.

從圖8可以看出,BNNT的軸向熱導率在低溫段(10-70 K)隨著溫度的升高近似呈線性增長趨勢;當溫度超過一定值(近似為100 K)時,其熱導率隨著溫度的升高逐漸降低.

在低溫段,當系統溫度遠小于德拜溫度時,由于溫度較低,晶格振動較弱,被激發出的聲子較少.因此在低溫段,材料的邊界對長波聲子的散射(也即邊界散射)是影響材料熱導率的主要因素.隨著系統溫度的降低,具有低頻率、長波長的聲子在聲子中的比重增加.根據聲子理論,長波聲子較易與邊界發生散射,這就導致了聲子平均自由程減小,由(6)式可知,熱導率也隨之減小,因此在低溫段,BNNT的熱導率隨著溫度的增加而增加.在高溫段,晶格振動劇烈,散發出大量的聲子,因此聲子與聲子間的相互碰撞即聲子間的散射取代邊界對長波聲子的散射成為影響材料熱導率的主要因素.隨著溫度的增高,聲子數增多,聲子間的碰撞概率增大,聲子平均自由程減小,因此在高溫段,BNNT的熱導率隨著溫度的增加而逐漸減小.

圖7 模擬溫度(TMD)與量子修正溫度(Tq)的關系Fig.7 Relationship between simulated temperature(TMD) and quantum corrected temperature(Tq)The slope of dotted line is 1 and can be used to contrast the change rate of corrected temperature with simulated temperature.

圖8 BNNT軸向熱輸運性能隨溫度的變化曲線Fig.8 Variation of axial thermal conductivity of BNNT under different temperatures

由圖8還可觀察到,在低溫段,修正后的熱導率與分子模擬直接獲得的熱導率有較大差別,隨著溫度的增加,這種差距逐漸減小.文獻23對比了量子修正前后硅晶體的熱導率變化曲線,當溫度在800 K左右時,兩者誤差約為10%.本文用兩種方法獲得的數據誤差為8%.文獻34用MD計算了1200-2000 K下(10,10)BNNT的熱導率,其熱導率隨著溫度升高逐漸降低,且1200 K時BNNT熱導率近似為48 W·m-1·K-1,1400 K時近似為43 W·m-1·K-1,1600 K時近似為38 W·m-1·K-1.計算得到上述三個溫度的熱導率分別為46、40、37 W·m-1·K-1,與文獻34結果較為接近.文獻35中通過理論分析BNNT的熱導率隨溫度變化規律,發現幾種直徑(3,3)(5,0)和幾種軸向長度(1,5,10 μm)的BNNT的熱導率發生轉折時系統的溫度均在70-100 K范圍內,而本文熱導率發生轉折的溫度在70-90 K之間,這可從一個角度印證本文結果,可證明本文方法的正確性.

本文通過理論求解Boltzmann-Peierls聲子傳輸方程得出了在低溫段BNNT的導熱性能與外界溫度的關系.

BTE聲子傳輸方程14可以表示為

左項描述聲子漂移部分,右項描述聲子散射部分,兩項可分別表示為36,37

其中N為聲子分布函數,N0為平衡態聲子分布函數, n可表示為n=N-N0,它用來描述實際聲子分布偏離平衡態的部分,τtot為總體弛豫時間,所以(7)式可以轉化為

結合公式(10),熱流可表示為

其中s為聲子模式,q為波數,ω為頻率.考慮熱流和溫度梯度的關系38

熱導率可表示為

其中V為體積.在低溫段,影響熱導率的主要因素是聲子與邊界的散射,近似認為聲子與邊界的散射時間是散射時間τtot的主要部分.邊界散射時間可表示為

其中p為描述邊界粗糙程度的參數,Λ為描述材料橫截面積的參數.

本文主要考慮邊界散射對熱導率的影響,可近似認為在低溫狀態下,τtot≈τb.并將公式(13)轉化為相對于波數q的積分形式

所以在低溫段納米管的熱導率與溫度符合κ∝T關系.

3.3 尺寸效應

根據聲子散射理論,模型軸向長度的變化將因為其逐漸接近或者遠離聲子平均自由程而對聲子在熱輸運過程中的散射產生影響,進而影響模型的熱輸運性能.本文建立了長度為10-120 nm的10_ 10扶手椅型BNNT,模擬溫度選擇為300 K,獲得結果如圖9所示.

如圖9可知,BNNT的軸向熱導率隨著軸向長度的增加而增加,且增加的速率逐漸減小.

圖9 BNNT軸向熱輸運性能與長度的關系Fig.9 Relationship between axial thermal conductivity and length

根據聲子散射理論,隨著BNNT軸向長度的增加,聲子在傳輸過程中與邊界的碰撞概率逐漸減小,這將增加聲子平均自由程,根據公式(6),BNNT的熱輸運性能逐漸增加;隨著BNNT軸向長度逐漸接近聲子平均自由程,聲子間相互碰撞概率逐漸增大,這對隨著長度增加而逐漸減小的聲子與邊界的碰撞概率起到一個抵消作用,根據上文分析,BNNT的熱輸運性能增加趨勢逐漸減緩.Mingo和Broido14通過計算CNT的BTE得到熱導率與長度呈κ∝Lα這一關系,且隨著L的增加,α值(描述模型彈道特性的參數)逐漸減小.本文通過MD獲得的結果如擬合為κ∝Lα形式,其α值大致介于0.1989-0.4899之間,且隨著長度的增加,α值逐漸減小,這與Mingo的理論結果相吻合.同時這也意味著隨著BNNT長度的增加,模型的彈道性能逐漸減弱,但其仍主要體現彈道特征.

4 結論

運用NEMD研究了BNNT的熱輸運性能,分別研究了長度、溫度和應變對BNNT的熱輸運性能的影響,并得到了以下結論:在拉伸應變作用下,隨著應變的增加,受聲子模式變化的影響,BNNT的熱導率呈逐漸降低的趨勢;在壓縮應變作用下,隨著應變的增加,受BNNT發生顯著塑性變形的影響,其熱導率逐漸減小.當外界溫度低于70 K時,受量子效應的影響,BNNT的熱導率隨著溫度的增加呈線性增加趨勢;當溫度高于100 K時,其熱導率隨著溫度的增加而逐漸減小.當BNNT軸向長度小于120 nm時,隨著長度的增加,其彈道性能逐漸減弱,但仍主要體現為彈道特征,其熱導率與長度基本滿足κ∝Lα這一關系.

(1) Iijima,S.Nature 1991,354,56.

(2) Rubio,A.;Corkill,J.L.;Cohen,M.L.Phys.Rev.B 1994,49 (7),5081.

(3) Loiseau,A.;Willaime,F.;Demoncy,N.;Hug,G.;Pascard,H. Phys.Rev.Lett.1996,76(25),4737.

(4) Wirtz,L.;Rubio,A.Phys.Rev.B 2003,68(4),045425.

(5) Zhang,G.;Li,B.W.The Journal of Physical Chemistry B 2005, 109(50),23823.

(6) Liew,K.M.;Yuan,J.H.Nanotechnology 2011,22,085701.

(7) Xiao,Y.;Yan,H.;Cao,J.X.;Ding,J.W.;Mao,Y.L.;Xiang,J. Phys.Rev.B 2004,69(20),205415.

(8) Chang,C.W.;Fennimore,A.M.;Afanasiev,A.;Okawa,D.; Ikuno,T.;Garcia,H.;Li,D.Y.;Majumdar,A.;Zettl,A.Phys. Rev.Lett.2006,97,085901.

(9) Chang,C.W.;Okawa,D.;Majumdar,A.;Zettl,A.Science 2006,314,1112.

(10) Zhi,C.Y.;Bando,Y.;Tang,C.C.;Golberg,D.Materials Science and Engineering R 2010,70,92.

(11) Barman,S.Europhysics Letters 2011,96,16004.

(12) Stewart,D.A.;Savic,N.;Mingo,N.Nano Lett.2009,9,81.

(13) Cho,H.B.;Tokoi,Y.;Tanaka,S.;Suematsu,H.;Suzuki,T.; Jiang,W.H.;Niihara,K.;Nakayama,T.Composites Science and Technology 2011,71,1046.

(14) Mingo,N.;Broido,D.A.Nano Lett.2005,5(7),1221.

(15) Chang,C.W.;Han,W.Q.;Zettl,A.Journal of Vacuum Science &Technology B 2005,23,1883.

(16) Tang,C.C.;Bando,Y.;Liu,C.H.;Fan,S.S.;Zhang,J.;Ding, X.X.;Golberg,D.The Journal of Physical Chemistry B 2006, 110(21),10354.

(17) Tersoff,J.Phys.Rev.B 1988,37(12),6991.

(18)Wang,X.W.;Huang,Z.;Wang,T.;Tang,Y.W.;Zeng,X.C. Physica B 2008,403,2021.

(19) Sandia National Laboratories.LAMMPS,Lammps-12Oct10; GNU:USA,2010.

(20) Schelling,P.K.;Phillpot,S.R.;Keblinski,P.Phys.Rev.B 2002, 65(14),144306.

(21) Lindsay,L.;Broido,D.A.;Mingo,N.Phys.Rev.B 2009,80, 125407.

(22) Wang,S.C.;Liang,X.G.;Xu,X.H.;Ohara,T.J.Appl.Phys. 2009,105,014316.

(23)Turney,J.E.;McGaughey,A.J.H.;Amon,C.H.Phys.Rev.B 2009,79(22),224305.

(24) Munoz,E.;Lu,J.X.;Yakobson,B.I.Nano Lett.2010,10(5), 1652.

(25) Hu,J.N.;Ruan,X.L.;Chen,Y.P.Nano Lett.2009,9(7),2730.

(26) Sese L.M.Mol.Phys.1993,78,1167.

(27) Tchouar,N.;Ould,K.F.;Levesque,D.J.Chem.Phys.2004, 121,7326.

(28) Goharshadi,E.K.;Abbaspour,M.;Kashani,H.;Baherololoom, M.Theor.Chem.Acc.2008,119,355.

(29) Yang,N.;Zhang,G.;Li,B.W.Nano Lett.2008,8(1),276.

(30) Kang,D.D.;Hou,Y.;Dai,J.Y.;Yuan,J.M.Phys.Rev.A 2009, 79,063202.

(31) Huang,M.J.;Weng,C.C.;Chang,T.M.International Journal of Thermal Sciences 2009,49,1095.

(32) Ren,C.L.;Zhang,W.;Xu,Z.J.;Zhu,Z.Y.;Huai,P. The Journal of Physical Chemistry C 2009,114(13),5786.

(33) Xu,Z.P.;Buehler,M.J.Nanotechnology 2009,20,185701.

(34) Shen,H.J.Computational Materials Science 2009,47,220.

(35) Jiang,J.W.;Wang,J.S.Phys.Rev.B 2011,84(8),085439.

(36) Goicochea,J.V.;Madrid,M.;Amon,C.ASME J.Heat Transfer 2010,132,012401.

(37)Chen,Y.F.;Li,D.Y.;Lukes,J.R.;Majumdar,A.ASME J.Heat Transfer 2005,127,1129.

(38) Nika,D.L.;Pokatilov,E.P.;Askerov,A.S.;Balandin,A.A. Phys.Rev.B 2009,79(15),155413.

November 30,2011;Revised:January 28,2012;Published on Web:February 27,2012.

Molecular Dynamics Simulation of Thermal Transport Properties for Boron Nitride Nanotubes

GAO Yu-Fei1MENG Qing-Yuan1ZHANG Lu2LIU Jia-Qiu3JING Yu-Hang1,*
(1School of Astronautics,Harbin Institute of Technology,Harbin 150001,P.R.China;2China Medical University, Shenyang 110001,Liaoning Province,P.R.China;3Harbin FRP Research Institute,Harbin 150001,P.R.China)

The Boltzmann-Peierls phonon transport equation(BTE)and non-equilibrium molecular dynamics simulation(NEMD)are used to investigate the thermal transport properties of boron nitride nanotubes(BNNTs).First,the thermal-mechanical coupling is explored using NEMD.Then,by combining BTE and NEMD,the influence of temperature and length is investigated.Quantum correction is used to extend the range over which NEMD can be used.The results demonstrate that under low-strain conditions,the thermal conductivity decreases with increasing tensile or compressive strain.Then the phonon density of state(PDOS)is used to analyze the trends in thermal transport properties theoretically; it is found that the variations in thermal transport properties under tension are caused by changes in the phonon modes,and that under compression changes are induced by the flection of the BNNT structure. The BNNT thermal conductivity increases linearly with increasing temperature because of the quantum effect at low temperatures,and it decreases significantly as the temperature reaches a certain value. When the BNNT length is less than 120 nm,the BNNT?s ballistic characteristics weaken with increasing length,but it also performs ballistic characteristic mainly,and thermal conductivity(κ)and length(L)obey the relationshipκ∝Lα.

Boron nitride nanotubes;Thermal transport property;Phonon Boltzmann-Peierls transport equation;Non-equilibrium molecular dynamics simulation;Phonon density of state;Quantum correction

10.3866/PKU.WHXB201202273

?Corresponding author.Email:jingyuhang@gmail.com;Tel:+86-451-86414143.

The project was supported by the National Natural Science Foundation of China(10772062).

國家自然科學基金(10772062)資助項目

O641;O642