大尺寸碳化硅晶體生長熱-質輸運過程建模及數值仿真

盧嘉錚,張 輝,鄭麗麗,馬 遠

(1.清華大學航天航空學院,北京 100084;2.清華大學工程物理系,北京 100084;3.中電化合物半導體有限公司,寧波 315336)

0 引 言

碳化硅(SiC)單晶襯底是制造新一代功率器件、射頻器件的重要基礎材料,但當前SiC襯底成本偏高,擴大SiC晶體直徑并提高晶體品質可增加晶圓利用率,達到降低器件成本目的。目前國內業界主要采用中頻感應加熱的物理氣相傳輸(physical vapor transport, PVT)法生產100 mm直徑SiC單晶[1],150 mm單晶具備批量供應能力[2-3],200 mm襯底已研發成功[3]。國外公司已大量供應200 mm晶片[4-5]。另一方面,由于電阻加熱方式可有效控制坩堝內局部溫場,電阻加熱式PVT系統逐漸嶄露頭角[6-7]。

PVT法的關鍵技術之一是設計合理的熱-質輸運環境,為生長低缺陷晶體提供穩定的熱場和均勻的流場。實驗表明:晶體生長界面中心低溫,晶體表面微凸,晶體質量較高[8];預燒原料可增大粉料顆粒度[9]、減小原料在長晶過程中的形狀變化,從而穩定坩堝內溫場[10];原料摻鈰(CeO2)能有效控制多型生長[11]。PVT坩堝是工作在2 000 ℃高溫的封閉結構,實驗手段難以監測坩堝內部熱-質輸運過程。計算機仿真是開展研究的重要辦法,在模擬原料多孔結構變化[12-13]、熱場設計[14]、預測成核點[15]、計算晶體熱應力和位錯密度[16-17],以及合成摻礬原料[18]等方面提供了理論指導。但針對電阻爐生長大直徑晶體的全過程模擬研究較少,原料消耗、晶體形貌變化等過程的數學模型不清楚,晶體生長過程中特征量之間的相互作用不夠清晰。

本團隊前期研究了PVT法SiC單晶生長熱場設計[14],為本研究提供了有力支撐。本研究針對電阻加熱式PVT爐生長150 mm的SiC單晶開展熱-質輸運過程建模和晶體生長過程的數值模擬研究。首先建立碳化硅原料分解及其多孔結構演變、氣相組分輸運、能量輸運、晶體生長的物理和數學模型,再研究不同原料溫度分布(側面高溫、底部高溫)對晶體生長形貌變化的影響規律,最后厘清晶體生長形貌變化與原料消耗、溫場變化等過程之間的關系。

1 熱-質輸運數學模型

PVT工藝生長SiC晶體的主要物理過程是:中頻感應線圈或電阻加熱器產生熱源,在坩堝內建立軸向和徑向溫度梯度,裝于坩堝下部的原料被加熱后分解,熱解產生的氣相組分在溫度梯度的驅動下被輸運至坩堝頂部的低溫籽晶面發生再結晶生成單晶體。如圖1所示,氣相組分輸運先后經過原料堆、坩堝上部的空腔區、籽晶面,每個區域內物理過程不同,因此將坩堝內部劃分為相應三個區域,并依次建立數學模型。電阻加熱式PVT工藝的詳細物理過程和熱力學過程參見前期工作[14]。

1.1 原料區控制方程

1.1.1 原料熱解、碳化和再結晶

原料區是坩堝內裝載SiC原料粉末的區域(見圖1),一般處于2 500～2 600 K,其最高溫度點位于側面或底部,最低溫度點一般在頂部。

高溫區域的SiC粉末受熱分解生成Si、SiC2、Si2C、SiC等氣相組分和疏松石墨,氣相組分在溫度梯度的驅動下輸運至低溫的SiC粉末表面和籽晶面再結晶,即:高溫的原料處于欠飽和熱力學狀態,被消耗;低溫的原料為過飽和狀態,SiC質量增加。

假設:

1)SiC原料為直徑0.24 mm的球形顆粒,SiC顆粒石墨化后直徑不變[19];

2)原料區為多孔介質區域;

