王彥成 邱吳劼 楊宏亮 席麗麗 楊炯 張文清

1)(中國(guó)科學(xué)院上海硅酸鹽研究所,高性能陶瓷和超微結(jié)構(gòu)國(guó)家重點(diǎn)實(shí)驗(yàn)室,上海 200050)

2)(中國(guó)科學(xué)院大學(xué),北京 100049)

3)(上海大學(xué)材料基因組工程研究院,上海 200444)

1 引 言

熱電材料需要優(yōu)異的電輸運(yùn)性質(zhì)和低熱導(dǎo)率.理解低晶格熱導(dǎo)率的形成機(jī)理,對(duì)于提高熱電材料的性能非常重要[1].方鈷礦材料具有良好的電輸運(yùn)性能,而且其內(nèi)部特殊的籠狀結(jié)構(gòu)能夠容納外來原子(填充原子),進(jìn)而降低材料的晶格熱導(dǎo)率.因此,填充方鈷礦體系是一類在中溫區(qū)熱電性能十分優(yōu)異的材料.再考慮到其不含鉛、碲等有毒元素的特點(diǎn),在過去的十多年里,方鈷礦材料受到了研究人員的廣泛關(guān)注和研究[2?4],其性能也得到了大幅提升,尤其是以三元填充化合物為代表[5?7],熱電優(yōu)值在800—850 K范圍達(dá)到了1.7—1.9.

然而在填充原子降低晶格熱導(dǎo)率的機(jī)理上,理論研究卻長(zhǎng)期存在著爭(zhēng)論.傳統(tǒng)的觀點(diǎn)認(rèn)為,填充原子與框架間較弱的化學(xué)鍵使其做類似于Einstein模式的大幅度獨(dú)立振動(dòng),并通過共振散射低頻的晶格聲子來降低晶格熱導(dǎo)率[8,9].一般來說,局域振動(dòng)模式對(duì)材料自身性質(zhì)所產(chǎn)生的影響,許多能夠在實(shí)驗(yàn)上測(cè)得并在理論計(jì)算中得到一定程度的驗(yàn)證[10].填充原子的這種特殊的振動(dòng)行為則在非彈性中子散射和熱容測(cè)量的實(shí)驗(yàn)研究中得到了確認(rèn)[11?13].但另一方面,第一性原理計(jì)算表明填充原子與晶格框架存在著非簡(jiǎn)諧的耦合[14],(La,Ce)Fe4Sb12的中子散射實(shí)驗(yàn)結(jié)果也支持這一觀點(diǎn),并認(rèn)為填充原子與框架之間的作用似乎可以用聲子與聲子作用來描述,即將填充原子的振動(dòng)描述為具有長(zhǎng)程相干性和周期性特征的聲子行為[15].Li和Mingo[16]用第一性原理方法在三聲子框架下計(jì)算了全填體系YbFe4Sb12的聲子輸運(yùn)性質(zhì),認(rèn)為填充方鈷礦和其他常規(guī)材料一樣,填充原子的引入增強(qiáng)了三聲子的散射截面,因此降低了晶格熱導(dǎo)率,但其晶格熱導(dǎo)率的計(jì)算結(jié)果卻遠(yuǎn)低于實(shí)驗(yàn)值.以上的研究結(jié)果都表明了填充原子與晶格之間相互作用的復(fù)雜性,而這種復(fù)雜性與填充原子偏離平衡位置的大幅度振動(dòng)行為密切相關(guān).最近,在具有化學(xué)鍵層級(jí)特征的Cu3SbSe3化合物的研究[17,18]中發(fā)現(xiàn),與周圍環(huán)境的成鍵較弱的部分,Cu原子會(huì)出現(xiàn)類似于非晶材料中協(xié)同特征明顯的低頻振動(dòng),其理論描述偏離了以晶態(tài)為基礎(chǔ)的“小參數(shù)”微擾近似.雖然在填充方鈷礦體系中填充原子處在封閉的籠子中,但是其同樣具有弱化學(xué)鍵下的動(dòng)態(tài)漲落的亞晶格特征.因此,我們初步認(rèn)為填充原子的振動(dòng)行為偏離了微擾近似的圖像,相應(yīng)地散射機(jī)理也可能同時(shí)存在非微擾的特征[18].本文要解決的關(guān)鍵問題是填充原子的運(yùn)動(dòng)是否偏離了微擾近似的范疇,并是否具有超出三聲子過程的散射機(jī)理,以及極低晶格熱導(dǎo)率的形成原因.

