黃桂芳，詹斯琦，黃維清

（湖南大學(xué) 物理與微電子科學(xué)學(xué)院，湖南 長(zhǎng)沙 410082）

近年來(lái)，一維非線性晶格體系的熱傳導(dǎo)性質(zhì)引起了廣泛的關(guān)注．人們?cè)噲D利用解析或非平衡分子動(dòng)力學(xué)方法從微觀動(dòng)力學(xué)機(jī)理的角度去解釋宏觀的熱傳導(dǎo)現(xiàn)象［1－6］．盡管至今沒(méi)能得到一個(gè)完善的關(guān)于低維材料的熱傳導(dǎo)理論，但在經(jīng)典一維系統(tǒng)的熱傳導(dǎo)研究過(guò)程中，發(fā)現(xiàn)許多新奇的現(xiàn)象，如熱整流效應(yīng)，界面熱阻，負(fù)微分熱阻等．這些現(xiàn)象的存在使得熱量的載體——聲子能像電子那樣得到控制與管理，因此開(kāi)辟了熱二極管［7］，熱晶體管［8］，熱邏輯門(mén)［9］等熱學(xué)器件的研究領(lǐng)域．其中熱二極管的研究最為廣泛，人們?cè)噲D建立各種模型來(lái)實(shí)現(xiàn)材料的熱整流性質(zhì)［10－12］．經(jīng)過(guò)大量的理論研究后，總結(jié)出熱二極管模型需要滿(mǎn)足兩個(gè)必要條件：系統(tǒng)的非線性（簡(jiǎn)諧）及不對(duì)稱(chēng)結(jié)構(gòu)．此外，實(shí)驗(yàn)方面也取得了不錯(cuò)的成績(jī)，借助納米科技的迅速發(fā)展，Chang等人［13］首次在實(shí)驗(yàn)中成功制備出碳納米管固體熱二極管．隨后不同形態(tài)及材料構(gòu)成的熱二極管也相繼在實(shí)驗(yàn)中成功制備［14－15］，這使熱二極管的實(shí)際應(yīng)用成為可能．

功能梯度材料在各領(lǐng)域中得到廣泛的應(yīng)用，它是一種在結(jié)構(gòu)或性能上連續(xù)／準(zhǔn)連續(xù)變化的非均質(zhì)復(fù)合材料．對(duì)于功能梯度材料的研究，目前已逐步擴(kuò)大到機(jī)械、電子、聲學(xué)、光學(xué)、核工程、生物等領(lǐng)域．相比這些性質(zhì)，人們對(duì)它的熱學(xué)性質(zhì)的研究并不夠多．文獻(xiàn)［16］最先將功能梯度材料概念引入一維非線性熱傳導(dǎo)領(lǐng)域，并提出一維線性質(zhì)量梯度下的Fermi－Pasta－Ulam（FPU）鏈，通過(guò)計(jì)算實(shí)現(xiàn)了熱整流效應(yīng)．隨后線性質(zhì)量梯度下的Frenkel－Kontorova（FK）鏈及指數(shù)型質(zhì)量梯度下的FPU 鏈也相繼被提出用以構(gòu)建熱二極管模型［17－18］．在實(shí)驗(yàn)方面，Chang等人［13］也是通過(guò)采用質(zhì)量梯度沉積制備出了不對(duì)稱(chēng)的無(wú)定型C6H16Pt固體熱二極管．由此可見(jiàn)，功能梯度材料在低維非線性熱傳導(dǎo)領(lǐng)域極具應(yīng)用價(jià)值．

本文利用非平衡分子動(dòng)力學(xué)模擬方法，系統(tǒng)地研究一維非線性功能梯度材料中不對(duì)稱(chēng)熱傳導(dǎo)性質(zhì)．在以往的研究中，改變系統(tǒng)的熱整流方向需要進(jìn)行結(jié)構(gòu)重建．然而，本文的模型能通過(guò)調(diào)節(jié)熱浴溫差進(jìn)行熱整流方向的控制，即系統(tǒng)出現(xiàn)了熱整流反轉(zhuǎn)現(xiàn)象．熱整流的反轉(zhuǎn)是由負(fù)微分熱阻引起的，并且它還與系統(tǒng)的體系大小及界面耦合強(qiáng)度緊密相關(guān)．鑒于該模型能在原地控制熱整流方向，文章進(jìn)一步討論了其實(shí)際應(yīng)用的可能性．