3)有限體積法的網格尺寸大于顆粒尺寸,同一網格體內的顆粒處于同樣的熱力學狀態;

4)球形顆粒呈緊密交錯排列;

5)氣相組分中Si、C摩爾比為1∶1,即總的SiC氣相組分摩爾流量JA=2JSiC2。

在顆粒表面,根據上述物理過程,可用Hertz-Knudsen公式描述原料粉末表面產生的熱解氣相流量:

(1)

(2)

(3)

式中:ΔA為控制體(網格)內所有SiC顆粒的表面積;ΔV為控制體(網格)體積。ΔA由下式計算:

ΔA=αNpAp

(4)

式中:Np為控制體內顆粒個數;Ap為單個顆粒的表面積;α為修正系數。根據假設4),控制體內顆粒個數為:

(5)

式中:ε0為原料的初始孔隙率;r為顆粒原始半徑。由SiC質量守恒與C質量守恒可得:

(6)

與

(7)

ε=1-γSiC-γC

(8)

定義原料消耗程度dg:

(9)

PSi∶PSi2C∶PSiC2≈4∶1∶2

(10)

根據上述組分所涉及反應的化學計量比和假設5),可近似認為:

(11)

(12)

通過式(1)～(7)可知,當顆粒表面SiC氣相組分的分壓P大于此處的平衡壓力P*時,顆粒處于過飽和狀態,氣相組分在顆粒表面再結晶,SiC體積分數增加;反之SiC顆粒被消耗,SiC體積分數降低,C體積分數增加。由于假設5),異相反應的質量不守恒,這是合理的。原因是:雖然SiC晶體生長一般處于富硅氣相,但晶體生長所需Si與C原子數是1∶1,在宏觀傳質層面,可忽略多余氣相Si對晶體生長的影響。

1.1.2 氣相組分輸運

根據假設1)和2),氣相組分在原料中的輸運可視作氣體在球形顆粒組成的多孔介質中流動和擴散。可用修正后的連續性方程、動量方程和組分方程描述:

(13)

(14)

(15)

式中:ρg為氣相平均密度;μ為黏度系數;K為Darcy定律滲透系數(permeability);Deff為有效擴散系數。分別由下式計算:

(16)

(17)

式中:τ為多孔介質的迂曲率[21],其大小與多孔介質基體形狀、孔隙率有關;Dij為氣體的二元擴散系數。

1.1.3 能量守恒

SiC原料粉末中存在能量輸運,顆粒間通過熱輻射、熱對流和固相導熱進行熱量傳遞(見圖1),原料區能量守恒方程如下:

(18)

式中:(ρcp)eff是等效熱容,為SiC原料、石墨和氣相組分熱容的體積平均值;keff是有效熱導系數,為SiC顆粒導熱、疏松石墨導熱、氣體導熱和孔隙間輻射傳熱的體積平均值。分別由下式計算:

(ρcp)eff=γSiC(ρcp)SiC+γC(ρcp)C+ε(ρcp)g

(19)

(20)

式中:εp為SiC顆粒表面的輻射發射系數;(ρcp)SiC、(ρcp)C和(ρcp)g分別為SiC原料、石墨和氣相組分的熱容;kSiC、kC和kg分別為SiC原料、石墨和氣相組分的熱導率。

1.2 生長室中的熱-質輸運

氣相組分從原料表面逸出進入坩堝上部的生長室空腔,在此區域內,氣相組分自由流動擴散。忽略氣相輻射,坩堝壁、晶體外表面和原料區上表面之間存在輻射傳熱。坩堝外的氣相區域也存在相似物理過程,如電阻加熱器與坩堝外壁間、坩堝與保溫棉間的輻射傳熱,坩堝外氬氣的流動擴散等。首先,坩堝內除原料區,對氣相有質量守恒、動量守恒和組分方程:

(21)

(22)

(23)

對包含保溫棉、坩堝和加熱器的整個系統(除原料區)有能量守恒:

(24)

(25)

式中:δij是Kronecker delta;Fji為視角系數,表示從j面發射的能量被i面攔截的份額,用下式計算:

(26)

