王甫 周毅 高士鑫 段振剛 孫志鵬 汪俊 鄒宇 付寶勤?

1)(四川大學(xué)原子核科學(xué)技術(shù)研究所,輻射物理及技術(shù)教育部重點(diǎn)實(shí)驗(yàn)室,成都 610064)

2)(中國核動力研究設(shè)計(jì)院,核反應(yīng)堆系統(tǒng)設(shè)計(jì)技術(shù)重點(diǎn)實(shí)驗(yàn)室,成都 610200)

碳化硅(SiC)由于性能優(yōu)異,已廣泛應(yīng)用于核技術(shù)領(lǐng)域.在輻照環(huán)境下,載能入射粒子可使材料中的原子偏離晶體格點(diǎn)位置,進(jìn)而產(chǎn)生過飽和的空位、間隙原子、錯位原子等點(diǎn)缺陷,這些缺陷將改變材料的熱物性能,劣化材料的服役性能.因此,本文利用平衡分子動力學(xué)方法(Green-Kubo 方法)采用Tersoff型勢函數(shù)研究了點(diǎn)缺陷對立方碳化硅(β-SiC 或 3C-SiC)熱傳導(dǎo)性能的影響規(guī)律.研究過程中考慮的點(diǎn)缺陷包括:Si 間隙原子(SiI)、Si 空位(SiV)、Si 錯位原子(SiC)、C 間隙原子(CI)、C 空位(CV)和C 錯位原子(CSi).研究結(jié)果表明,熱導(dǎo)率(λ)隨點(diǎn)缺陷濃度(c)的增加而減小.在研究的點(diǎn)缺陷濃度范圍(點(diǎn)缺陷與格點(diǎn)的比例范圍為0.2%—1.6%),額外熱阻率(ΔR=Rdefect—Rperfect,R=1/λ,Rdefect 為含缺陷材料的熱阻率,Rperfect 為不含缺陷材料的熱阻率)與點(diǎn)缺陷的濃度呈線性關(guān)系,其斜率為熱阻率系數(shù).研究表明:空位和間隙原子的熱阻率系數(shù)高于錯位原子的熱阻率系數(shù);高溫下點(diǎn)缺陷的熱阻率系數(shù)高于低溫下點(diǎn)缺陷的熱阻率系數(shù);Si 空位和Si 間隙原子的熱阻率系數(shù)高于C 空位和C 間隙原子的熱阻率系數(shù).這些結(jié)果有助于預(yù)測及調(diào)控輻照條件下碳化硅的熱傳導(dǎo)性能.

1 引言

碳化硅(SiC)因具有優(yōu)異的性能,如抗輻照能力強(qiáng)、化學(xué)穩(wěn)定性好、熱導(dǎo)率高[1]等,在核技術(shù)領(lǐng)域備受關(guān)注[2-4].由于SiC 堆垛層錯的形成能低[5],所以具有多型體現(xiàn)象,即成分相同,但構(gòu)型和物理性質(zhì)有差異[6].目前,已經(jīng)發(fā)現(xiàn)的SiC 多型體有250 多種[7].其中,受關(guān)注比較多的是立方結(jié)構(gòu)(β-SiC 或3C-SiC),六方結(jié)構(gòu)(如4H-SiC 和6H-SiC).在高溫氣冷堆[8]中,3C-SiC 作為TRISO (tristructural-isotropic)型包覆燃料顆粒(簡稱TRISO 顆粒)中的包覆層[9],可阻擋裂變產(chǎn)物的釋放,且可高效導(dǎo)出核芯的熱能,避免堆芯溫度過高.然而,在輻照環(huán)境下,載能入射粒子(如中子和裂變碎片等)可使材料中的原子偏離晶體格點(diǎn)位置,引起離位損傷,進(jìn)而產(chǎn)生過飽和的空位、間隙原子、錯位原子等點(diǎn)缺陷.3C-SiC 中的基本點(diǎn)缺陷類型包括:Si 間隙原子(SiI)、Si 空位(SiV)、Si 錯位原子(SiC)、C 間隙原子(CI)、C 空位(CV)和C 錯位原子(CSi),這些缺陷是構(gòu)成復(fù)合缺陷的基本單元.Katoh 等[10]研究了中子輻照下3C-SiC 的微觀結(jié)構(gòu)變化,結(jié)果表明,在相對較低的溫度(≤800 ℃)下,輻照后材料里的缺陷主要為點(diǎn)缺陷及小間隙團(tuán)簇.Devanathan 等[11]和Ran 等[12]基于分子動力學(xué)研究表明,基本撞擊原子輻照產(chǎn)生的剩余點(diǎn)缺陷主要為Frenkel 對和錯位原子.這些缺陷會影響材料的熱物性能,在輻照條件下,SiC 材料的熱導(dǎo)率會降低[13].這可能會導(dǎo)致堆芯過熱進(jìn)而造成TRISO 顆粒失效,影響反應(yīng)堆的安全運(yùn)行.因而研究點(diǎn)缺陷對3C-SiC 熱傳導(dǎo)性能的影響規(guī)律有重要意義.