在半導(dǎo)體和絕緣體中,材料本征的晶格熱導(dǎo)率取決于其內(nèi)部的聲子及其散射過程,可由聲子的玻爾茲曼輸運(yùn)方程(Boltzmann transport equation)描述.基于零溫的第一性原理方法計(jì)算并預(yù)測(cè)材料晶格熱導(dǎo)率的準(zhǔn)確程度,受到原子間作用力的計(jì)算和非彈性聲子散射機(jī)理選取的雙重影響[19].不同于利用凍聲子法[20,21]得到簡(jiǎn)諧及非諧力常數(shù)的常規(guī)手段,本文采用考慮了有限溫度效應(yīng)的第一性原理分子動(dòng)力學(xué)(ab initiomolecular dynamics,AIMD)的方法,覆蓋了填充原子全部的相互作用,即充分考慮“全”相互作用的力常數(shù).本文以YbFe4Sb12為例模擬填充原子這種大幅偏離平衡位置的振動(dòng)行為,并捕捉得到其更為接近真實(shí)情況的軌跡與受力信息.AIMD的結(jié)果結(jié)合Hellman等[22,23]最近發(fā)展的溫度相關(guān)有效勢(shì)場(chǎng)(temperature-dependent effective potential,TDEP)方法,統(tǒng)計(jì)其動(dòng)力學(xué)軌跡和相應(yīng)的原子受力信息,擬合得到有效的二階、三階力常數(shù).這實(shí)質(zhì)上是將填充原子的作用等效為非諧作用進(jìn)行考慮.在聲子非彈性散射方面,仍以三聲子過程為主要散射機(jī)理,計(jì)算聲子自能與壽命,以及晶格熱導(dǎo)率.通常情況下,這種做法在非諧作用較強(qiáng)的材料中符合較好,但是對(duì)于可能存在的非微擾的散射機(jī)理卻難以給予恰當(dāng)?shù)倪^程描述和機(jī)理解釋.從反證的角度來講,將計(jì)算得到的晶格熱導(dǎo)率和實(shí)驗(yàn)結(jié)果相比較,進(jìn)而檢驗(yàn)填充原子引入的散射機(jī)理是否具有非微擾的特征,則可以為解釋相關(guān)體系熱導(dǎo)率降低的機(jī)理提供清晰的線索.因此,為了突出填充原子在熱輸運(yùn)微觀過程的影響,選用Yb原子對(duì)方鈷礦進(jìn)行填充.由于滿填體系具有最好的空間對(duì)稱性和最少的計(jì)算量,且CoSb3存在填充量上限我們選用YbFe4Sb12進(jìn)行研究.

2 計(jì)算方法

非諧作用是聲子與聲子碰撞散射的主要原因,它會(huì)使聲子的頻率發(fā)生變化,而且峰也會(huì)產(chǎn)生一定的展寬,后者意味著聲子存在著有限的壽命.晶體的勢(shì)能展開形式可以寫為

其中,Φ0為勢(shì)能的常數(shù)項(xiàng);ui分別為原子i的位移;αβγ為笛卡爾分量的指標(biāo);表示二階與三階力常數(shù);(1)式中右邊的第二、三項(xiàng)分別為簡(jiǎn)諧勢(shì)與立方非諧勢(shì),記為U2與U3.其中立方非諧勢(shì)的二次量子化形式可以寫作三聲子碰撞的形式,

這里λ=(q,j),?λ=(?q,j),其中q,j分別為聲子模式的波矢和所在格波數(shù),?a?λ和?aλ為聲子模式λ的增加和湮沒算符.Φλλ′λ′′為參與散射的三個(gè)聲子λ,λ′與λ′′間的作用強(qiáng)度,其表達(dá)式[24,25]為

這里的?(q+q′+q′′)表達(dá)的是準(zhǔn)動(dòng)量守恒條件,即當(dāng)q+q′+q′′為0或者等于倒格矢G的整數(shù)倍時(shí)取1,否則取0.在微擾理論下,考慮三聲子過程為一階近似,利用Fermi’s golden rule[24]得到聲子模式λ的自能虛部項(xiàng)Γλ(ω),

其中,nλ為平衡態(tài)時(shí)聲子的占據(jù)數(shù).利用聲子壽命和聲子展寬2Γλ(ω)的關(guān)系,τλ(ω)=1/2Γλ(ω),可以計(jì)算得到晶格熱導(dǎo)率其中,Cλ為摩爾熱容.

