空位缺陷和原子摻雜對氮化鎵熱導率影響的分子動力學模擬

高志樂 ,王繼芬

(1.上海海事大學商船學院,上海 201306;2.上海第二工業大學 資源與環境工程學院,上海 201209)

近年來,由于氮化鎵(GaN)具有電子飽和率高、擊穿電壓大、抗腐蝕性強、穩定性好等優點[1-2],基于GaN 材料的光電子器件[3]、高溫大功率器件和高頻微波器件[4]得到了廣泛應用。目前設備中的絕大多數GaN 晶體是通過在基片上異質外延生長獲得的[5],這些襯底與GaN 之間的晶格常數和熱膨脹系數存在著不同程度的不匹配,導致了實際制得的GaN 存在很多缺陷[6],如空位缺陷和原子摻雜。缺陷的存在對GaN 的光學、電學和熱學等相關性能有非常重要的影響[7],嚴重影響了GaN 基半導體設備的使用壽命、發光效率、電子遷移率以及熱傳導等各方面性能[8]。

實際應用中,GaN 往往存在空位缺陷及原子摻雜的情況,這是研究人員在研究GaN 時比較關注的。Cai 等[9]采用第一性原理計算研究了缺鎵空位(VGa)和缺氮空位(VN)缺陷對GaN 電子結構和光學性質的影響,研究發現VGa 中價帶向更高能量區域移動,VN 中價帶和導帶明顯向低能區移動。空位缺陷的存在導致GaN 的帶隙增大,使得GaN 的光學性質發生了改變。Abdalla 等[10]利用第一性原理系統研究了Fe 摻雜對纖鋅礦型GaN 電子和磁性的影響。研究發現,當Fe 原子在Ga 位置摻雜時,GaN 由半導體狀態轉變為磁性半導體狀態,從而導致Ga 位Fe 摻雜的GaN 比N位Fe 摻雜的GaN 更穩定。Beechem 等[11]通過時域熱反射率研究了n 型和p 型摻雜的GaN 的熱導率,結果顯示,當GaN 厚度增大時,n 型GaN 的熱導率幾乎恒定,而p 型GaN 的熱導率會出現下降的趨勢。Jezowski 等[12]通過實驗研究了摻雜的自由載流子對GaN 晶體中聲子散射的貢獻,實驗結果表明,通過電場影響了聲子的散射率、鍵硬度和離子性以及晶體對稱性,從而明顯改變了GaN 的熱導率。為了解決空位缺陷和原子摻雜給GaN 材料性能帶來的影響,需要進一步了解GaN 空位缺陷和原子摻雜的基本特性及其變化規律。

本文利用非平衡分子動力學(NEMD)的方法研究了空位缺陷和原子摻雜等缺陷對GaN 熱導率的影響,介紹了GaN 勢函數的選取以及模擬計算過程主要的參數設置;對GaN 的完整結構、空位缺陷結構以及原子摻雜結構的計算結果進行了分析與討論。該研究可為GaN 材料器件的熱管理應用提供一定的理論參考。

1 模擬模型和計算方法

1.1 GaN 及其空位缺陷和原子摻雜結構

GaN 晶體通常有纖鋅礦(WZ)、閃鋅礦(ZB)和巖鹽礦(RS)三種結構,如圖1 所示。因為WZ 結構的GaN 是穩定的,ZB 結構的GaN 是易變的,而RS 結構GaN 只在高壓情況下才會出現,所以GaN 晶體通常為WZ 結構,其單元格對應的空間群結構為P63mc(C46V)。本文中,以WZ 結構的GaN 作為研究對象,討論了空位缺陷和原子摻雜對GaN 的C 面(0001)方向熱導率影響。

圖1 (a) 纖鋅礦結構;(b) 閃鋅礦結構;(c) 巖鹽礦結構Fig.1 (a) WZ structure;(b) ZB structure;(c) RS structure