在實(shí)驗(yàn)中獲得單種點(diǎn)缺陷類型的試樣和準(zhǔn)確調(diào)控缺陷濃度是極其困難的[14],因而有必要借助計(jì)算機(jī)模擬技術(shù)來研究缺陷對材料熱導(dǎo)率的影響.分子動力學(xué)(molecular dynamics,MD)是一種合適研究晶體缺陷對熱導(dǎo)率影響的方法.只需針對待研究體系建立適當(dāng)?shù)膭莺瘮?shù),就可以模擬其主要的熱力學(xué)性質(zhì)和動力學(xué)性質(zhì).分子動力學(xué)計(jì)算材料聲子熱導(dǎo)率的方法主要包括平衡分子動力學(xué)(equilibrium molecular dynamics,EMD)和非平衡分子動力學(xué)(non-equilibrium molecular dynamics,NEMD)兩種方法[15,16].Crocombette 等[17]利用NEMD 方法研究了輻照對β-SiC 熱導(dǎo)率的影響.相似的方法被用來計(jì)算SiC 納米管、納米線的熱導(dǎo)率[18].這些研究表明了分子動力學(xué)研究SiC 中點(diǎn)缺陷對熱傳導(dǎo)性能影響的可行性.然而NEMD 方法的不足之處是存在明顯的尺寸效應(yīng)[19],計(jì)算過程中存在較大的溫度梯度引起非線性的人為現(xiàn)象,且僅能獲得給定熱流方向的熱導(dǎo)率.EMD 方法避免了NEMD的主要不足,在計(jì)算過程中系統(tǒng)長時間保持恒定溫度,始終處于線性響應(yīng)狀態(tài).因而本文采用EMD方法研究點(diǎn)缺陷對3C-SiC 熱導(dǎo)率的影響規(guī)律.研究表明,含有點(diǎn)缺陷的SiC 材料的熱阻率與點(diǎn)缺陷的濃度(考慮的點(diǎn)缺陷與格點(diǎn)比范圍為0.2%—1.6%)呈線性關(guān)系,對應(yīng)的斜率為熱阻率系數(shù);空位和間隙原子的熱阻率系數(shù)高于錯位原子的熱阻率系數(shù);高溫下點(diǎn)缺陷的熱阻率系數(shù)高于低溫下點(diǎn)缺陷的熱阻率系數(shù);Si 空位和Si 間隙原子的熱阻率系數(shù)高于C 空位和C 間隙原子的熱阻率系數(shù).

本文余下部分組織如下:第2 節(jié)給出熱導(dǎo)率計(jì)算所用的方法和具體的模擬細(xì)節(jié);第3 節(jié)結(jié)果與討論,包括熱流密度的時間自相關(guān)函數(shù)收斂性分析和不同點(diǎn)缺陷及其濃度對熱傳導(dǎo)性能的影響;最后是全文的結(jié)論.

2 方法

2.1 平衡分子動力學(xué)

在分子模擬中,計(jì)算熱導(dǎo)率的EMD 方法一般是指利用Green-Kubo 公式[20]計(jì)算熱流密度的時間自相關(guān)函數(shù)(HFACF)來求解體系聲子熱導(dǎo)率的方法.計(jì)算SiC 熱導(dǎo)率(λ)的Green-Kubo 公式為