1 模型和計(jì)算方法

圖1為系統(tǒng)結(jié)構(gòu)示意圖，該系統(tǒng)由兩段質(zhì)量連續(xù)變化的FK鏈組成，并由一彈性系數(shù)為kint的簡(jiǎn)諧彈簧連接，我們稱(chēng)該彈簧為界面彈簧．系統(tǒng)通過(guò)引入線性遞減分布的質(zhì)量梯度來(lái)實(shí)現(xiàn)不對(duì)稱(chēng)性（功能梯度），為方便觀察，質(zhì)量大小由晶格的大小表示．最左端與最右端晶格分別與溫度為T(mén)L和TR的熱浴接觸．

系統(tǒng)的哈密頓量為：

式中：N為系統(tǒng)總粒子數(shù)；HLFK和HRFK分別為左右兩段FK 鏈的哈密頓量．

圖1 系統(tǒng)結(jié)構(gòu)示意圖Fig．1 Schematic picture of the system

式中：pi和xi分別為第i個(gè)晶格的動(dòng)量和位置．取k＝a＝1，分別表示晶格間的簡(jiǎn)諧作用強(qiáng)度及晶格常數(shù)．第i個(gè)粒子的質(zhì)量為：

固定Mmax＝20，Mmin＝1，Mmax為最左端晶格的質(zhì)量，Mmin為最右端晶格的質(zhì)量．基底勢(shì)強(qiáng)度為V＝5．

模擬計(jì)算中，采用固定邊界條件，并使兩端的晶格與NOSE－HOOVER 熱浴耦合，保持溫度分別為T(mén)L和TR．為簡(jiǎn)單起見(jiàn)，定義系統(tǒng)的平均溫度為T(mén)0＝（TL＋TR）／2，晶格鏈兩端的溫度差為ΔT＝TL－TR＝2T0Δ．當(dāng)Δ＞0時(shí)，最左端晶格與溫度較高的熱浴耦合，最右端的晶格與溫度較低的熱浴耦合；當(dāng)Δ＜0時(shí)，情況則相反．系統(tǒng)的局域化溫度定義為格點(diǎn)i處晶格的平均動(dòng)能：Ti＝m〈˙xi2〉．系統(tǒng)的局域熱流定義為Ji＝k〈˙xi（xi－xi－1）〉，（〈〉均表示對(duì)時(shí)間的平均）．計(jì)算中，使用改進(jìn)后的龍格庫(kù)塔算法對(duì)晶格的運(yùn)動(dòng)方程進(jìn)行積分，設(shè)置步長(zhǎng)為0．1，并有足夠長(zhǎng)的計(jì)算時(shí)間（107）使系統(tǒng)最終達(dá)到非平衡穩(wěn)態(tài)．當(dāng)系統(tǒng)的熱傳導(dǎo)過(guò)程達(dá)到非平衡穩(wěn)態(tài)以后，系統(tǒng)中流過(guò)不同格點(diǎn)的局域熱流量相同且為一恒定值．我們定義從左到右流過(guò)系統(tǒng)的熱流量的絕對(duì)值為J＋，當(dāng)熱浴反轉(zhuǎn)，從右到左的熱流量的絕對(duì)值則為J－，r＝J＋／J－ 為系統(tǒng)的熱整流效率．

2 結(jié)果與討論

2．1 不對(duì)稱(chēng)熱傳導(dǎo)