式中:A為微元面i的面積;R為i面與j面心連線長度;θi為R與i面外法向量夾角。參與輻射的面包括晶體表面、生長室坩堝內壁、原料區頂部邊界在內的所有氣-固邊界。

1.3 晶體生長界面

在晶體生長界面上存在Stefan流,組分邊界條件采用Hertz-Knudsen模型,對SiC氣相組分和氬氣分別有:

(27)

JAr=0

(28)

式中:J為摩爾流;ci為組分i的摩爾濃度;vn為法向速度。忽略晶體的徑向生長,晶體在法向的生長速率G為:

(29)

式中:Mcrystal為SiC晶體的摩爾質量;ρcrystal為SiC晶體密度。采用動網格技術,在每個時間步前根據上一步計算得到的G值調整晶體生長界面的節點位置,模擬晶體生長界面形狀變化和晶體增厚,并重構長晶界面兩側的體網格。

2 數值模擬研究

本節展開對電阻加熱式150 mm晶體生長過程的數值模擬研究,探究不同原料溫度分布條件對晶體生長的影響規律,梳理晶體生長與原料和熱場變化之間的耦合關系。

2.1 算例說明

2.1.1 幾何模型與計算說明

計算所用幾何模型是對前期工作[14]的改進,如圖2所示,系統結構呈圓柱形,半徑約570 mm,高約1 600 mm。布置側面和底部共2個電阻式加熱器。

圖2 計算域幾何結構示意圖Fig.2 Schematic diagram of computational domain geometry

為簡化計算,采用軸對稱計算域,對稱軸邊界條件為溫度0梯度、壓力0梯度、軸向速度0梯度和徑向速度為0,其余外邊界為320 K定溫、速度為0、壓力0梯度。坩堝內氬氣壓力300 Pa。采用商業軟件ICEM對計算域劃分非結構網格,網格尺度為0.5 mm。采用商業軟件Fluent對方程進行離散和求解,使用Dynamic Mesh對晶體生長界面網格及其附近體網格進行重構。物性參數和幾何尺寸參見表1。

表1 物性和參數[22-24]Table 1 Material properties and parameters for computation[22-24]

2.1.2 算例設置

150 mm晶體生長爐尺寸較大,原料內溫度、氣流在徑向上的不均勻分布影響晶體生長形貌。因此設置2個算例,模擬在典型原料溫度分布下的晶體生長過程,探究原料溫度對晶體生長形貌變化的影響規律,分析坩堝內關鍵參量間的耦合關系。算例1僅打開側面加熱器,算例2同時打開側面和底部加熱器,算例模擬時長90 h,算例設計說明如表2所示。

表2 算例設計說明Table 2 Case description

2.2 原料溫度分布對晶體形貌影響

2.2.1 側面高溫

算例1模擬的工況是只用側面加熱器對坩堝進行加熱,初始時刻坩堝內溫場如圖3(a)所示,原料最高溫度約2 556 K,位于側面中部。原料頂部和底部是低溫區,溫度2 475～2 485 K。晶體生長界面中心點溫度約2 410 K,邊緣約2 440 K。圖3(b)展示初始時刻坩堝內的流場,高溫區原料熱解生成的氣相產物在溫度梯度的驅動下被輸運至低溫的原料底部、原料頂部和晶體生長界面。對比流場和多孔結構變化,與文獻[19,25]符合度較高。

圖3 算例1熱/流場結果圖。初始時刻坩堝內溫度云圖(單位:K)(a)和流線與速度矢量圖(b);30 h(c)、60 h(d)和90 h(e)的晶體形貌、流場、原料含量和孔隙率,左半圖展示流線與SiC原料體積變化,右半圖呈現速度矢量與孔隙率(初始值0.65)Fig.3 Heat/flow field results of case 1. At the initial moment, the temperature distribution (unit: K) (a), and streamline and velocity inside the crucible (b); the crystal morphology, flow field, feedstock content and porosity after 30 h (c), 60 h (d) and 90 h (e), the left half of the graph shows the streamlines and the volume change of the SiC source material, and the right half of the graph presents the velocity vector and porosity (initial value 0.65)