式中kB為玻爾茲曼常數(shù);T為系統(tǒng)溫度;V表示β-SiC 超胞體積;J為系統(tǒng)瞬時微觀熱流密度;〈·〉表示HFACF 的系綜統(tǒng)計(jì)平均,在MD 計(jì)算過程中通過時間平均來計(jì)算系綜平均.微觀熱流密度的計(jì)算公式為

式中t為時間,εi表示原子i能量,vi表示原子i的速度,Si表示原子i的virial 應(yīng)力.

2.2 模擬細(xì)節(jié)

采用LAMMPS 軟件包進(jìn)行MD 模擬[21].采用Tersoff型[22]勢函數(shù)[23]描述SiC 中原子間的相互作用,該勢函數(shù)可較好地描述3C-SiC 的點(diǎn)缺陷相對穩(wěn)定性和彈性性能及熱學(xué)性能.在計(jì)算過程中考慮了不同尺寸的超胞,包括4a×4a×4a(512),5a×5a×5a(1000),6a×6a×6a(1728)和7a×7a×7a(2744),其中a表示立方SiC 結(jié)構(gòu)胞的晶胞參數(shù),括號中的數(shù)字表示相應(yīng)尺寸完美超胞中的原子數(shù),超胞的3 個軸取向分別為[100],[010]和[001].為了在超胞中引入點(diǎn)缺陷,通過隨機(jī)刪除一定數(shù)量的C(Si)格點(diǎn)原子構(gòu)建一定濃度的C(Si)空位,通過隨機(jī)選擇一定數(shù)量的Si(C)格點(diǎn)改變?yōu)镃(Si)原子構(gòu)成C(Si)錯位原子,通過在超胞的隨機(jī)位置引入一定數(shù)量的C(Si)原子構(gòu)成C(Si)間隙原子.所有超胞x,y和z方向均采用周期性邊界條件.在MD 運(yùn)行前,先對含有點(diǎn)缺陷的超胞進(jìn)行原子位置最優(yōu)化處理,采用共軛梯度算法,使各原子優(yōu)化至局部最穩(wěn)定位置.MD 運(yùn)行的時間步長是0.3 fs,經(jīng)過測試,發(fā)現(xiàn)這個時間步長足夠小,具有足夠精度獲得可靠的計(jì)算結(jié)果.針對優(yōu)化后的超胞,對原子進(jìn)行熱化以接近目標(biāo)溫度,本文考慮的溫度包括600 和1500 K.然后通過Nose-Hoover調(diào)溫技術(shù)[24]和Parrinello-Rahman 調(diào)壓技術(shù)[25]產(chǎn)生等溫等壓系綜(NPT 系綜),使用時間可逆速度Verlet 算法[26]運(yùn)行約160000 步(48 ps).正式采樣時改為微正則系綜(NVE 系綜)運(yùn)行MD,每隔20 步(即6 fs)取當(dāng)前的瞬態(tài)構(gòu)型(包括當(dāng)前的超胞軸取向、各個原子的坐標(biāo)、速度和受力)以計(jì)算當(dāng)前的瞬間微觀熱流密度(如(2)式所述),然后通過時間平均來獲得超胞的HFACF,最后根據(jù)(1)式對HFACF 進(jìn)行積分計(jì)算熱導(dǎo)率.MD 計(jì)算的總步數(shù)至少32000000 步(9600 ps),關(guān)聯(lián)時間至少大于15 ps,以保證計(jì)算的收斂和可靠性.另外為了防止計(jì)算的偶然性,又分別進(jìn)行了至少20 次獨(dú)立的計(jì)算,以獲得熱導(dǎo)率的統(tǒng)計(jì)平均值.

3 結(jié)果與討論

3.1 熱流密度的時間自相關(guān)函數(shù)(HFACF)

