弋英民，張潼，劉丹

(西安理工大學 自動化與信息工程學院，陜西 西安 710048)

大尺寸、高品質(zhì)硅單晶是制造微細納米級集成電路芯片的關(guān)鍵材料，目前我國全部依賴進口。硅單晶爐作為生長電子級硅單晶材料的主要設(shè)備，其工藝復雜，技術(shù)集成度高，實現(xiàn)困難，在集成電路生產(chǎn)流程中處于首要地位。因此，針對制備單晶硅的各種基礎(chǔ)研究和仿真模擬成為了學者們研究的熱點。單晶硅按生長方法分有直拉法和區(qū)熔法等[1]，直拉法也稱CZ法，目前廣泛應用于晶體制造領(lǐng)域。在CZ法制備單晶硅的過程中，影響晶體質(zhì)量的因素有坩堝內(nèi)的硅熔體流動[2]、熔體內(nèi)雜質(zhì)分布[3]以及晶體生長面彎月面的形狀[4]等。使用計算機模擬晶體的生長環(huán)境在晶體制備的仿真中起到了重要的作用，為了衡量仿真算法的正確性，文獻[5]提出了Wheeler標準問題，Wheeler標準問題是驗證仿真算法正確與否的檢驗標準。

為了能夠模擬晶體的生長情況，一般采用專用軟件或者編程實現(xiàn)。文獻[6]和[7]采用FEMAG-CZ對單晶爐內(nèi)的熱場和氧濃度分布進行了模擬；文獻[8]采用CGSIM對晶體中的微缺陷進行了模擬研究；文獻[9]使用Fortran編程，利用Petrov-Galerkin有限元法實現(xiàn)了不同物性參數(shù)對晶體生長的影響的模擬；文獻[10]也使用通用編程軟件，模擬了晶體的生長過程，并通過了Wheeler標準的驗證。然而，專用軟件需要大量的授權(quán)費用，使用通用軟件編程需要較高的編程功底和理論水平，因此本文研究使用通用的有限元軟件來實現(xiàn)晶體生長相關(guān)的數(shù)值模擬。

本文采用Fluent 6.3，利用有限容積法結(jié)合三階迎風QUICK離散格式模擬了單晶爐中坩堝內(nèi)熔體的流場和熱場，并用Wheeler標準檢驗了算法的有效性。對坩堝直徑為600 mm的熔體進行了自然對流下的流場、熱場以及雜質(zhì)場的數(shù)值模擬。

1 數(shù)值模型和物性參數(shù)

直拉法晶體生長的簡化結(jié)構(gòu)如圖1所示。

圖1 晶體生長結(jié)構(gòu)簡化圖