從圖3可知,在原料表面附近,速度大小沿徑向(A至B)逐漸變大,這是由于原料側面(B至D)溫度梯度大,氣流驅動力大,A至C點溫度梯度小,驅動力小,此現象在后續過程中進一步加劇。SiC氣相組分在A點附近再結晶,該區域SiC顆粒體積在90 h增幅約30%,孔隙率從0.65降至0.55,氣流通道縮小。而D至B點區域原料不斷消耗,孔隙率增加,氣流通道逐漸打開,導致A點附近氣流速度不斷減小,B點附近氣流速度不斷增大,氣流不均勻性凸顯。原料底部再結晶程度比頂部更嚴重,90 h后該區域SiC固相體積增加了70%,造成氣源浪費,原料使用率低。

徑向分布不均勻的氣流量導致晶體中心區域的生長速率比5/8半徑處慢,由圖3(d)可見生長界面呈“W”形。根據圖4所示,平均晶體生長速率先增加后降低,拐點位于40 h附近;A點的生長速率從0.22 mm/h不斷降低至0.15 mm/h;前期B點生長速率高于A點,先從0.24 mm/h增至0.25 mm/h,在40 h后迅速降低,最終減至0。因此,中前期B點晶體厚度比A點大,即界面呈“W”形狀,后期A點晶體厚度增加較快,最終A/B點晶體厚度相近,生長界面形狀平整。

圖4 算例1晶體生長速率曲線。(a)晶體生長界面上的A和B點(A點位于界面中心,B點在處);(b)A/B點的晶體生長速率(GA、GB)、晶體厚度(TA、TB)和長晶界面平均生長速率隨時間變化圖Fig.4 Crystal growth rate curve of case 1. The crystal growth interface taking point A at the center of the interface and point B at 5/8 R; (b) plots of the crystal growth rate (GA, GB), crystal thickness (TA, TB) at point A/B and average growth rate as a function of time

綜上,對于原料內最高溫度點位于側面的情況,原料側面被消耗使氣流通道打開,同時原料頂部再結晶使氣流通道縮小,造成貼壁附近氣流量大、中心區域氣流量小,導致中前期籽晶面中心的生長速率低于邊緣的生長速率。后期晶體中心生長速率大于晶體邊緣生長速率,中心的凹陷被補齊,此現象是PVT法生長晶體的共性問題,在后文2.3小節中會進一步討論。側面氣流通道的擴張為晶體生長提供較穩定的氣源,長晶界面平均生長速率變化較平緩。

2.2.2 底部高溫

算例2用側面加熱器和底部加熱器同時加熱坩堝。通過圖5可知,初始時刻,原料底部是高溫,頂部是低溫,且原料區徑向溫度梯度小。在此原料溫度分布情況下,底部的SiC原料首先熱解,氣相組分沿負溫度梯度方向輸運至籽晶面,氣流速度在徑向上的變化小,氣流均勻。值得注意,算例2與算例1的總加熱功率相等,但算例2中坩堝內溫度更高,因為本研究使用的幾何模型中,底部加熱器離坩堝更近,且系統底部的保溫棉更厚,使坩堝底部受熱更充分。

由圖5可知,前期晶體從中心區域開始生長,長晶界面中心凸起,邊緣平緩。中后期長晶界面形狀幾乎無變化,各處生長速率相近,晶體表面近似做“平移”運動。該現象與算例1后期的晶體中心區域補長類似。

原料底部的SiC原料持續消耗,無再結晶現象,氣相組分從原料區下方流向上方的低溫區,并在原料頂部低溫區域再結晶。氣相組分均勻流出原料表面,氣流速度在徑向上變化小。但隨時間增加,未熱解原料減少,且原料頂部氣流通道在90 h縮小60%以上,導致氣流速度逐漸減小,氣流量降低。理論上算例2的原料使用率更高,晶體更厚,晶體形貌更符合需求。