采用基于密度泛函理論的平面波基矢軟件VASP軟件包進(jìn)行第一性原理分子動(dòng)力學(xué)的模擬,采用綴加平面波的方法,交換關(guān)聯(lián)勢(shì)采取廣義梯度近似(GGA)下的Perdew-Burker-Ernzerhof形式,平面波截?cái)嗄転?00 eV,在正則系綜下模擬了YbFe4Sb12的2×2×2(272個(gè)原子)超胞在300,400,500和600 K等溫度下的動(dòng)力學(xué)過程,其中溫度控制采取了Nosé thermostat的熱浴,布里淵區(qū)積分僅限于Γ點(diǎn),在系統(tǒng)達(dá)到平衡狀態(tài)后模擬時(shí)間在30 ps以上,模擬時(shí)間步長(zhǎng)為2 fs.為考慮熱膨脹的因素,不同溫度下YbFe4Sb12的晶胞常數(shù)采用了quasi-harmonic approximation(QHA)方法進(jìn)行估算,所用的軟件為phonopy[20].使用TDEP軟件[22,23]對(duì)AIMD中原子的軌跡與受力進(jìn)行最小二乘法擬合,進(jìn)而得到二階與三階力常數(shù).該方法中考慮了平移不變性和晶格自身的對(duì)稱性,以減少力常數(shù)矩陣中所需擬合的獨(dú)立矩陣元數(shù)目,進(jìn)而減少了計(jì)算量.在計(jì)算中,取二階力常數(shù)的截?cái)喟霃酱笾聻?0 ?,三階力常數(shù)的截?cái)喟霃綖?.0 ?.這是由于通常情況下三階力常數(shù)比較小[26],只需考慮較為鄰近的三原子對(duì)的受力情況即可.在二階和三階力常數(shù)的基礎(chǔ)上,采用11×11×11的q點(diǎn)網(wǎng)格來計(jì)算YbFe4Sb12的聲子壽命以及晶格熱導(dǎo)率等熱輸運(yùn)性質(zhì).

3 結(jié)果與討論

所用AIMD模擬結(jié)合TDEP的方法的特點(diǎn),能夠從更趨實(shí)際的原子運(yùn)動(dòng)軌跡中統(tǒng)計(jì)分析出材料的力常數(shù)等內(nèi)在性質(zhì).高溫下材料內(nèi)部原子大幅振動(dòng),并帶來顯著的非諧效應(yīng)等,是這一方法的長(zhǎng)處所在[23].而在低溫或者材料內(nèi)原子位移參數(shù)(atomic displacement parameters,ADP)比較小的情況下,該方法和凍聲子法的圖像趨向一致.

方鈷礦材料最顯著的特點(diǎn)是框架部分與填充原子在振動(dòng)幅度存在著較大的差別,框架原子主要是圍繞著平衡位置附近做較小的振動(dòng),而填充原子的振動(dòng)幅度是前者的3—5倍之多.在討論填充體系YbFe4Sb12之前,用AIMD結(jié)合TDEP方法計(jì)算純方鈷礦Co4Sb12的熱輸運(yùn)性質(zhì),來檢驗(yàn)這一方法的有效性與準(zhǔn)確性.

3.1 純方鈷礦材料Co4Sb12的聲子輸運(yùn)性質(zhì)

圖1 使用TDEP方法計(jì)算的Co4Sb12(a)在300 K的聲子壽命隨頻率變化的關(guān)系,(b)晶格熱導(dǎo)率 κL的溫度依賴關(guān)系(紅線),其中低溫和高溫實(shí)驗(yàn)晶格熱導(dǎo)率分別來自Nolas等[29]和Chen等[28]的研究,作為對(duì)照的凍聲子法計(jì)算的κL(藍(lán)色虛線)來自于Guo等[30]的研究.Fig.1.(a)Phonon lifetimes of Co4Sb12at 300 K vs.frequency and(b)its lattice thermal conductivity vs.temperature,calculated by TDEP(red solid line).The experimental data of κLat low and high temperature is from references of Nolas et al.[29]and Chen et al.[28]respectively.The referent calculated κLby fi nite-displacement approach(blue dash line)is from reference of Guo et al.[30].

用AIMD模擬300 K時(shí)Co4Sb12的動(dòng)力學(xué)行為,并用TDEP擬合得到其二階與三階力常數(shù),計(jì)算得到了聲子壽命和晶格熱導(dǎo)率的結(jié)果,如圖1所示.從圖1(a)中可以看出Co4Sb12的聲子壽命τλ(ω)正比于ω?2.通常來說,非諧作用所帶來的聲子非彈性散射機(jī)理主要為三聲子U過程,其聲子壽命的表達(dá)式具有τU∝ω?2的關(guān)系[27].這一結(jié)果表明了在純方鈷礦Co4Sb12中,三聲子過程是其主要的散射機(jī)理.因此,在三聲子框架下計(jì)算得到的晶格熱導(dǎo)率同實(shí)驗(yàn)結(jié)果[28,29]符合得比較好,而且該結(jié)果與Guo等[30]用凍聲子方法計(jì)算的結(jié)果(GGA potential的結(jié)果)幾乎一致,見圖1(b).這說明用AIMD模擬結(jié)合TDEP方法計(jì)算材料熱輸運(yùn)性質(zhì),確實(shí)是有效可行的.

3.2 全填體系YbFe4Sb12的聲子輸運(yùn)性質(zhì)