研究表明,空位缺陷是GaN 的本征缺陷。VGa 是在禁帶附近引入自旋極化缺陷態的受主缺陷;VN 是一種施主缺陷,它導致費米能級的升高。GaN 結構中的本征空位點缺陷將引起低能區附近光學性質的顯著變化[9]。有研究顯示GaN 中摻雜Fe,將在GaN 結構中形成一個深受體中心,和一些電子生成的內在缺陷補償通過機動的孔,從而賦予GaN 材料高電阻特性[13]。除了Fe 的摻雜,Mg 和Si 也是GaN 晶體中的重要摻雜劑。Mg 是實現p 型GaN 導電性的唯一摻雜劑[14],Si 是實現n 型GaN 摻雜的重要摻雜劑,能有效提高發光效率[15]。GaN 的p 型和n 型摻雜是制備光電子器件的關鍵。

GaN 含有空位缺陷的結構如圖2 所示,包含了16 個Ga 原子和16 個N 原子的2×2×2 超晶胞模型。圖2 (a)、(b)和(c)分別是完整GaN、VGa-GaN 和VN-GaN,相應的系統分別是Ga16N16、Ga15N16 和Ga16N15。

圖2 (a) 無空位;(b) 缺鎵空位;(c) 缺氮空位Fig.2 (a) No vacancies;(b) VGa;(c) VN

GaN 含有原子摻雜的結構如圖3 所示,包含了16 個Ga 原子和16 個N 原子的2×2×2 超晶胞模型。Fe、Mg、Si 原子以3.13%的摻雜率直接取代超胞中的一個Ga 原子。相應的系統分別是Ga15Fe1N16、Ga15Mg1N16和Ga15Si1N16。

圖3 (a)無摻雜;(b) Fe 摻雜;(c) Mg 摻雜;(d) Si 摻雜Fig.3 (a) No doping;(b) Fe-doping;(c) Mg-doping;(d) Si-doping

1.2 非平衡分子動力學模擬

因NEMD 的原理與實驗測量相似,本文利用NEMD 的方法模擬計算GaN 的傳熱性能。計算時,通過熱浴的方法建立了熱源與熱匯兩個區域,在兩個區域之間有熱流通過并且產生穩定的溫度梯度。根據傅立葉的熱傳導定律計算GaN 的熱導率。

式中:Q,A,和k分別是熱流速率、橫截面積、溫度梯度和熱導率。能量從熱源流向熱匯的橫截面積表示如下:

式中:W是模型在X方向的長度;δ是模型在Y方向的長度。本文在計算過程中利用定熱流的方法來產生穩定的溫度梯度。

1.3 計算模型

圖4 為采用NEMD 方法計算GaN 熱導率的模型示意圖。Z軸方向設定為計算熱導率的方向與GaN 初始晶胞的(0001)方向相對應。X和Y方向設定為周期性邊界條件。傳熱的Z方向設定為可收縮性邊界條件。在兩端各設置一個1.4 nm 的絕熱區域,以防模擬過程中的原子溢出。熱浴的寬度為4 nm,熱傳導區域的長度為52 nm,在GaN 的聲子自由程范圍內[16]。

圖4 計算模型區域劃分Fig.4 Area divisions of Calculation model

模擬計算中,時間步長設置為1 fs,熱源與熱匯之間設定一個恒定的熱流(q=0.5 eV/ps)。其中,系統中的原子在某一時刻的溫度可以表示如下:

式中:T是原子的溫度;m是原子的質量;Vi是i原子的速度;Nl是原子的總數;KB是玻爾茲曼常數。模擬計算時沿著傳熱方向將傳導區域平均分成40 塊,統計每塊內原子的平均溫度,從而得到GaN 晶體的熱源和熱匯之間沿傳熱方向的溫度分布。

聲子態密度(PDOS)和聲子參與率(PPR) 都是表征材料中聲子活動的有效方法。通過對PDOS 和PPR 的研究,分析空位缺陷和原子摻雜對GaN 的聲子輸運特性的影響。PDOS 是通過速度自相關函數(VACF)的傅立葉變換得到[17]。

式中:PDOS(ω) 表示為振動頻率為ω時的總PDOS,其VACF 由公式(5)表示[18]。