圖5 算例2熱/流場結果圖。初始時刻的坩堝內溫場(a)和流場(b);30 h(c)、60 h(d)和90 h(e)的形貌,左半圖展示坩堝內SiC原料體積變化量、氣流的流線,右半圖呈現孔隙率變化量、速度矢量Fig.5 Heat/flow field results of case 2. The temperature field (a) and the flow field (b) in the crucible at the initial moment of case 2. Crystal morphology after 30 h (c), 60 h (d) and 90 h (e), the left half of the figures show the volume change of the SiC source material in the crucible and the streamline of the gas flow, and the right half of the figures show the porosity change and the velocity vector

2.3 晶體生長、原料演變與熱場、流場變化等過程間的耦合關系

晶體生長界面的溫度分布是影響界面形狀的重要因素,由2.2小節已知晶體形貌受原料溫度分布、氣流均勻性影響,長晶界面形狀動態變化是多因素綜合作用的結果。本節分別對原料和晶體等特征區域溫度、形狀結構變化帶來的影響展開討論。

2.3.1 原料區溫度和孔隙結構變化的影響

原料區溫度是決定晶體生長界面溫度的主要因素之一,也是晶體生長的氣源。電阻加熱器通過熱輻射加熱坩堝外壁,再由坩堝導熱將熱量傳遞給SiC原料,因此原料內最高溫度點總是貼近坩堝壁。高溫區的原料熱解后不斷生成疏松石墨,形成較大熱阻,以算例2為例,原料區平均溫度從初始時刻的2 551 K降低至90 h時的2 522 K,如圖6(d)所示。原料溫度下降導致晶體生長界面溫度降低,加劇了原料低溫區的再結晶速度,使原料頂部孔隙率降低、氣流通道縮小(見圖6(e)),抵達晶體生長界面的氣流速度減小(見圖5)。

原料頂部再結晶使SiC固相體積分數從0.35增至0.57,增幅60%。頂部區域導熱性能增強,A點附近溫度梯度減小。原料底部a點附近由于SiC消耗生成疏松石墨,導熱性能減弱,該處溫度梯度變大。但由于加熱器功率、坩堝外環境無變化,a點處溫度穩定在2 598 K左右。若用類似算例2溫場進行實際長晶試驗,上述現象可能造成測量誤差,即坩堝底部的溫度測量值變化較小,但原料內實際溫度可能已低于設計值。

2.3.2 晶體生長界面溫度與形貌動態變化

長晶界面與原料表面存在較強輻射傳熱,面對面輻射傳熱強度主要受各表面溫度和面與面距離影響。晶體增厚,晶體表面與原料表面間距縮短,晶體生長界面溫度與原料表面溫度、晶體厚度成正相關。以算例2結果為例(見圖6(d)),在0～20 h,原料溫度略降(2 551～2 548 K),長晶界面均溫升高(2 438～2 445 K),即該時段內,晶體增厚對長晶界面溫度影響更大。20～40 h,原料溫度持續降低,二者對晶體表面溫度的影響相互抵消,長晶界面均溫保持在2 443 K左右。40 h后,原料溫度降低導致長晶界面均溫從40 h時的2 443 K降至90 h時的2 424 K。

圖6 算例2典型時刻溫場和原料體積分數變化。30 h(a)、60 h(b)和90 h(c)坩堝內溫度分布;(d)原料區平均溫度(Tc)、晶體生長界面平均溫度(Ts)和原料底部a點(見圖(a))溫度(Tc_max)隨時間變化情況;(e)原料頂部A點(見圖(a))SiC原料的體積分數、孔隙率隨時間變化圖Fig.6 Typical temperature field and source material volume fraction changes in case 2. The temperature distributions in the crucible at 30 h (a), 60 h (b) and 90 h (c); (d) average temperature in the feed zone (Tc), average temperature of the growth interface (Ts) and the temperature (Tc_max) at point ‘a’ (see Fig(a)) changes with time; (e) volume fraction and porosity of the SiC source material at point A (see Fig(a)) change with time

長晶界面形狀主要由界面溫度控制[8]。如圖7所示,晶體從中心開始生長,中心生長速率大,邊緣生長速率小,長晶界面向等溫線形狀趨近。長晶界面徑向溫度梯度不斷減小,90 h晶體表面形狀與2 426 K等溫線近乎平行。晶體附近等溫線形狀是徑向溫度梯度的體現,主要由系統頂部的散熱與加熱設計決定。