不同于純方鈷礦材料Co4Sb12僅用三聲子相互作用就能計(jì)算得到較為準(zhǔn)確的熱導(dǎo)率結(jié)果,填充原子的引入使得體系的熱導(dǎo)計(jì)算變得復(fù)雜起來.為研究填充原子對(duì)于方鈷礦材料熱輸運(yùn)性質(zhì)的影響,計(jì)算了YbFe4Sb12在300 K時(shí)的聲子壽命的頻率變化關(guān)系,并分析了填充原子振動(dòng)模式的局域性特征,以及這些模式在壽命上所具有的特點(diǎn).從圖2(b)中可以看到,相比于純的體系Co4Sb12,YbFe4Sb12的聲子壽命整體降低了許多,這表明后者的非諧作用比純的體系增強(qiáng)了很多.其中低頻的聲學(xué)支模式和部分光學(xué)支模式的壽命滿足ω?2的關(guān)系,但是在填充原子振動(dòng)頻率附近的模式,其壽命相比τU∝ω?2趨勢(shì)有著顯著的降低.考慮到TDEP方法仍是在微擾近似的前提下,以三聲子過程為主要框架來計(jì)算聲子輸運(yùn)性質(zhì)的這一特點(diǎn),我們認(rèn)為這是由于填充原子的引入帶來了非常規(guī)的聲子散射機(jī)理所致.

圖2 YbFe4Sb12的聲子色散關(guān)系(a),(c)和聲子壽命的頻率變化關(guān)系(b),以及不同頻率的原子參與率APR分布(d),其中(a)和(b)中的藍(lán)色部分為聲子參與率PR＜0.2部分,紅色部分為PR＞0.2部分,(b)中粉紅色點(diǎn)為Co4Sb12的聲子壽命,(b)中虛線為擬合得到的各自顏色部分的頻率變化趨勢(shì),均為ω?2關(guān)系,(c)中藍(lán)色部分為填充原子APR(Yb)＞0.1的部分Fig.2.Calculated phonon dispersions of YbFe4Sb12(left),its phonon lifetimes(b)and atomic participation ratio for phonon modes of YbFe4Sb12(d).Blue parts of top panels(a),(b)mean modes with low PR(PR＜0.2).In panel(b)pink dots are phonon lifetimes of Co4Sb12,while red(PR＞0.2)and blue dots(PR＜0.2)are lifetimes of YbFe4Sb12,and dash lines are fi tted ω?2trends of lifetimes for U process.In panel(c),blue part means phonon modes with high APR of Yb,which is larger than 0.1.

為分析這部分振動(dòng)模式的性質(zhì)及其與Yb的關(guān)系,計(jì)算了聲子參與率p(ωλ)(participation ratio,PR)和原子參與率pi(ωλ)(atomic participation ratio,APR),它們分別描述了所有原子在聲子模式λ中的參與程度[31]與不同原子i對(duì)于聲子模式λ的貢獻(xiàn)[32],其公式如下:

圖2(b)中偏離τU∝ω?2關(guān)系的振動(dòng)模式,其PR值恰好小于0.2,可知這部分模式為局域振動(dòng)模式.這些低頻局域模式的波矢分布可參見圖2(a)中 0—2 THz范圍內(nèi)的藍(lán)色部分.而原子參與率APR的分析能夠幫助我們找到這些局域模式中的主要參與原子.圖2(d)中0.8—1.4 THz范圍內(nèi)Yb的參與率要遠(yuǎn)高于框架原子Fe和Sb,這一頻率范圍又和圖2(b)中的低頻局域振動(dòng)模重合,其波矢分布(圖2(c)中藍(lán)色部分)與圖2(a)中局域振動(dòng)模式的低頻部分一致.這些表明填充原子Yb的振動(dòng)屬于分布在較大的波矢區(qū)間內(nèi)的低頻局域振動(dòng)模.因此可以認(rèn)為聲子壽命偏離ω?2關(guān)系所對(duì)應(yīng)的散射機(jī)理,實(shí)際是由填充原子Yb主導(dǎo)的局域振動(dòng)模式參與的聲子散射過程.這一物理過程的實(shí)質(zhì)是在材料內(nèi)部傳播的晶格聲子,借助填充原子與局域振動(dòng)模發(fā)生不停的轉(zhuǎn)換,而這種振動(dòng)特征的變換,正是除去較強(qiáng)非諧作用之外填充方鈷礦材料熱阻的主要來源.這種相似的結(jié)論在Pailhès等[32]的研究中也有所體現(xiàn).因此,鑒于填充原子Yb的局域振動(dòng)特點(diǎn),其對(duì)于晶格框架的影響原則上可以近似認(rèn)為具有累加的特點(diǎn),即當(dāng)Yb的填充量減小時(shí),填充體系的非諧作用以及非常規(guī)的聲子散射均會(huì)減弱.實(shí)驗(yàn)上填充體系晶格熱導(dǎo)率隨填充量增加而降低[34],也印證了這一點(diǎn).

從以上討論中,我們對(duì)填充原子帶來的散射機(jī)理有了較為初步的整體了解.接下來從細(xì)節(jié)的角度來驗(yàn)證這一論斷,并分析填充原子在低熱導(dǎo)率成因中所起的具體作用.通常來講,非諧作用對(duì)于材料聲子譜的影響,主要表現(xiàn)在聲子的峰位和半高寬,峰位的變化即為頻率的位移,展寬則表明存在著聲子的阻尼,意味著其有限的聲子壽命.這分別對(duì)應(yīng)著聲子自能Σ(?)=?(?)+iΓ(?)中的實(shí)部?與虛部Γ,其中?=E/?,E為中子的探測(cè)能量,這些都可以在聲子的中子散射截面σλ(?)中得到比較直觀的反映.其公式[35]寫作