本文模型中，當(dāng)界面耦合強(qiáng)度為kint＝1．0時(shí)，系統(tǒng)為一完整的FK 鏈．線性質(zhì)量梯度下的FK 鏈已被證明具有熱整流效應(yīng)［17］．因此，我們先研究界面耦合強(qiáng)度較小時(shí)系統(tǒng)的熱傳導(dǎo)性質(zhì)．圖2為當(dāng)N＝100，kint＝0．05時(shí)，系統(tǒng)在不同平均溫度下的溫度分布圖．由圖2 可見(jiàn)，系統(tǒng)的溫度分布是不對(duì)稱(chēng)的，并且由于界面耦合彈簧的存在，所有情況中，系統(tǒng)都存在由界面熱阻（ITR）引起的溫度躍變．而且熱浴反轉(zhuǎn)后，界面熱阻是不對(duì)稱(chēng)的，正是因?yàn)榻缑鏌嶙璧牟粚?duì)稱(chēng)性導(dǎo)致了系統(tǒng)的不對(duì)稱(chēng)熱傳導(dǎo)現(xiàn)象．此外，該現(xiàn)象還與系統(tǒng)的平均溫度有關(guān)，溫度過(guò)高或過(guò)低都不利于實(shí)現(xiàn)不對(duì)稱(chēng)熱傳導(dǎo)，因此我們固定系統(tǒng)的平均溫度為T(mén)0＝0．1．

圖2 不同平均溫度下系統(tǒng)的歸一化溫度分布圖Fig．2 Normalized temperature profile of the system for different average temperature

2．2 熱整流反轉(zhuǎn)

質(zhì)量梯度下的FK 鏈，沿晶格質(zhì)量增加方向的熱流量要大于沿質(zhì)量減小方向的熱流量［17］．然而，當(dāng)FK 鏈中引入一彈性系數(shù)較小的界面彈簧后，這種單向的不對(duì)稱(chēng)熱傳導(dǎo)現(xiàn)象將被打破．圖3 為當(dāng)kint＝0．05時(shí)，3種不同體系大小下熱流量與熱浴溫差因子的關(guān)系圖．由圖3（a）（N＝100）可以看到，沿晶格質(zhì)量減小方向的熱流量J＋隨｜Δ｜的增加幾乎呈線性增長(zhǎng)，而沿質(zhì)量增加方向的熱流在時(shí)達(dá)到最大值后隨著｜Δ｜的增加而減小，導(dǎo)致了時(shí)J＋＜J－到J＋＞J－的熱整流反轉(zhuǎn)現(xiàn)象．因此系統(tǒng)僅在條件下，沿著晶格質(zhì)量增加方向的熱流量大于沿質(zhì)量減小方向的熱流量，而當(dāng)｜Δ｜＞0．4時(shí)，則轉(zhuǎn)為沿質(zhì)量減小方向的熱流量大于沿質(zhì)量增加方向的熱流量．可見(jiàn)，熱整流方向的改變僅僅是通過(guò)控制熱浴兩端溫差實(shí)現(xiàn)的．一般而言，熱流量會(huì)隨著熱浴溫差的增加而增大，然而在圖3（a）中，當(dāng)｜Δ｜＞0．5時(shí)出現(xiàn)了負(fù)導(dǎo)熱現(xiàn)象，正因?yàn)樗某霈F(xiàn)直接導(dǎo)致了熱整流的反轉(zhuǎn)．并且負(fù)導(dǎo)熱現(xiàn)象出現(xiàn)的區(qū)域隨著體系的增大逐步移向更高的熱浴溫差區(qū)．因此，當(dāng)體系大小為N＝1 000時(shí)，從J＋＜J－到J＋＞J－的熱整流的反轉(zhuǎn)現(xiàn)象出現(xiàn)在｜Δ｜≈0．7，如圖3（b）所示．隨著體系進(jìn)一步增大，負(fù)導(dǎo)熱現(xiàn)象最終消失，熱整流反轉(zhuǎn)也隨之消失，因此，當(dāng)N＝2 000時(shí)，系統(tǒng)中沿著晶格質(zhì)量增加方向的熱流量大于沿質(zhì)量減小方向的熱流量，即J＋＜J－ ．為了進(jìn)一步研究界面耦合強(qiáng)度及體系大小對(duì)熱傳導(dǎo)的影響，圖4給出了系統(tǒng)在以上3種體系下熱流量J±與界面耦合強(qiáng)度kint的關(guān)系，固定熱浴溫差因子．所有的情況下，隨著界面耦合強(qiáng)度的加強(qiáng)，J±迅速增加，并且在一定范圍內(nèi)呈現(xiàn)J±～k2int的關(guān)系，這種關(guān)系普遍存在于聲子通過(guò)一較弱耦合界面的情況中．隨著kint進(jìn)一步增加，J±將逐漸達(dá)到飽和，二者的關(guān)系也趨近穩(wěn)定．特別是，在圖4（a）和（b）中還存在J＋＞J－到J＋＜J－的交叉現(xiàn)象，即當(dāng)界面耦合強(qiáng)度較小時(shí)，系統(tǒng)沿晶格質(zhì)量減小方向的熱流量大于沿質(zhì)量增加方向的熱流量，隨著耦合強(qiáng)度的增加，變?yōu)檠刭|(zhì)量增加方向的熱流量大于沿質(zhì)量減小方向的熱流量．而當(dāng)體系較大時(shí)（圖4（c）），系統(tǒng)的熱整流方向不受界面耦合強(qiáng)度的影響，沿晶格質(zhì)量增加方向的熱流量始終大于沿質(zhì)量減小方向的熱流量，即J＋＜J－．可見(jiàn)，無(wú)論體系多大，當(dāng)界面彈簧彈性系數(shù)較大時(shí)，均有J＋＜J－，這是因?yàn)楫?dāng)系統(tǒng)處于強(qiáng)耦合狀態(tài)時(shí)，兩段FK 鏈組成一個(gè)整體，成為一條完整的FK 鏈，與文獻(xiàn)［17］中得出的熱整流方向一致，且系統(tǒng)的熱整流方向不隨體系的改變而改變．而當(dāng)系統(tǒng)處于弱耦合狀態(tài)時(shí)，系統(tǒng)的熱整流方向會(huì)隨著熱浴溫差，體系大小的改變而改變，出現(xiàn)熱整流反轉(zhuǎn)現(xiàn)象．