式中:Vi(t) 是粒子i在t時刻的速度矢量;N是原子總數;〈 〉是尖括號,表示集合平均;VACF(t) 是原子速度自相關函數。

計算GaN 的PPR 時,定義為[19]:

式中:N是原子總數;PDOSi(ω) 表示通過計算單粒子的VACF,得到的單粒子的PDOS。

1.4 勢函數

本文分別利用Stillinger-Weber (SW)勢函數[20]和Tersoff 勢函數[21-22]計算了完整GaN 的PDOS,計算結果如圖5 所示。結果表明利用SW 勢和Tersoff 勢得到的PDOS 在低頻段0～10 THz 范圍內吻合得比較好,在高頻段20～30 THz 范圍內,SW 勢態密度頻率比Tersoff 勢態密度頻率稍微偏大。通過圖5 可以發現,完整GaN 在SW 和Tersoff 兩種勢函數下得到的PDOS 相差不大。考慮到本文涉及了GaN 的空位缺陷和原子摻雜等晶體缺陷,SW 勢函數更適合GaN 晶體缺陷的情況,因此本文采用SW 勢函數作為計算過程中原子間的作用力函數,其參數是采用Béré 和Serra[23]校正后的值。

圖5 SW 和Tersoff 勢函數下完整GaN 的PDOS 對比Fig.5 PDOS of intact GaN obtained by SW potential and Tersoff potential

2 結果與討論

2.1 模型收斂性

對空位缺陷和原子摻雜GaN 結構進行模擬計算后,通過VACF[18]得到了各個GaN 模型的收斂性,如圖6 (a)～(f)所示分別為完整GaN、VGa-GaN、VNGaN 和GaN 摻雜Fe、Mg、Si 原子模型的收斂性圖。

圖6 (a) 完整GaN 收斂性;(b) VGa-GaN 收斂性;(c) VN-GaN 收斂性;(d) Fe 摻雜GaN 收斂性;(e) Mg 摻雜GaN 收斂性;(f) Si 摻雜GaN 收斂性Fig.6 (a) Convergence of intact GaN;(b) Convergence of VGa;(c) Convergence of VN;(d) Convergence of Fe-doped GaN;(e) Convergence of Mg-doped GaN;(f) Convergence of Si-doped GaN

通過圖6(a)～(f)可以發現,VACF 在初始階段都會有劇烈的波動,當后期時間趨近于3 ps 時,可以發現VACF 的振蕩減緩,趨近于0。根據圖6(a)～(f)可知,不同結構GaN 模型的收斂性在3 ps 附近的振蕩幅度不同。這是因為它們的結構存在著差異,從而導致他們收斂性的不同展現,這并不影響熱導率的計算[24-25]。通過圖6 可以證明在本文的模擬計算條件下,這六種結構的GaN 都能夠得到相當精確的熱導率值。

2.2 GaN 晶體的熱導率

利用NEMD 方法計算GaN 熱導率,得到如圖7 所示的熱源和熱匯之間良好的溫度梯度分布。通過圖7可以發現,忽略熱源和熱匯附近的非線性區域,GaN體系的溫度趨于線性分布。通過擬合熱源和熱匯之間的線性區域確定了GaN 體系的溫度梯度。根據公式(1)和(2)可得到GaN 傳熱區長度為52 nm 時的熱導率。

圖7 GaN 體系溫度分布圖Fig.7 GaN system temperature distribution chart

考慮垂直傳熱方向的尺寸效應,當GaN 沿傳熱方向的體系長度都為52 nm 時,計算得到不同截面積的GaN 熱導率如表1 所示。

表1 GaN 不同截面尺寸時的熱導率Tab.1 Thermal conductivity of GaN at different section sizes

通過表1 發現,四個GaN 算例的熱導率屬于同一數量級。當X方向尺寸不變,Y方向尺寸分別為1.2865,1.9298,2.5730 和3.8596 nm 時,可以發現GaN 熱導率的變化很小,最大值45.1 W·m-1·K-1與最小值42.3 W·m-1·K-1之間相差約2.8 W·m-1·K-1。按樣品重復誤差分析,則平均值約為43.7 W·m-1·K-1,體系誤差為±1.5 W·m-1·K-1,即3.4%。因此當使用NEMD 方法計算GaN (0001)方向的熱導率時,X和Y方向尺寸對結果有較小的影響,在精度要求不高時,可以不予考慮。