假設(shè)坩堝和晶體為同軸圓柱，坩堝半徑為Rc，坩堝角速度為Ωc，坩堝內(nèi)熔體高度為H，晶體半徑為Rx(Rx

模擬中假定：①晶體等徑生長；②熔體流動和換熱是軸對稱的；③熔體與晶體界面的自由表面為平面；④忽略熔體自由表面張力和旋轉(zhuǎn)時界面的變化；⑤熔體為軸對稱的不可壓縮的牛頓流體；⑥不考慮磁場；⑦熔體內(nèi)的流動為層流；⑧坩堝和生長界面處的流動滿足無滑移條件，浮升力滿足Bousinesq假設(shè)。除密度外，其它物性參數(shù)均為常數(shù)，熔體的流動和傳熱滿足N-S方程和能量方程。

該模型引入無量綱變量：

(1)

式中帶*標的變量為有量綱變量，(r,z)為質(zhì)點坐標；u、v、w分別為質(zhì)點的軸向速度、徑向速度和周向速度；p為壓力，T為溫度；μ為熔體的粘性系數(shù)；ρ為熔體密度；g為重力加速度。

基于之前的假設(shè)，無量綱的控制方程為：

(2)

(3)

(4)

(5)

(6)

本文所采用熔體模型的邊界條件為：

(7)

相應的，模型中的無量綱參數(shù)定義如下：

(8)

式中，α為熔體高度與半徑之比；β為晶體與熔體半徑之比；Pr為普朗特數(shù)；Gr為格拉斯霍夫數(shù)；Rex和Rec分別是因晶體和坩堝旋轉(zhuǎn)產(chǎn)生的雷諾數(shù)。

仿真中采用的硅熔體的物性參數(shù)如表1所示[10]。

表1 物性參數(shù)表

2 Wheeler標準化問題

在1989年召開的北約晶體生長計算模擬大會上，大會主席Wheeler提出了四個測試問題，即Wheeler標準問題。這些標準問題包括四組情況：A組只考慮晶體旋轉(zhuǎn)；B組考慮晶體和坩堝共同旋轉(zhuǎn)；C組只考慮自然對流，晶體坩堝不旋轉(zhuǎn)；D組考慮自然對流和晶體旋轉(zhuǎn)，坩堝不旋轉(zhuǎn)。每組又分三個不同的算例，Wheeler標準問題的參數(shù)如表2所示。本文計算了表2中的12組算例。

表2 Wheeler標準問題參數(shù)

3 數(shù)值模擬及分析

本文使用Fluent 6.3[11]進行數(shù)值模擬，采用層流模型，壓力和速度的解耦合用Fluent 6.3中的SIMPLE算法，擴散項差分采用二階截差中心差分格式，對流項差分采用三階截差QUICK格式，壓力插值采用PRESTO！格式。對Wheeler標準問題的12個算例進行了數(shù)值模擬。

算例A結(jié)果的流函數(shù)如圖2所示。在這組算例中，熔體內(nèi)部只有一個由晶體旋轉(zhuǎn)產(chǎn)生的一個順時針的渦，熔體沿軸線向上運動，直至晶體邊界，然后流向坩堝側(cè)壁，晶體附近的流體由于粘性的作用和晶體一起旋轉(zhuǎn)，產(chǎn)生從軸心指向坩堝壁的離心力，在離心力的作用下，熔體流向坩堝壁，遇到坩堝壁阻擋，被迫下流。轉(zhuǎn)速將渦流的中心移向坩堝側(cè)壁，偏向右上角，并發(fā)生形變。同時，速度最快的區(qū)域從左上角移至右下角，速度最小的區(qū)域在左下角擴大。

圖2 算例A的流函數(shù)圖

在本組算例中，隨著晶體轉(zhuǎn)動速度的增加，流函數(shù)變大，這說明渦流強度增加。算例A1、A2、A3的流函數(shù)最大值分別為8.79E-5、1.18E-4、1.47E-4。

算例B結(jié)果的流函數(shù)如圖3所示。在算例B中，熔體內(nèi)部出現(xiàn)兩個方向相反的渦流，左上角的順時針渦流是由晶體旋轉(zhuǎn)形成的，右下角的逆時針渦流是由坩堝旋轉(zhuǎn)造成的；隨著晶體旋轉(zhuǎn)速度和坩堝旋轉(zhuǎn)速度的同時增大，雖然晶體轉(zhuǎn)速遠大于坩堝轉(zhuǎn)速，但是右下角的渦流仍主導著整個流場；晶體生長的理想條件是熔體區(qū)域只有一個渦流，顯然，這樣的情況是不利于晶體生長的。

圖3 算例B的流函數(shù)圖

在本組算例中，隨著晶體轉(zhuǎn)動速度和坩堝轉(zhuǎn)動速度的增加，流函數(shù)變大，說明渦流強度顯著增加。算例B1、B2、B3的流函數(shù)最大值分別為8.25E-7、1.86E-6、3.92E-5。加入晶體和坩堝旋轉(zhuǎn)是為了平衡自然對流產(chǎn)生的渦流，但晶體和坩堝的轉(zhuǎn)速過大，產(chǎn)生的渦流強度過大，也對晶體等徑生長帶來擾動，因此晶體和坩堝的轉(zhuǎn)速應控制在一定范圍之內(nèi)。

算例C中只有浮升力引起的自然對流，熔體內(nèi)部由一個逆時針的渦占據(jù)，隨著Gr數(shù)的增大，渦流也在變大，最大速度移向坩堝側(cè)壁，渦流的中心幾乎沒有移動位置。算例C結(jié)果的流函數(shù)圖如圖4所示。Gr數(shù)在108時，層流模型無法得到收斂的解，證明層流模型不適用，Gr數(shù)在108以上的數(shù)量級時，應該使用湍流模型。

圖4 算例C的流函數(shù)圖

算例D結(jié)果的流函數(shù)圖如圖5所示。算例D相比算例C加上了晶體旋轉(zhuǎn)，結(jié)果加上晶體旋轉(zhuǎn)的自然對流和只有自然對流的情況，流函數(shù)不管是在形狀上，還是在數(shù)值上都沒有什么變化，和C組的幾乎一樣。

圖5 算例D的流函數(shù)圖

通過Wheeler標準問題檢驗后，證明此方法可用，于是，本文將此方法應用于大直徑硅晶體，對坩堝半徑為300 mm內(nèi)的硅熔體，進行自然對流下的流場、熱場和氧濃度的數(shù)值模擬。熔體高度與半徑之比α=0.6；晶體與熔體半徑之比β=0.4。

圖6為坩堝半徑為300 mm時，硅熔體在自然對流下的流函數(shù)圖?？梢钥闯?，流函數(shù)的基本形狀和層流模型的基本相同，是充滿整個熔體區(qū)域的一個逆時針大渦，但是渦心更偏向坩堝側(cè)壁，此渦是由于坩堝底部和側(cè)壁的溫度比熔體高，坩堝壁附近區(qū)域的熔體受熱后，密度變小，在浮升力的作用下，沿著側(cè)壁上升，然后在自由表面和生長界面處受冷，密度變大，下沉至坩堝底部，從而形成完整的渦流。

圖6 硅熔體流函數(shù)圖

圖7為坩堝半徑為300 mm時，硅熔體在自然對流下的等溫圖。低缺陷晶體生長的溫度條件是在整個晶體生長過程中固液交界面保持平坦，這要求等溫線近似垂直于生長方向，即等溫線比較平坦。圖7中，靠近晶體生長面的等溫線比較平坦，隨著遠離晶體生長界面，等溫線變得更為陡峭。自然對流的湍流會加劇溫度的波動，而溫度的波動是產(chǎn)生生長條紋的重要原因，因此要獲得高品質(zhì)單晶硅，控制其溫度波動是關(guān)鍵。

圖7 硅熔體等溫圖

圖8顯示的是在自然對流作用下坩堝半徑為300 mm的硅熔體中的氧原子濃度分布圖?？梢钥闯觯谌垠w內(nèi)部的大部分區(qū)域，氧原子濃度都很小。并且，由于渦流在自然對流的作用下，沿坩堝壁經(jīng)由自由表面再流向?qū)ΨQ軸中心，因此一部分氧雜質(zhì)也會從自由表面揮發(fā)出去，在圖中就可以看出大部分的氧雜質(zhì)聚集在坩堝底部和側(cè)壁。如果熔體內(nèi)的流動加強，那么氧原子雜質(zhì)有可能會來不及揮發(fā)而被帶入熔體，從而形成熔體內(nèi)的雜質(zhì)，因此抑制自然對流的流動強度對降低熔體內(nèi)的氧雜質(zhì)至關(guān)重要。