其中ωλ為聲子模式λ的本征頻率.

圖3(a)為Γ-P方向中子散射截面的計(jì)算結(jié)果.聲子模式振動(dòng)特征的分析,顯示出YbFe4Sb12的低頻聲子具有類似于clathrate的avoid crossing[36]的特征,如圖3(a)中標(biāo)識(shí)的橫光學(xué)支TO/TO′,與縱聲學(xué)支LA/LA′,在其相互臨近的區(qū)域(Q3點(diǎn)附近)產(chǎn)生作用,然后彼此分開(即為TO′與LA′部分).可以看出TO/TO′的展寬都比較大,相應(yīng)地壽命都比較小.從圖2的結(jié)論中,我們了解到這部分振動(dòng)模式中Yb的貢獻(xiàn)很大,通過分析振動(dòng)模式中Yb原子的權(quán)重定量地得到了在avoid crossing前后部分TO/TO′,以及LA/LA′中填充原子Yb在聲子展寬上的具體貢獻(xiàn),如圖3(b)和圖3(c)所示.可以看到,Q1,Q2,Q4和Q5處TO/TO′中Yb所貢獻(xiàn)的展寬較為一致,且都很高,這說明Yb的振動(dòng)模式壽命比較低.而在Q3處,由于相互作用,TO中的Yb的振動(dòng)有一部分轉(zhuǎn)移到了LA中,提高了LA的展寬.如果將Q3處TO與LA展寬中Yb的貢獻(xiàn)量加起來,得到的Yb振動(dòng)模式的展寬同avoid crossing前后Q1,Q2,Q4和Q5處也保持相近的水平.這說明填充原子Yb的振動(dòng)模式自身展寬比較大,壽命很低,而與之發(fā)生相互作用的聲學(xué)支聲子由于混入了低弛豫時(shí)間的填充原子Yb的振動(dòng),其壽命也隨之降低.這種圖像與填充原子Yb局域振動(dòng)模參與的非常規(guī)的三聲子過程,令這部分模式的壽命比ω?2關(guān)系更低的結(jié)果是相近的,也能夠部分地解釋晶格熱導(dǎo)率降低的原因.

圖3 YbFe4Sb12在Γ-P方向的(a)中子散射截面σ的計(jì)算結(jié)果,其中顏色表示峰的強(qiáng)度;(b)橫光學(xué)支TO/TO′;(c)縱聲學(xué)支LA/LA′在avoid crossing前后部分的聲子自能虛部Γ;(b)和(c)中Q1—Q5對(duì)應(yīng)(a)中白色虛線標(biāo)識(shí)的5個(gè)波矢,藍(lán)色和橘色分別對(duì)應(yīng)Yb和框架原子Fe/Sb各自的貢獻(xiàn)量Fig.3. (a)Calculated scattering cross section of YbFe4Sb12along Γ-P direction.Intensities is plotted on a logarithmic scale.Imaginary self energy Γ (ω)of TO/TO′branch(b)and LA/LA′branch(c)at Q1–Q5 wave vector points around the avoid crossing part,which are labeled in panel(a).In panels(b)and(c),blue and orange parts mean contributions of Yb and Fe/Sb atoms to Γ (ω)respectively.

3.3 動(dòng)力學(xué)上填充原子的局域振動(dòng)特征

以上討論均是從晶格動(dòng)力學(xué)的角度出發(fā),以聲子以及聲子與聲子相互作用為前提,即將晶體材料中原子的振動(dòng)簡(jiǎn)化為聲子這一集體振動(dòng)的元激發(fā)行為.雖然前面的討論中,確定了填充原子Yb的振動(dòng)模式具有局域振動(dòng)的特征,但是填充原子之間的在振動(dòng)行為上的關(guān)聯(lián)性質(zhì),以及是否適用于聲子的圖像,仍舊沒有確切的結(jié)論.而這又是回答填充原子以何種方式與框架原子或晶格聲子作用的物理圖像上關(guān)鍵的基礎(chǔ)問題.

我們從動(dòng)力學(xué)的角度出發(fā),分析填充原子與周圍框架原子在微觀運(yùn)動(dòng)過程中的關(guān)聯(lián)性質(zhì),力求為理解其相互作用構(gòu)建一個(gè)清晰而準(zhǔn)確的物理圖像.研究中引入了原子速度的相干性這一物理量,它實(shí)際上是歸一化的原子速度的互關(guān)聯(lián)譜.可以將原子對(duì)中的兩個(gè)原子的速度理解為兩個(gè)或隨機(jī)或存在著某種程度線性關(guān)系的信號(hào),它們?cè)谔囟l率上的關(guān)聯(lián)性質(zhì)可以用互關(guān)聯(lián)譜