2．3 熱整流機(jī)理

針對(duì)以上不對(duì)稱(chēng)熱傳導(dǎo)現(xiàn)象，我們根據(jù)有效聲子理論研究了其熱整流機(jī)理：當(dāng)FK 鏈中的晶格處于高溫極限時(shí)，其動(dòng)能足以克服基底勢(shì)能的束縛，因此基底勢(shì)可忽略不計(jì)并且系統(tǒng)可近似地看做為簡(jiǎn)諧晶格鏈，晶格中聲子的振動(dòng)頻率范圍為0＜ω／2π＜；當(dāng)晶格處于低溫極限時(shí)則被困在基底勢(shì)壘中，通過(guò)線性化運(yùn)動(dòng)方程可得到晶格中聲子的振動(dòng)頻率范圍為

因此，當(dāng)系統(tǒng)處于強(qiáng)耦合狀態(tài)時(shí)，系統(tǒng)回歸為一條完整的FK 鏈，不存在界面熱阻，我們主要研究熱浴兩端晶格的聲子頻率：當(dāng)系統(tǒng)左端質(zhì)量較大的晶格與高溫?zé)嵩〗佑|時(shí)，晶格處于高溫極限中，聲子的振動(dòng)頻率為0＜ω／2π＜0．071．而與低溫?zé)嵩〗佑|的右端質(zhì)量較小的晶格中，聲子的頻率為0．35＜ω／2π＜0．47．可見(jiàn)左右兩端晶格中聲子的頻率不重合，不能匹配，這就意味著聲子不能較順利地從左往右輸運(yùn)熱量，因此系統(tǒng)中的熱流量較?。?/p>