圖8 硅熔體氧濃度圖

根據(jù)本研究結(jié)論，通過在直徑300 mm硅單晶爐修改拉晶的SOP表有關(guān)參數(shù)，進行不同SOP溫度參數(shù)下的拉晶實驗如圖9所示。

圖9 在不同拉晶SOP參數(shù)下的晶棒

由圖9可以看出，在晶體生長界面處溫度波動大，會導致生長處的晶體無法成晶，甚至斷裂。而當控制晶體生長界面(彎月面)溫度為硅熔體成晶溫度1 420℃±0.1℃時，可以長出完美的硅單晶。

4 結(jié) 論

本文通過Fluent6.3對Wheeler標準問題中的12個算例進行了數(shù)值計算，結(jié)果表明，使用Fluent計算的結(jié)果與文獻[10]中的基本一致。本文的仿真實驗結(jié)果表明，對于直徑600 mm坩堝必須控制合適的晶轉(zhuǎn)和堝轉(zhuǎn)范圍，以保證晶體生長，控制溫度的波動可以抑制生長條紋，抑制自然對流的可以降低熔體內(nèi)的氧雜質(zhì)含量。本研究下一步將在自然對流下增加晶體和坩堝旋轉(zhuǎn)，并增加間壁對氧濃度的抑制方面的研究。