來表示.借助于和卷積相似的數(shù)學(xué)手段,除下原子i或j的速度自關(guān)聯(lián)函數(shù)的功率譜sii(ω)或sjj(ω),得到歸一化后的互關(guān)聯(lián)譜,

這里稱cij(ω)為原子i和j速度在頻率ω上的相干性;cij(ω)是一無(wú)量綱的數(shù),取值范圍為0≤cij(ω)≤ 1. 當(dāng)cij(ω)≈1.0時(shí),vi(t)和vj(t)是線性相關(guān)的;當(dāng)cij(ω)≈0時(shí),vi(t)和vj(t)彼此之間是相互獨(dú)立和隨機(jī)的關(guān)系.

圖4是分析YbFe4Sb12在600 K時(shí)AIMD的數(shù)據(jù)得到的填充原子Yb和框架原子Fe/Sb之間的速度關(guān)聯(lián)性質(zhì)的結(jié)果,其中原子對(duì)中原子的距離選為9.2 ?,即為晶胞Yb2Fe8Sb24的周期長(zhǎng)度.從YbFe4Sb12的聲子分波態(tài)密度的結(jié)果中可知填充原子Yb振動(dòng)模式的頻率范圍在1—1.3 THz,Sb原子的振動(dòng)頻率主要在1—5 THz,Fe的振動(dòng)頻率主要為6—7.5 THz.由圖4可以看到,在1—1.3 THz的特征振動(dòng)頻率范圍內(nèi),原子對(duì)Yb-Yb的相干性小于0.1,可知其關(guān)聯(lián)性非常弱.而框架原子Sb-Sb和Fe-Fe的在其特征振動(dòng)頻率范圍的相干性均在0.35左右,其關(guān)聯(lián)性比Yb要強(qiáng)很多.鑒于框架原子Fe/Sb這種具有關(guān)聯(lián)性的運(yùn)動(dòng),通常被解釋為具有周期性和集體振動(dòng)等元激發(fā)特征的晶格聲子,基本可以認(rèn)為填充原子Yb這種關(guān)聯(lián)性較弱且周期性不明顯的運(yùn)動(dòng)偏離了聲子的圖像,更傾向于彼此間相互獨(dú)立的局域振動(dòng).有了這一認(rèn)識(shí),再將填充原子和框架部分的晶格聲子間的作用解釋為三聲子過程則是不合適的.前面部分中用TDEP的三聲子框架來計(jì)算聲子輸運(yùn)性質(zhì),可以視作一種等效的處理,而等效的效果如何,則可以從其晶格熱導(dǎo)率的計(jì)算結(jié)果直觀地看出來.

圖4 在600 K時(shí)YbFe4Sb12中不同種類原子對(duì)(a)Yb-Yb,(b)Sb-Sb,(c)Fe-Fe的原子速度相關(guān)性.為方便查看,各自原子主要貢獻(xiàn)的振動(dòng)頻率范圍由紅色填充作為標(biāo)識(shí),且原子對(duì)中的原子編號(hào)在圖(d)中均有標(biāo)明,其彼此間的距離均為9.2 ?Fig.4.Coherence of velocities between pairs of atoms in YbFe4Sb12at 600 K:(a)Yb-Yb,(b)Sb-Sb,(c)Fe-Fe,and numbers behind atom name are index of these atoms in the supercell.As a guide to see,the portion under curves within the range of vibration frequency dominated by each element is fi lled with red color.Distances between pairs of atoms in(a)–(c)labeled in panel(d)are both 9.2 ?.

3.4 YbFe4Sb12晶格熱導(dǎo)率的計(jì)算與修正

如圖5所示,在YbFe4Sb12的晶格熱導(dǎo)率的計(jì)算上做了以下的幾點(diǎn)嘗試.首先用300 K時(shí)的二階和三階力常數(shù),計(jì)算得到的晶格熱導(dǎo)率隨溫度變化的關(guān)系,不出意料地滿足T?1的規(guī)律,這是其散射機(jī)理以三聲子過程為主所決定的.然后,考慮填充原子作用可能受溫度的影響,我們通過AIMD和TDEP方法得到了300—600 K幾個(gè)不同溫度下的力常數(shù),在三聲子框架下分別計(jì)算得到各自溫度下的晶格熱導(dǎo)率.可以看到其溫度依賴關(guān)系不再為T?1的關(guān)系,變得和實(shí)驗(yàn)數(shù)據(jù)的趨勢(shì)有些接近,這說明考慮溫度對(duì)于非諧或者填充原子作用的影響的思路是對(duì)的.但是在數(shù)值上,采用不同溫度力常數(shù)得到的結(jié)果與實(shí)驗(yàn)結(jié)果仍然存在著一定的差距.結(jié)合動(dòng)力學(xué)上原子關(guān)聯(lián)性討論的結(jié)果,我們認(rèn)為這種差距的原因在于三聲子散射過程的圖像不適用于填充原子與晶格聲子的相互作用.