基于瞬間構(gòu)型提煉超胞的體積和各原子的動能、勢能、Virial 應(yīng)力和速度等物理量,體系的微觀熱流密度可根據(jù)(2)式獲得,熱流密度隨時間的變化如圖1(a)所示,這里以完美超胞(6a×6a×6a)在1500 K 時的模擬結(jié)果為例進(jìn)行分析,可以發(fā)現(xiàn)超胞不同方向的微觀熱流密度(Jx,Jy和Jz)并無明顯差別,這是由于3C-SiC 的3 個軸向是完全等效的.微觀熱流密度數(shù)值均在0 eV/(ps·?2)上下變動,變動范圍約為—0.2—0.2 eV/(ps·?2).該超胞在1500 K 時HFACF 與時間的關(guān)系如圖1(b)所示,可以看到HFACF 在初始時變化較大,在時間達(dá)到10 ps 后減小至接近為0,因此可以認(rèn)為這時熱流密度不再時間相關(guān)或相關(guān)較小,因而如前文所述,關(guān)聯(lián)時間設(shè)置為15 ps.因此,可用有限時間的MD 模擬,通過HFACF 與時間的關(guān)系,運(yùn)用(1)式對HFACF 進(jìn)行一定時間的積分可獲得熱導(dǎo)率.從圖1(b)還可以發(fā)現(xiàn),超胞的取向?qū)FACF 無影響,同樣是由于3C-SiC 的高對稱性.除了關(guān)聯(lián)時間,使用Green-Kubo 公式計(jì)算體系的熱導(dǎo)率時,還需足夠長的時間進(jìn)行累積平均,相當(dāng)于系綜平均的時間,本文要求正式采樣的MD 運(yùn)行總時間不少于9600 ps,可保證計(jì)算結(jié)果的收斂.本文測試了不同尺寸完美超胞的熱導(dǎo)率,如圖2 所示.可以看出,在尺寸為6a×6a×6a時熱導(dǎo)率收斂至67.8 W/mK(1500 K),是與實(shí)驗(yàn)值(62.6—68.6 W/mK[27])接近.而600 K 時完美超胞的熱導(dǎo)率為325 W/mK,若進(jìn)行量子修正[28]也與實(shí)驗(yàn)值接近.

圖1 (a)微觀熱流密度(J)隨時間的變化關(guān)系;(b)微觀熱流密度的時間自相關(guān)函數(shù)(HFACF)隨時間的變化關(guān)系.溫度為1500 K,完美超胞的尺寸為6a×6a×6aFig.1.(a)Relationship of microscopic heat flux (J)with time;(b)relationship of heat flux autocorrelation function(HFACF)with time.The temperature is 1500 K,the size of the perfect supercell is 6a×6a×6a.

圖2 尺寸對熱導(dǎo)率的影響Fig.2.Effect of size on the thermal conductivity.

3.2 點(diǎn)缺陷對熱導(dǎo)率的影響

如前文所述,SiC 材料在輻照環(huán)境下服役,將受到中子和裂變碎片等載能粒子轟擊,引起材料中格點(diǎn)原子離位,進(jìn)而產(chǎn)生過飽和的空位、間隙原子、錯位原子等點(diǎn)缺陷.點(diǎn)缺陷的種類和濃度受到載能粒子(能量、通量和劑量等)、SiC 材料狀態(tài)(結(jié)構(gòu))和其他服役條件(溫度和壓力等)等因素的影響[10-12].為了研究點(diǎn)缺陷類型及其濃度對3CSiC 熱傳導(dǎo)性能的影響,本文針對不同的點(diǎn)缺陷類型,構(gòu)建了不同點(diǎn)缺陷濃度的超胞.根據(jù)前文的參數(shù)設(shè)置及計(jì)算流程,計(jì)算了不同類型點(diǎn)缺陷在不同濃度下的熱導(dǎo)率,結(jié)果如圖3 所示.可以發(fā)現(xiàn)隨著點(diǎn)缺陷濃度的增加,SiC 材料的熱導(dǎo)率(λ)下降.當(dāng)然不同類型的點(diǎn)缺陷,對應(yīng)超胞的熱導(dǎo)率數(shù)值不同,這是由不同類型點(diǎn)缺陷對聲子的散射行為的差異造成的.另外從圖3 還可以發(fā)現(xiàn),錯位原子(尤其是CSi)對應(yīng)超胞的熱導(dǎo)率較大,即對聲子的散射較小;而Si型點(diǎn)缺陷(尤其是SiI)對應(yīng)超胞的熱導(dǎo)率較小,即對聲子的散射較大.對比圖3(a)和圖3(b)可以發(fā)現(xiàn),溫度不影響上述點(diǎn)缺陷類型及濃度與熱導(dǎo)率的主要規(guī)律,僅影響熱導(dǎo)率數(shù)值的大小.另外研究過程中還發(fā)現(xiàn),間隙原子構(gòu)型(對于本文使用的Tersoff 勢函數(shù),C 間隙的穩(wěn)定構(gòu)型為C+-C〈100〉,Si 間隙的最穩(wěn)定構(gòu)型為Si+-C〈100〉)將影響超胞的熱傳導(dǎo)行為,dumbbell 間隙軸向(即其中一個〈100〉取向)的熱導(dǎo)率要小于橫向(即另外兩個〈100〉取向)的熱導(dǎo)率.