圖3 不同體系下熱流量J 與｜Δ｜的關(guān)系圖Fig．3 Dependence of the heat flux Jon｜Δ｜for different system size

圖4 不同體系下熱流量J 與kint的關(guān)系圖Fig．4 Dependence of the heat flux J on kintfor different system size

當(dāng)系統(tǒng)左端晶格與低溫?zé)嵩〗佑|時(shí)，其聲子頻率為0．079＜ω／2π＜0．106．而右端晶格處于高溫極限，其聲子頻率為0＜ω／2π＜0．308．可見(jiàn)左右兩端晶格中聲子的頻率能較好地匹配，更有利于聲子的熱輸運(yùn)，因此，J＋＜J－，這與我們的結(jié)果吻合．然而，當(dāng)系統(tǒng)處于弱耦合狀態(tài)時(shí)，以上結(jié)論無(wú)法解釋熱整流的反轉(zhuǎn)．但可以很直觀地看出，負(fù)導(dǎo)熱現(xiàn)象直接導(dǎo)致了熱整流反轉(zhuǎn)，而負(fù)導(dǎo)熱現(xiàn)象被認(rèn)為是由負(fù)微分熱阻所引起，因此需要考慮界面效應(yīng)．正是由于彈性系數(shù)較小的界面彈簧存在，當(dāng)系統(tǒng)右端質(zhì)量較小粒子與低溫?zé)嵩〗佑|時(shí)，導(dǎo)致系統(tǒng)出現(xiàn)了負(fù)微分熱阻，即隨著系統(tǒng)右端的溫度進(jìn)一步減小，界面兩端晶格中聲子頻率呈現(xiàn)由匹配到分離的現(xiàn)象．于是J－隨著熱浴溫差的增加而減小，熱整流方向由最初的J＋＜J－變?yōu)镴＋＞J－，因此熱整流反轉(zhuǎn)了．由此可見(jiàn)，只有界面耦合強(qiáng)度較小情況下，才能通過(guò)調(diào)節(jié)熱浴溫差控制熱整流的方向．

此外，熱整流反轉(zhuǎn)除了與耦合強(qiáng)度有關(guān)外，還與體系的大小密切相關(guān)．由圖3和圖4可知，只有在體系較小的條件下才能出現(xiàn)熱整流反轉(zhuǎn)，當(dāng)體系大于一定的范圍后，無(wú)論界面耦合強(qiáng)度如何，系統(tǒng)都呈現(xiàn)出和強(qiáng)耦合狀態(tài)下一樣的熱整流性質(zhì)，沒(méi)有熱整流反轉(zhuǎn)現(xiàn)象．因此，圖3（c）中熱整流反轉(zhuǎn)的消失可理解如下：界面彈簧對(duì)系統(tǒng)的影響隨著體系的增大而減小，在較大的體系中，界面彈簧對(duì)系統(tǒng)的作用可以忽略不計(jì)，系統(tǒng)的兩部分仍可被看做強(qiáng)耦合在一起的整體．因此，可以總結(jié)得出，熱整流反轉(zhuǎn)的根本原因是界面彈簧的引入，并且只有在特定的體系下，當(dāng)系統(tǒng)為弱耦合狀態(tài)時(shí)，才能通過(guò)調(diào)節(jié)熱浴溫差控制熱整流方向．

2．4 熱整流效率

值得一提的是，對(duì)于特定的體系（如N＝100），系統(tǒng)的熱整流效率r＝J＋／J－ 可從r→0．94變?yōu)閞→2．1，與目前實(shí)驗(yàn)中制備的固體熱二極管的熱整流效率較為一致．這意味著，該模型在熱整流方向控制領(lǐng)域具有潛在的應(yīng)用價(jià)值．此外，由于負(fù)微分熱阻的存在以及功能梯度材料的廣泛應(yīng)用，該模型也是構(gòu)建熱晶體管，熱邏輯門(mén)等熱學(xué)器件不錯(cuò)的候選者．

3 結(jié) 論