圖5 用TDEP方法計(jì)算得到的YbFe4Sb12晶格熱導(dǎo)率κL結(jié)果.黑色虛線為300 K的單一力常數(shù)的計(jì)算結(jié)果,藍(lán)色方塊為不同溫度下力常數(shù)計(jì)算得到的κL,紅色圓點(diǎn)為Yb振動(dòng)模式壽命經(jīng)過共振散射修正后所得到的κL,紫色五角星為YbFe4Sb12的實(shí)驗(yàn)晶格熱導(dǎo)率結(jié)果,藍(lán)色、紅色線為各自數(shù)據(jù)擬合所得趨勢(shì)線Fig.5.Calculated lattice thermal conductivity κLof YbFe4Sb12vs.temperature through TDEP method.Black dash line is the calculated κLwith 2nd/3rdforce constants at 300 K and is proportional to 1/T.Blue squares are κLof 300–600 K calculated by 2nd/3rd force constants at each temperature.Red dots are κL got by modi fi ed lifetimes of Yb-dominant modes by resonant scattering form.Blue dash line and red solid line are fi tted by the data of blue squares and red dots respectively.Purple pentagrams are experimental κL of YbFe4Sb12from Ref.[38].

前面關(guān)于填充原子Yb的振動(dòng)與聲學(xué)支聲子耦合,并降低其壽命的表述,也適合于共振散射的解釋.用共振散射的壽命形式[37],即

作為由Yb所主導(dǎo)的局域振動(dòng)模式的壽命,其他模式的聲子壽命保持不變,代入到熱導(dǎo)率的計(jì)算中,即

則可以將晶格熱導(dǎo)率大致降低到了和實(shí)驗(yàn)結(jié)果[38]相近的水平.這說明了填充原子降低晶格熱導(dǎo)率的原因,需要用超出三聲子過程的非微擾的散射機(jī)理來解釋,而共振散射是其中的一種可能.

4 結(jié) 論

以全填充方鈷礦YbFe4Sb12為例,采用第一性原理分子動(dòng)力學(xué)模擬(AIMD),并結(jié)合溫度相關(guān)有效勢(shì)場(chǎng)(TDEP)的方法研究了填充原子在熱輸運(yùn)過程中的作用,發(fā)現(xiàn)填充原子Yb的振動(dòng)模式具有強(qiáng)烈的局域特征,其降低晶格熱導(dǎo)率的機(jī)理超出了三聲子倒逆過程的范疇.由于AIMD充分考慮了有限溫度下填充原子與框架的相互作用,而TDEP方法將之納入到了等效的非諧力常數(shù)內(nèi).在此基礎(chǔ)上計(jì)算得到的聲子輸運(yùn)性質(zhì)表明,Yb局域振動(dòng)模式通過與聲學(xué)支耦合,確實(shí)大幅降低了聲學(xué)支模式的壽命,但是將之作為三聲子倒逆過程考慮,理論計(jì)算的晶格熱導(dǎo)和實(shí)驗(yàn)結(jié)果存在著嚴(yán)重的偏離.而且動(dòng)力學(xué)過程中原子運(yùn)動(dòng)的關(guān)聯(lián)性質(zhì)分析表明,填充原子Yb的運(yùn)動(dòng)存在明確的局域和獨(dú)立特征,不符合聲子周期性的集體振動(dòng)圖像,所以填充原子散射晶格聲子、貢獻(xiàn)熱阻的機(jī)理具有超出傳統(tǒng)三聲子框架的特征.將填充原子Yb局域振動(dòng)模式的壽命做共振散射的修正,可以將晶格熱導(dǎo)率的計(jì)算結(jié)果降到和實(shí)驗(yàn)相符的水平.因此,在類似于Yb原子填充體系中,即在填充方鈷礦體系中存在強(qiáng)烈且孤立的局域振動(dòng)模時(shí),填充原子的局域振動(dòng)通過與聲學(xué)支聲子耦合,以共振散射為主、聲子間相互作用為輔的形式大幅降低傳熱聲子的壽命,進(jìn)而降低了材料的晶格熱導(dǎo)率.

[1]Shi X,Xi L L,Yang J,Zhang W Q,Chen L D 2011Physics40 710(in Chinese)[史迅,席麗麗,楊炯,張文清,陳立東2011物理40 710]

[2]Nolas G S,Morelli D T,Tritt T M 1999Annu.Rev.Mater.Sci.29 89

[3]Shi X,Bai S,Xi L,Yang J,Zhang W,Chen L,Yang J 2011J.Mater.Res.26 1745

[4]Rull-BravoM,MoureA,FernándezJF,Martín-González M 2015RSC Adv.5 41653

[5]Shi X,Yang J,Salvador J R,Chi M,Cho J Y,Wang H,Bai S,Yang J,Zhang W,Chen L 2011J.Am.Chem.Soc.133 7837