圖3 含有點(diǎn)缺陷的SiC 超胞熱導(dǎo)率(λ)與點(diǎn)缺陷濃度(c,相當(dāng)于點(diǎn)缺陷與格點(diǎn)的比例)的關(guān)系 (a)600 K;(b)1500 KFig.3.Relationship between the thermal conductivity (λ)of SiC supercells containing point defects and the concentration of point defects (c,equivalent to the ratio of point defects to lattice points):(a)600 K;(b)1500 K.

SiC 熱傳導(dǎo)的主要載體是聲子,在溫度高的區(qū)域晶格振動具有更大的振幅和更多的模式,即聲子數(shù)更多,這些聲子將傳遞至低溫區(qū)域,然而聲子間存在相互作用,傳遞過程將發(fā)生碰撞,在材料中引入點(diǎn)缺陷,也將影響聲子的傳遞,聲子會與缺陷發(fā)生碰撞,降低聲子壽命.根據(jù)聲子的玻爾茲曼輸運(yùn)理論[29],(其中Cμ為對應(yīng)聲子模μ比熱容,vαμ是α方向聲子模μ的群速度,τμ為聲子模μ的馳豫時間),晶格熱導(dǎo)率主要取決于聲子的比熱容、群速度和壽命(馳豫時間).點(diǎn)缺陷降低了聲子壽命,從而使得熱導(dǎo)率減小.也就是說引入的點(diǎn)缺陷將增加材料的熱阻率(R=1/λ),我們將點(diǎn)缺陷引起的熱阻率增加部分稱為額外熱阻率(ΔR=Rdefect—Rperfect,其中Rdefect為含有點(diǎn)缺陷材料的熱阻率,Rperfect為不含點(diǎn)缺陷材料的熱阻率),顯然額外熱阻率受到點(diǎn)缺陷濃度的影響,點(diǎn)缺陷的濃度越高,額外熱阻率越大,在較小的濃度范圍內(nèi),可以認(rèn)為額外熱阻率隨點(diǎn)缺陷濃度的升高而線性升高.如圖4 所示,額外熱阻率與點(diǎn)缺陷濃度確實(shí)呈現(xiàn)較好的線性關(guān)系,該結(jié)果與Crocombette和Proville[30]的結(jié)論一致.那么圖4 中ΔR-c關(guān)系可用以下公式來描述:

式中l(wèi)定義為熱阻率系數(shù);ΔR0為擬合系數(shù),相當(dāng)于點(diǎn)缺陷濃度為0 時的額外熱阻率,該值應(yīng)該為0,事實(shí)上,擬合結(jié)果發(fā)現(xiàn)ΔR0基本上接近于0.另外應(yīng)該說明的是,這種ΔR-c的線性關(guān)系僅存在于較低的濃度范圍內(nèi),對于高的點(diǎn)缺陷濃度,點(diǎn)缺陷之間將存在相互作用,甚至相互反應(yīng)生成復(fù)雜的點(diǎn)缺陷團(tuán)簇,進(jìn)而影響點(diǎn)缺陷對聲子的散射過程.