參考文獻：

[1] 蔣榮華,肖順珍. 半導體硅材料的進展與發(fā)展趨勢[J].四川有色金屬, 2000,3:1-7.

Jiang Ronghua,Xiao Shunzhen. Progress and development trend of semiconductor Si materials[J]. Sichuan Nonferrous Metals,2000,3:1-7.

[2] 宇慧平,隋允康,張峰翊,等. 紊流模型模擬分析旋轉(zhuǎn)對提拉大直徑單晶硅的影響[J]. 人工晶體學報, 2004,33(5): 835-840.

Yu Huiping,Sui Yunkang,Zhang Fengyi, et al. Numerical simulation on the effect of rotation in a czochralski silicon crystal growth with a turbulence model[J]. Journal of Synthetic Crystal, 2004,33(5): 835-840.

[3] Ren Bingyan,Zhao Long,Zhao Xiuling, et al. Effects of argon gas flow rate and guide shell on oxygen concentration in czochralski silicon growth [J]. Rare Metals, 2006, 25(1): 7-10.

[4] 滕冉,戴小林,肖清華,等. 大直徑硅單晶生長過程中固/液界面形狀及熔體流動的數(shù)值分析[J]. 人工晶體學報,2013,42(4): 611-615.

Teng Ran,Dai Xiaolin,Xiao Qinghua, et al. Numerical simulation of melt/crystal interface and melt flow during large diameter single silicon crystal growth[J]. Journal of Synthetic Crystal, 2013,42(4): 611-615.

[5] Wheeler A A. Four test problems for the numerical simulation of flow in czochraski crystal growth [J]. Journal of Crystal Growth,1990,102(4):691-695.

[6] 王慶,張婷曼. 基于FEMAG-CZ熱場仿真軟件的單晶爐熱場分析[J]. 計算機光盤軟件與應用, 2011,17: 65-66.

Wang Qing,Zhang Tingman. Single crystal furnace thermal field analysis on FEMAG-CZ thermal field simulation software[J]. Computer CD Software and Applications, 2011,17: 65-66.

[7] 常麟,周旗鋼,戴小林,等. CUSP磁場對直拉硅單晶氧濃度分布影響的數(shù)值模擬[J].稀有金屬,2011,35(6): 909-915.

Chang Lin,Zhou Qigang,Dai Xiaolin, et al. Numerical simulation of CUSP magnetic field on oxygen concentration distribution in CZ-Si crystal growth[J]. Chinese Journal of Rare Metals, 2011,35(6): 909-915.

[8] 曾慶凱,關(guān)小軍,潘忠奔,等. Ф400mm直拉硅單晶生長過程中氧濃度對微缺陷影響的數(shù)值模擬[J]. 人工晶體學報,2011,40(5):1150-1156.

Zeng Qingkai,Guan Xiaojun,Pan Zhongben, et al. Simulation of microdefect with different oxygen concentration in Ф400 mm CZ silicon crystal growth[J]. Journal of Synthetic Crystal,2011,40(5): 1150-1156.

[9] 李友榮,魏東海,余長軍,等. 物性參數(shù)對硅單晶Czochralski生長過程的影響[J]. 熱科學與技術(shù),2006,5(4): 351-355.

Li Yourong,Wei Donghai,Yu Changjun, et al. Effects of thermophysical properties on silicon single crystal czochraski growth processes[J]. Journal of Thermal Science and Technology,2006,5(4): 351-355.

[10] 宇慧平. MCZ大直徑單晶體生長的數(shù)值模擬及控制參數(shù)優(yōu)化[D]. 北京:北京工業(yè)大學,2005.

Yu Huiping. The numerical simulation of the large diameter MCZ crystal growth and the optimization of the operational parameters[D]. Beijing:Beijing University of Technology,2005.

[11] 李國棟,哈岸英,鐘小彥,等. 基于FLUENT 的滲流場數(shù)值模擬分析[J]. 西安理工大學學報,2011,27(3): 317-320.

Li Gudong,Ha Anying,Zhong Xiaoyan, et al. Numerical simulation of seepage field based on FLUENT[J]. Journal of Xi’an University of Technology, 2011,27(3): 317-320.