2.3 空位缺陷對GaN 熱導率的影響

圖8 所示是完整GaN、VGa-GaN 和VN-GaN 的熱導率在300～600 K 溫度范圍內的變化趨勢。通過圖8 可以發現,當溫度從300 K 變化到600 K 時,三種結構GaN 的熱導率都隨著溫度的升高呈下降趨勢。完整GaN 熱導率的下降速率快于VGa-GaN 和VN-GaN 熱導率的下降速率。這表明完整GaN 的熱導率和含有空位缺陷結構的GaN 熱導率對于溫度依賴性的一般趨勢都是隨著溫度的升高而減小,且完整GaN 的熱導率對溫度的依賴性更大。根據圖8 可知,當溫度為300 K時,完整GaN 的熱導率最大為43.06 W·m-1·K-1。VGa-GaN 和VN-GaN 的熱導率與完整GaN 的熱導率相比會出現較大的下降,但VGa 缺陷與VN 缺陷的GaN 熱導率相差不大,分別為10.83 W·m-1·K-1和8.66 W·m-1·K-1,與完整GaN 的熱導率相比分別下降74.85%和79.89%。這表明GaN 空位缺陷的存在對GaN 熱導率的影響較大。

圖8 完整GaN 與含空位缺陷GaN 的熱導率對比Fig.8 Thermal conductivity of intact GaN compared to GaN with vacancy defects

為了解空位缺陷對GaN 熱傳輸機理的影響,計算了300 K 時完整GaN、VGa-GaN 和VN-GaN 相關的PDOS,如圖9 (a)～(d) 所示。

圖9 (a) 300 K 時完整GaN 各原子的PDOS;(b) 300 K 時完整GaN 與兩種空位缺陷GaN 的Ga 原子PDOS 對比;(c) 300 K 時完整GaN 與兩種空位缺陷GaN 的N 原子PDOS 對比;(d) 300 K時完整GaN 與兩種空位缺陷GaN 的PDOS 對比Fig.9 (a)PDOS for each atom of intact GaN at 300 K;(b) Comparison of Ga atomic PDOS between intact GaN and two vacancy-deficient GaN at 300 K;(c) Comparison of N atomic PDOS between intact GaN and two vacancy-deficient GaN at 300 K;(d) Comparison of PDOS between intact GaN and two vacancy-deficient GaN at 300 K

通過圖9(a)可以發現,Ga 原子的振動峰在低頻區域(5～10 THz)具有較高的峰值,而N 原子的振動峰在高頻區域(20～30 THz)具有較高的峰值。圖9(a)的現象說明GaN 結構中的低頻聲子主要是Ga 原子起主要作用,高頻聲子主要是N 原子起主要作用,這和原子質量的相互作用是一致的。通過圖9(b)和(c)可以發現,在相同的條件下,完整GaN、VGa-GaN 和VN-GaN 中的Ga 原子與N 原子在低頻區域和高頻區域都存在著振動峰,且Ga 原子的低頻振動峰值高于高頻振動峰值,而N 原子的高頻振動峰值高于低頻振動峰值。根據圖9(d)可知,300 K 時,三種結構GaN的PDOS 曲線都在低頻區域和高頻區域存在著振動峰,但VGa-GaN 的PDOS 有一點特殊,它在15.2 THz 有一個振動峰。由于GaN 的熱導率主要受低頻區域聲子的影響[26],所以重點分析了低頻區域的聲子態密度。通過分析圖9(d)發現,在5.3 THz 處存在著一個振動峰,在這個振動峰處完整GaN 的PDOS 值＞VGa-GaN的PDOS 值＞VN-GaN 的PDOS 值。出現這種變化趨勢是因為當GaN 結構中出現空位缺陷時,GaN 結構中參與熱輸運的聲子數量減少了,根據經典晶格傳熱理論,聲子數量的減少導致了GaN 結構中聲子群速度的下降[27-28],同時VN-GaN 的費米能級高于VGa-GaN的費米能級。因此同一條件下,完整GaN 的熱導率高于VGa-GaN 和VN-GaN 的熱導率,VGa-GaN 的熱導率高于VN-GaN 的熱導率。