本文研究了質(zhì)量梯度下一維非線性功能梯度材料的熱整流性質(zhì)．該系統(tǒng)在特定的界面耦合強(qiáng)度及特定的體系大小范圍內(nèi)，可以?xún)H通過(guò)控制熱浴溫差實(shí)現(xiàn)由J＋＜J－到J＋＞J－的熱整流反轉(zhuǎn)．引起熱整流反轉(zhuǎn)的根本原因是系統(tǒng)中存在的界面彈簧．鑒于該系統(tǒng)無(wú)需進(jìn)行結(jié)構(gòu)重建便能在原地進(jìn)行熱整流方向的控制，討論了它在實(shí)驗(yàn)制備及應(yīng)用方面的可能性．

［1］LEPRI S，LIVI R，POLITI A．Heat conduction in chains of nonlinear oscillators［J］．Physical Review Letters，1997，78（10）：1896－1899．

［2］HU Bam－bi，LI Bao－wen，ZHAO Hong．Heat conduction in one－dimensional chains［J］．Physical Review E，1998，57（3）：2992－2995．

［3］LEPRI S，LIVI R，POLITI A．Thermal conduction in classical low dimensional lattices［J］．Physics Reports－Review Section of Physics Letters，2003，377（1）：1－80．

［4］PROSEN T，CAMPBELL D K．Normal and anomalous heat transport in one－dimensional classical lattices ［J］．Chaos，2005，15（1）：015117．

［5］GAUL C，BUETTNER H．Quantum mechanical heat transport in disordered harmonic chains［J］．Physical Review E，2007，76（1）：011111．

［6］LI Nian－bei，REN Jie，WANG Lei，etal．Colloquium：Phononics：Manipulating heat flow with electronic analogs and beyond［J］．Reviews of Modern Physics，2012，84（3）：1045－1066．

［7］TERRANEO M，PEYRARD M，CASATI G．Controlling the energy flow in nonlinear lattices：A model for a thermal rectifier［J］．Physical Review Letters，2002，88（9）：094302．

［8］LI Bao－wen，WANG Lei，CASATI G．Negative differential thermal resistance and thermal transistor［J］．Applied Physics Letters，2006，88（14）：143501．

［9］WANG Lei，LI Bao－wen．Thermal logic gates：Computation with phonons［J］．Physical Review Letters，2007，99（17）：177208．

［10］LI Bao－wen，WANG Lei，CASATI G．Thermal diode：Rectification of heat flux［J］．Physical Review Letters，2004，93（18）：184301．

［11］YANG Nuo，ZHANG Gang，LI Bao－wen．Thermal rectification in asymmetric graphene ribbons［J］．Applied Physics Letters，2009，95（3）：033107．

［12］YANG Nuo，ZHANG Gang，LI Bao－wen．Carbon nanocone：A promising thermal rectifier［J］．Applied Physics Letters，2008，93（24）：243111．

［13］CHANG C W，OKAWA D，MAJUMDAR A，etal．Solidstate thermal rectifier［J］．Science，2006，314（5802）：1121－1124．

［14］KOBAYASHI W，TERAOKA Y，TERASAKI I，An oxide thermal rectifier［J］．Applied Physics Letters，2009，95（17）：171905．

［15］TIAN He，XIE Dan，YANG Yi，etal．A novel solid－state thermal rectifier based on reduced graphene oxide［J］．Scientific Reports，2012，2：523．

［16］YANG Nuo，LI Nian－bei，WANG Lei，etal．Thermal rectification and negative differential thermal resistance in lattices with mass gradient［J］．Physical Review B，2007，76（2）：020301．

［17］HU Tao，HU Ke，TANG Yi．Anharmonic on－site potential implies asymmetry thermal transport in 1D mass－graded harmonic lattice［J］．Physica B－Condensed Matter，2010，405（21）：4407－4412．

［18］SHAH TEJAL N，GAJJAR P N．Study of thermal conductivity and thermal rectification in exponential mass graded lattices［J］．Physics Letters A，2012，376（8）：438－441．