長晶界面形狀在一定程度上受氣流均勻性的影響。如圖3、圖4所示,算例1中氣流不均勻,5/8半徑處生長速率大于晶體中心,晶體生長界面先呈“W”形。后期溫度變化調節生長速率,中心厚度被補上,90 h界面平坦,已消除“W”形。若繼續模擬,晶體形貌將與算例2類似。

原料區為晶體生長提供氣源。一是原料溫度高則熱解速率快,晶體生長速率大;二是剩余SiC原料越多,能提供的氣源越多,晶體生長越快。由圖8可知,晶體生長速率與原料區平均溫度、剩余SiC原料質量接近成正比。因此在估算生長速率時,應考慮原料溫度和剩余質量帶來的影響。

圖7 算例2長晶界面上的氣相組分濃度與生長速率變化情況。(a)0、20、40、60和80 h晶體生長界面溫度分布;(b) 20、40、60和80 h的晶體生長界面形狀,與80/90 h時的2 426 K等溫線形狀;(c)晶體生長界面中心點A與邊緣點B位置;(d)A/B點晶體生長速率(GA、GB)、平均生長速率溫度(TA、TB)隨時間變化折圖線;(e)A/B點SiC2氣相組分濃度(cA、cB)和平衡濃度隨時間變化線圖(平衡濃度由平衡分壓換算獲得,晶體生長速率與SiC2濃度和平衡濃度的差值Δc成正比)Fig.7 Variations of gas species concentration and growth rate at the growth interface in case 2. (a) Temperature distribution at the interface of crystal growth at 0, 20, 40, 60 and 80 h; (b) shape of the growth front at 20, 40, 60 and 80 h, and the shape of the isotherm of 2 426 K; (c) position of the center point A and the edge point B on the crystal growth interface; (d) growth rate (GA, GB) at point A/B, the average growth rate and the temperature(TA, TB) versus time; (e) species concentration (cA, cB) of SiC2 and equilibrium concentrations versus time at point A/B. Equilibrium concentration is obtained by conversion of equilibrium pressure, the crystal growth rate is proportional to the difference (Δc) between the SiC2 concentration and the equilibrium concentration

圖8 原料溫度及其剩余量與長晶速率關系。(a)算例1、2晶體平均生長速率(G1、G2)隨原料平均溫度的變化;(b)平均生長速率隨剩余SiC原料質量的變化Fig.8 Graph of source temperature and its residual amount versus growth rate. Variation of the average crystal growth rate (G1, G2) with the average temperature of source materials (a) and quality of the remaining SiC source materials (b)

3 結 論

對150 mm晶體生長全過程的數值仿真結果表明,本研究提出的電阻加熱PVT法SiC晶體生長熱-質輸運模型能較準確地模擬原料消耗、再結晶、晶體生長界面形貌和熱場變化,流場和晶體生長的仿真趨勢與文獻吻合。若進一步利用實驗數據對此模型進行修正,能廣泛用于預測電阻加熱式PVT爐中SiC晶體生長形貌、原料和熱場變化情況。

僅當側面加熱器工作時,原料區徑向溫度梯度較大使氣流不均勻,氣相組分主要沿坩堝壁附近流出原料區,導致長晶界面四周的過飽和度比界面中心大,晶體表面生長成“W”形。但生長后期,晶體中心生長速率變大,凹陷被補上。當底部、側面加熱器同時開啟時,能獲得軸向溫梯大、徑向溫梯小的原料區溫度分布狀態,有利于氣相組分均勻流動,晶體生長呈表面微凸的理想形狀。

晶體表面形狀主要受表面的溫度控制,晶體表面上各處溫度變化引起氣相組分平衡濃度變化,進而改變組分的過飽和度,影響生長速率。長晶界面的低溫區域逐步升溫,生長速率從快到慢,高溫區域逐漸降溫,生長速率從慢到快,晶體生長界面最終形成等溫線形狀,此時各處生長速率相近。工程實際中需估算晶體生長速率,晶體平均生長速率大小和原料區溫度、剩余原料質量均成正相關。