如圖10 所示,是計算得到的完整GaN、VGa-GaN和VN-GaN 在300 K 時的PPR。為實現聲子局部化的概念,將局域狀態定義為參與率小于0.2 的狀態,而大于0.2 的部分設定為離域狀態,這個定義將在整個研究中保持一致。通過圖10 可以發現,在PPR 譜的左端(0～10 THz),聲子處于離域狀態,右端(55～65 THz)的聲子處于局域狀態。對GaN 熱傳輸有主要貢獻的聲子來自離域狀態[17]。其中低頻區域(0～10 THz)的聲子主要是聲學模式,高頻區域(55～65 THz)的聲子則是光學模式。根據圖10 可知,無論是在低頻區域還是在高頻區域,PPR 譜都遵循倒U 形,不同結構GaN 的離域狀態大都落在相同的區域中。這說明空位缺陷對GaN 熱傳輸時PPR 的影響并不大,即PPR并不是導致GaN 熱導率隨著空位缺陷出現而下降的主因。

圖10 300 K 時完整GaN 與含空位缺陷GaN 的PPR 對比Fig.10 PPR of intact GaN compared to GaN with vacancy defects at 300 K

2.4 原子摻雜對GaN 熱導率的影響

如圖11 所示,是完整GaN、Fe 摻雜GaN、Mg 摻雜GaN 和Si 摻雜GaN 的熱導率在300～600 K 溫度范圍內的變化趨勢。根據圖11 可知,當溫度從300 K 變化到600 K 時,四種結構GaN 的熱導率隨著溫度的升高均呈下降趨勢。完整GaN 的熱導率下降速率快于GaN 摻雜結構的熱導率下降速率,Fe 摻雜GaN 熱導率的下降速率快于其他兩種GaN 摻雜結構熱導率的下降速率。這說明含有原子摻雜結構的GaN 熱導率對于溫度依賴性的一般趨勢是隨著溫度的升高而減小。當溫度為300 K 時,完整GaN 的熱導率最大,為43.06 W·m-1·K-1。Fe 摻雜GaN、Mg 摻雜GaN 和Si 摻雜GaN 的熱導率與完整GaN 的熱導率相比會出現較大的下降。Fe 摻雜GaN、Mg 摻雜GaN 和Si 摻雜GaN 的熱導率分別為10.15,8.22 和4.70 W·m-1·K-1,與完整GaN 的熱導率相比分別下降76.43%,80.91%和89.08%。這表明完整GaN 原子經摻雜后,它的熱傳導性能會發生較大的改變。

圖11 完整GaN 與原子摻雜GaN 的熱導率對比Fig.11 Comparison of the thermal conductivity between intact GaN and atomically doped GaN

為了解原子摻雜對GaN 熱傳輸機理的影響,計算了完整GaN 和原子摻雜GaN 結構在300 K 時的相關PDOS,如圖12 (a)～(d) 所示。

圖12 (a) 300 K 時完整GaN 與三種摻雜體系GaN 的PDOS 對比;(b) 300 K 時完整GaN 與三種摻雜體系GaN 的Ga 原子PDOS 對比;(c) 300 K 時完整GaN 與三種摻雜體系GaN 的N 原子PDOS 對比;(d) 300 K 時三種摻雜原子的PDOS 對比Fig.12 (a) Comparison of PDOS between intact GaN and three doping systems of GaN at 300 K;(b) Comparison of Ga atomic PDOS between intact GaN and three doping systems of GaN at 300 K;(c) Comparison of N atomic PDOS between intact GaN and three doping systems of GaN at 300 K;(d) Comparison of PDOS for three doped atoms at 300 K