顯然根據(jù)(3)式,在確定了某溫度下各類型點(diǎn)缺陷的熱阻率系數(shù)后,可估算低點(diǎn)缺陷濃度下SiC 材料的熱導(dǎo)率,預(yù)測計(jì)算的公式為

(4)式假定了各種類型點(diǎn)缺陷的熱阻率系數(shù)是相互獨(dú)立的,顯然僅在低點(diǎn)缺陷濃度下成立.另外該預(yù)測公式未考慮實(shí)際材料中位錯、晶界、層錯和相界等缺陷結(jié)構(gòu)的影響.

根據(jù)圖4 中的擬合曲線可得到在溫度為600和1500 K 條件下各種點(diǎn)缺陷的熱阻率系數(shù),結(jié)果如圖5 所示.可以看出,空位和間隙型點(diǎn)缺陷引起的熱阻率系數(shù)一般比錯位型點(diǎn)缺陷的熱阻率系數(shù)高;高溫下點(diǎn)缺陷的熱阻率系數(shù)一般比低溫下點(diǎn)缺陷的熱阻率系數(shù)高;Si 類型點(diǎn)缺陷(Si 空位和Si間隙)引起的熱阻率系數(shù)一般比C 類型(C 空位和C 間隙)的熱阻率系數(shù)高.這是因?yàn)樵跓o機(jī)非金屬中,熱傳導(dǎo)的主要載體是聲子,即晶格振動的能量量子[31].因此,熱阻(L/λ,L為熱流方向材料長度)主要是由聲子之間的相互作用和晶體缺陷散射聲子引起的[32].空位或間隙引起了最大的質(zhì)量差和半徑差,這將導(dǎo)致更強(qiáng)的點(diǎn)缺陷-聲子散射[33,34],從而導(dǎo)致更高的熱阻率系數(shù).

圖4 含有點(diǎn)缺陷的SiC 超胞的額外熱阻率(ΔR)與點(diǎn)缺陷濃度(c,相當(dāng)于點(diǎn)缺陷與格點(diǎn)的比例)的關(guān)系 (a)600 K;(b)1500 KFig.4.Relationship between the excess thermal resistance(ΔR)of SiC supercells containing point defects and the concentration of point defects (c,equivalent to the ratio of point defects to lattice points):(a)600 K;(b)1500 K.

圖5 各種點(diǎn)缺陷不同溫度下的熱阻率系數(shù)l (單位為mK/W)Fig.5.Thermal resistivity coefficient l (in mK/W)of various point defects at different temperatures.

4 結(jié)論

本文利用EMD 方法(Green-Kubo 公式),研究了6 種點(diǎn)缺陷(Si 間隙原子(SiI)、Si 空位(SiV)、Si 錯位原子(SiC)、C 間隙原子(CI)、C 空位(CV)和C 錯位原子(CSi))在不同濃度(點(diǎn)缺陷與格點(diǎn)比為0.2%—1.6%)下對立方碳化硅(β-SiC 或3CSiC)聲子熱導(dǎo)率的影響規(guī)律.研究發(fā)現(xiàn)熱導(dǎo)率(λ)隨點(diǎn)缺陷濃度(c)的增加而減小.在研究的點(diǎn)缺陷濃度范圍內(nèi),額外熱阻率(即點(diǎn)缺陷引起的熱阻率的增加值)與點(diǎn)缺陷的濃度呈線性關(guān)系,間隙原子構(gòu)型取向使熱傳導(dǎo)行為呈現(xiàn)各向異性.以熱阻率系數(shù)作為衡量點(diǎn)缺陷對材料熱阻率(熱導(dǎo)率)影響程度大小的影響因子,空位和間隙引起的熱阻率系數(shù)高于錯位原子的熱阻率系數(shù);高溫下點(diǎn)缺陷的熱阻率系數(shù)高于低溫下點(diǎn)缺陷的熱阻率系數(shù);Si 空位和Si 間隙比C 空位和C 間隙引起更高熱阻率系數(shù).這些結(jié)果可為預(yù)測輻照條件下SiC 熱導(dǎo)率以及調(diào)控SiC 熱傳導(dǎo)性能提供參考.