[6]Rogl G,Aabdin Z,Scha fl er E,Horky J,Setman D,Zehetbauer M,Kriegisch M,Eibl O,Grytsiv A,Bauer E 2012J.Alloys Compd.537 183

[7]Xi L L,Yang J,Shi X,Zhang W Q,Chen L D,Yang J H 2011Sci.China:Phys.Mech.Astron.41 706(in Chinese)[席麗麗,楊炯,史迅,張文清,陳立東,楊繼輝2011中國(guó)科學(xué):物理學(xué)力學(xué)天文學(xué)41 706]

[8]Slack G A,Tsoukala V G 1994J.Appl.Phys.76 1665

[9]Nolas G,Cohn J,Slack G 1998Phys.Rev.B58 164

[10]Huang L F,Li Y L,Ni M Y,Wang X L,Zhang G R,Zeng Z 2009Acta Phys.Sin.58 306(in Chinese)[黃良鋒,李延齡,倪美燕,王賢龍,張國(guó)仁,曾雉 2009物理學(xué)報(bào)58 306]

[11]Keppens V,Mandrus D,Sales B C,Chakoumakos B C,Dai P,Coldea R,Maple M B,Gajewski D A,Freeman E J,Bennington S 1998Nature395 876

[12]Hermann R P,Jin R,Schweika W,Grandjean F,Mandrus D,Sales B C,Long G J 2003Phys.Rev.Lett.90 135505

[13]Dimitrov I K,Manley M E,Shapiro S M,Yang J,Zhang W,Chen L D,Jie Q,Ehlers G,Podlesnyak A,Camacho J,Li Q 2010Phys.Rev.B82 174301

[14]Feldman J L,Singh D J,Mazin I I,Mandrus D,Sales B C 2000Phys.Rev.B61 R9209

[15]Koza M M,Johnson M R,Viennois R,Mutka H,Girard L,Ravot D 2008Nat.Mater.7 805

[16]Li W,Mingo N 2015Phys.Rev.B91 144304

[17]Qiu W,Xi L,Wei P,Ke X,Yang J,Zhang W 2014Proc.Natl.Acad.Sci.USA111 15031

[18]Qiu W,Ke X,Xi L,Wu L,Yang J,Zhang W 2016Sci.China:Phys.Mech.Astron.59 627001

[19]Broido D A,Malorny M,Birner G,Mingo N,Stewart D A 2007Appl.Phys.Lett.91 231922

[20]Togo A,Oba F,Tanaka I 2008Phys.Rev.B78 134106

[21]Li W,Carrete J,A.Katcho N,Mingo N 2014Comput.Phys.Commun.185 1747

[22]Hellman O,Steneteg P,Abrikosov I A,Simak S I 2013Phys.Rev.B87 104111

[23]Hellman O,Abrikosov I A 2013Phys.Rev.B88 144301

[24]Srivastava G P 1990The Physics of Phonons(Boca Raton:CRC press)p88

[25]Hellman O,Broido D A 2014Phys.Rev.B90 134309

[26]Li C W,Hellman O,Ma J,May A F,Cao H B,Chen X,Christianson A D,Ehlers G,Singh D J,Sales B C,Delaire O 2014Phys.Rev.Lett.112 175501

[27]Slack G A,Galginaitis S 1964Phys.Rev.133 A253

[28]Chen L D,Kawahara T,Tang X F,Goto T,Hirai T,Dyck J S,Chen W,Uher C 2001J.Appl.Phys.90 1864

[29]Nolas G S,Fowler G,Yang J 2006J.Appl.Phys.100 043705

[30]Guo R,Wang X,Huang B 2015Sci.Rep.5 7806

[31]Hafner J,Krajci M 1993J.Phys.:Condens.Matter5 2489

[32]Pailhes S,Euchner H,Giordano V M,Debord R,Assy A,Gomes S,Bosak A,Machon D,Paschen S,de Boissieu M 2014Phys.Rev.Lett.113 025506

[33]Euchner H,Pailhès S,Nguyen L T K,Assmus W,Ritter F,Haghighirad A,Grin Y,Paschen S,de Boissieu M 2012Phys.Rev.B86 224303

[34]Zhao X Y,Shi X,Chen L D,Zhang W Q,Bai S Q,Pei Y Z,Li X Y,Goto T 2006Appl.Phys.Lett.89 092121

[35]Cowley R A 1968Rep.Prog.Phys.31 123

[36]Christensen M,Abrahamsen A B,Christensen N B,Juranyi F,Andersen N H,Lefmann K,Andreasson J,Bahl C R,Iversen B B 2008Nat.Mater.7 811

[37]Pohl R 1962Phys.Rev.Lett.8 481

[38]Qiu P F,Yang J,Liu R H,Shi X,Huang X Y,Snyder G J,Zhang W,Chen L D 2011J.Appl.Phys.109 063713