通過圖12(a)可以發現,300 K 時,這四種結構GaN 的PDOS 曲線具有相似的變化趨勢,在低頻區域(5～10 THz)和高頻區域(20～30 THz)都存在著振動峰。通過分析圖12(a)發現,在5.3 THz 和24.2 THz處存在振動峰,在兩個振動峰處完整GaN 的PDOS 值＞Fe 摻雜GaN 的PDOS 值＞Mg 摻雜GaN 的PDOS 值＞Si 摻雜GaN 的PDOS 值。根據圖12(b)和(c)可知,相同的條件下,完整GaN、Fe 摻雜GaN、Mg 摻雜GaN 和Si 摻雜GaN 中的Ga 原子與N 原子在低頻區域和高頻區域都存在著振動峰。Ga 原子的低頻振動峰值高于高頻振動峰值,而N 原子的高頻振動峰值高于低頻振動峰值。根據圖12(d)可知,Fe 摻雜GaN、Mg摻雜GaN 和Si 摻雜GaN 中摻雜的原子在低頻區域和高頻區域都存在著振動峰,且摻雜原子的低頻振動峰值高于高頻振動峰值。出現上述的變化趨勢是因為當完整GaN 結構中出現原子摻雜時,GaN 熱輸運時的聲子散射現象增強了,參與熱輸運的聲子數量減少,導致了GaN 結構中聲子群速度的下降,從而導致PDOS峰值的下降。當完整GaN 摻雜Fe、Mg 和Si 后,聲子在熱傳輸過程中更容易散射,導致聲子壽命降低。同時Fe 是磁性原子[29],在Fe 取代GaN 中的一部分Ga原子后,形成了比Mg 和Si 摻雜GaN 更穩定的磁性材料。另外由于Mg、Si 摻雜具有較強的聲子非諧性,增強了聲子-聲子散射速率,而Si 摻雜相較于Mg 摻雜對于聲子弛豫時間有更高的要求[30]。因此在相同條件下,完整GaN 的熱導率＞Fe 摻雜GaN 的熱導率＞Mg 摻雜GaN 的熱導率＞Si 摻雜GaN 的熱導率。

如圖13 所示,對完整GaN 和原子摻雜GaN 結構中的聲子性質進行了深入研究,計算了四種結構在300 K 時的PPR。為了實現聲子局部化的概念,同樣將局域狀態定義為參與率小于0.2 的狀態,而大于0.2的部分設定為離域狀態。通過圖13 可以發現,PPR 譜的左端(0～10 THz)聲子處于離域狀態,右端(55～65 THz)的聲子處于局域狀態。對GaN 熱傳輸有主要貢獻的聲子來自離域狀態[17]。其中低頻區域(0～10 THz)聲子主要是聲學模式,高頻區域(55～65 THz)的聲子則是光學模式。無論是在低頻區域還是在高頻區域,PPR 譜都遵循倒U 形,不同結構的GaN 的離域狀態大都落在相同的區域中。這說明完整的GaN 結構中出現原子摻雜時,PPR 對于GaN 的熱導率并沒有產生重要的影響,即PPR 并不是導致原子摻雜GaN 的熱導率小于完整GaN 熱導率的主因。

圖13 300 K 時完整GaN 與原子摻雜GaN 的PPR 對比Fig.13 PPR of intact GaN compared to atomically doped GaN at 300 K

3 結論

本文基于NEMD 方法,利用SW 勢函數研究了空位缺陷和原子摻雜對GaN 熱導率的影響。得出以下結論:

(1) 當GaN 含有空位缺陷結構時,相同溫度下,完整GaN 的熱導率高于含有空位缺陷的GaN 熱導率,且VGa-GaN 的熱導率總高于VN-GaN 的熱導率。

(2) 對于含有原子摻雜的GaN 結構,相同條件下,GaN 摻雜體系的熱導率都低于完整GaN 的熱導率。Fe 摻雜的GaN 熱導率為GaN 三種摻雜結構中最高,Si 摻雜的GaN 熱導率為GaN 三種摻雜結構中最低。

(3) 當溫度從300 K 上升到600 K 時,完整GaN、含有空位缺陷GaN 與原子摻雜結構GaN 的熱導率均會受到溫度的影響,且隨溫度的上升呈下降趨勢。