Modeling of the Growth Rate during Top Seeded Solution Growth of

Dec 5, 2011 - A numerical model of the SiC solution growth process has been developed, with the main aim to give quantitative outcomes in addition to ...
0 downloads 0 Views 4MB Size
Article pubs.acs.org/crystal

Modeling of the Growth Rate during Top Seeded Solution Growth of SiC Using Pure Silicon as a Solvent Julien Lefebure, Jean-Marc Dedulle, Thierry Ouisse, and Didier Chaussende* LMGP, CNRS UMR 5628, Grenoble INP, Minatec, 3 parvis Louis Neel, BP257, 38016 Grenoble, France ABSTRACT: To assist the development of high quality single crystalline SiC ingot using the top seeded solution growth process, we have implemented a numerical model with the aim of giving quantitative outcomes in addition to qualitative information. The major role of the convection patterns on the carbon flux is demonstrated. We also evidence that the carbon solubility in liquid silicon is the actual limiting parameter of the SiC solution growth process. A good agreement between computed and experimental growth rates is obtained as a function of temperature, making simulation an adapted predictive tool for the further development of the process.



INTRODUCTION Ten years ago, top seeded solution growth (TSSG) was reconsidered as a possible process for the growth of large size and high quality single crystalline ingots of silicon carbide.1 Despite many attempts (see, for instance, Chaussende et al.2), it is only recently that major achievements have been reported with the real demonstration of a 2 in. diameter single crystalline 6H-SiC ingot, having 10 mm in length.3 The well-known problem of the low growth rate related to solution growth has even been ruled out, as a rate of 1 mm/h could be reached, using high temperature together with chromium−silicon solvent.4 This latter result raises many questions about a technically realistic long lasting process under such operating conditions, but it also brought a new breath of optimism. Similarly to sublimation growth, an effective breakthrough could be expected only if the simulation tool was able to open step by step the black box of the process and to gain insight into the intricate mixture of phenomena. Through a coupled modeling of heat transfers and fluid dynamics, we have demonstrated in a previous work that the different convective flow contributions have a very important effect on the stabilization/destabilization of the growth front.5 For example, we have pointed out the importance of controlling both the electromagnetic and the Marangoni convections, which can easily give rise to turbulent flows in specific areas of the crucible. Very recently, Mercier et al.6 have emphasized through modeling the benefit of using an additional low frequency magnetic field to improve the carbon transport and homogenization of the liquid. After this first phenomenological approach, the present work aims at developing a much more robust and complete numerical model, which could be used on a quantitative footing, i.e. able to predict the growth rate. Note that, to develop the model, we have drawn a specific crucible with very small temperature gradients in order to reduce the strong convective flows as much as possible. The resulting low growth rates obtained are thus far from being optimized and only constitute a first step toward optimization. Finally, the © 2011 American Chemical Society

validity of the model will be discussed with respect to the corresponding experiments.



EXPERIMENTAL SETUP

The experimental setup is a modified Metal Research Czochralski puller with medium frequency induction heating (around 60 kHz). Details about the growth chamber and crucible elements have already been given elsewhere.5 Briefly, pure liquid silicon is used as solvent and placed in a high density graphite crucible, which acts both as a container and as the carbon source. The liquid volume is typically 3 cm in diameter and 1.5 cm high. We use top seeded solution growth (TSSG) geometry, meaning that a SiC seed of 10 mm in diameter is mounted on a graphite rod and dipped into the liquid. The crucible is presented is Figure 1. The temperature is measured with an optical pyrometer on the top of the inner crucible.

Figure 1. Schematic representation of the TSSG crucible. The external graphite crucible (called furnace) is directly heated by induction. The inner crucible (called crucible) acts both as a container for the liquid and as the carbon source. It is mainly heated by radiation from the furnace. Note that the red parallelepiped corresponds to the liquid region for which the simulations are represented in the paper. The reactor atmosphere is argon with a slight overpressure of 0.5 bar. Three experiments have been performed for three growth Received: October 10, 2011 Revised: November 17, 2011 Published: December 5, 2011 909

dx.doi.org/10.1021/cg201343w | Cryst. Growth Des. 2012, 12, 909−913

Crystal Growth & Design

Article

temperatures, 1700, 1800, and 1900 °C, in order to compare with the calculated results. The growth procedure consists of 60 min of temperature ramp up to the growth temperature. Once the growth temperature is reached, we wait 15 min as a thermal stabilization step. Then, the SiC substrates are dipped at a depth of 3 mm in the liquid for 8 h with a constant rotation rate of 40 rpm. On-axis (0001) 6H-SiC substrates are stuck on the graphite rod using a graphite glue. At the end of the growth, the crystals are pulled out of the liquid and naturally cooled down. As the growth surfaces are mirror-like, they are perfectly dewetted by the liquid and can be directly observed after removal from the pulling rod. Experimental growth rates are measured by optical microscopy on cross sections, with the formed epilayers having a different color compared to the substrates due to their higher nitrogen content. For these experiments, the carbon polarity is used for the growth, as it is less sensitive to polytype switch during epitaxial growth.7

order to compute the SiC growth rate in the diffusion boundary layer. To evidence the dissolution and crystallization areas, we define the supersaturation S as a simple deviation from equilibrium concentration:

S = C − Ceq

(4)

The growth rate is calculated from the flux of carbon arriving normal to the crystal surface (eq 5). In this case, we assume that all carbon atoms that reach the surface contribute to the growth without any kinetic limitation.

Vg = (3.6 × 109)



MSiC ( − D∇C + Cu) · n ρSiC

(5)

with n being the unit vector perpendicular to the crystal surface. The first term before the parentheses is a constant necessary to change the carbon flux (in mol/m2·s) used in calculations, to growth rate (in μm/h) directly comparable with the experimental results. It is worth noticing that our model does not contain any adjustable parameter. The results are only dependent on the input data accuracy such as carbon solubility, diffusion coefficient, liquid silicon viscosity, and so on.

CALCULATION DETAILS The simulations are conducted by the finite element method with the COMSOL multiphysics software. Most of the calculation details, regarding induction heating, heat transfers, and fluid dynamics have been presented in a previous paper which included both models and input parameters used.5 Buoyancy, forced convection, Marangoni convection, and electromagnetic convection have been considered in the present study. Solutal convections are not taken into account, as the carbon concentration in liquid silicon is lower than 1% over all the investigated temperature range. To get rid of thermocapillary convection and more generally of any turbulent flow, which appears to be very harmful in the TSSG process of SiC,5 we first designed a totally new crucible geometry with special care paid on the reduction of the temperature gradients. This is practically done by: • increasing the crucible size and insulation thickness with respect to the liquid volume, • increasing the induction coil height so as to spread the heating area and to avoid any tight hot finger, • closing as much as possible the crucible to avoid any strong radiative heat sink. The carbon concentration (molar concentration C) is computed using the classical convection−diffusion equation for the steady state (eq 1), with D being the diffusion coefficient of carbon in liquid silicon (1.7 × 10−8 m2/s given in ref 6b) and u the fluid velocity vector.



RESULTS AND DISCUSSION Hereafter will be displayed only the liquid area, as no additional information can be deduced from the full representation of the crucible. Figure 2 gives an example of the temperature

(1) − D∇2C + u · ∇C = 0 A carbon equilibrium concentration (Ceq) is set as boundary condition for both the crucible−liquid and the crystal−liquid interfaces. It is obtained from the carbon molar fraction xCeq, calculated according to eq 2: xC eq ρ Ceq = Si MSi 1 − xC eq (2)

with ρSi being the liquid silicon density and MSi the molar weight of silicon. The temperature dependence of carbon solubility in liquid Si, given by Durand et al.,8 is accounted for by eq 3. (3)

Figure 2. Temperature distribution in the liquid zone for three different ranges of temperature: (a) 1700 °C, (b) 1800 °C, and (c) 1900 °C. A half crucible only is represented, with the dotted line being the rotation axis.

Special attention is paid to the discretization of the liquid area. Triangular meshing is used for the whole computation domain, except at the crystal surface, where rectangular meshing has proven to yield a better accuracy and stability in

distribution for three different ranges of temperature, namely 1700, 1800, and 1900 °C. The temperature scales are adapted each time for a better visualization. The hot point is the bottom

xCeq = exp(6.249 − 24460/T )

910

dx.doi.org/10.1021/cg201343w | Cryst. Growth Des. 2012, 12, 909−913

Crystal Growth & Design

Article

the crystal, a strong convection loop close to the liquid surface is formed for every condition. This loop is more packed against the liquid surface with the increase in the rotation rate of the crystal. As can be expected, the carbon concentration distribution (shape of the isovalues) is directly correlated to the convection patterns (convection loops). The maximum carbon concentration is thus encountered along the symmetric axis of the crucible, i.e. at the intersection between the two dominating convection loops (Figure 4a). The supersaturation S is

corner of the crucible, and the cold point is the crystal. Note that, for the three cases, the axial temperature difference is about 3−5 °C and the radial one close to 1 °C gives rise to very low temperature gradients, typically about 1 °C/cm. Obviously, the increase of absolute temperature is associated with a slight increase of the axial temperature gradient. For a fixed crucible geometry, it is practically not possible to work with a constant temperature gradient while changing the temperature range. As shown in Figure 2, increasing the temperature from 1700 to 1900 °C does not affect significantly the distribution, and the isotherms keep exactly the same shape. As a first approximation, the temperature gradient will be considered as constant in the forthcoming discussions. The effect of the forced convection is evaluated by representing the fluid flow patterns in the crucible as a function of the crystal rotation rate for 10, 20, and 40 rpm (Figure 3).

Figure 4. Simulation of (a) the distribution of carbon concentration (gray scale) and (b) the distribution of supersaturation (color scale) in the liquid. Blue arrows indicate the normalized carbon flux. Calculations are done for a temperature of 1900 °C and a crystal rotation rate of about 40 rpm. S is calculated from eq 3. A half crucible only is represented, with the dotted line being the rotation axis.

consequently mainly driven by the convective phenomena in the liquid (Figure 4b). Dark blue corresponds to supersaturated area (S > 0), meaning that crystallization can occur. Dark red corresponds to the undersaturated area (S < 0), where dissolution of carbon from the crucible should occur. As expected, the hottest points of the crucible act as the carbon source, and the coldest parts constitute the crystallization areas. The strong lateral convection loop should not be problematic, as the temperature gradient along the liquid surface is close to 0. The supersaturation is thus much higher at the crystal surface than on the lateral side of the crystal. This means that no parasitic crystallization should occur at the liquid surface. This is experimentally confirmed after 8 h or growth. The effect of temperature on the growth rate is presented in Figure 5. Values are taken on the symmetry axis of the crucible, which correspond to the center of the crystal. As can be seen, the growth rate increase with temperature is correctly depicted by the model. The activation energies related to the growth rate are 396 and 226 kJ/mol for the calculated and experimental values, respectively. The activation energy related to the carbon solubility in pure liquid silicon (from eq 3) is 203 kJ/mol. From the very close values of the activation energies of both the measured growth rate and the carbon solubility, we can reasonably speculate that the carbon solubility in liquid silicon is the limiting parameter for the growth rate. A transport limited process should have led to a much lower value. Nasir Khan et al. have also measured the activation energy for both SiC polarities (C-face and Si-face) for two different growth

Figure 3. Representation of the fluid flow patterns in the crucible as a function of the crystal rotation rate, for (a) 10 rpm, (b) 20 rpm, and (c) 40 rpm. The color scale represents the magnitude of the fluid velocity vectors. The red arrows are normalized velocity vectors, added for a better visualization of the different convection loops. Calculations done at 1900 °C. A half crucible only is represented, with the dotted line being the rotation axis.

The mixture of the different convective flows gives rise to complex flow patterns, which can be controlled by increasing the forced convection contribution.5 At a low crystal rotation rate (10 rpm), several convection loops form. As a consequence, the carbon transport associated with the convective flows is not entirely directed to the crystal surface. By increasing the crystal rotation rate to 20 rpm, one loop is developing to the detriment of the others. Further increase of the rotation rate to 40 rpm leads to a single dominant convection loop below the crystal surface. This can be considered as an optimal situation, as all the carbon transport through convection is now directed toward the crystal. Besides 911

dx.doi.org/10.1021/cg201343w | Cryst. Growth Des. 2012, 12, 909−913

Crystal Growth & Design

Article

crystal where a drastic increase of the growth rate is observed. This effect is directly related to an increase of the supersaturation gradient normal to the surface at the edge of the crystal (see Figure 4b). In this area, the calculated growth rate is clearly underestimated compared to the measured one. From optical microscopy observation, we do not observe any difference in the growth mechanism between the edge and the center, which can explain the difference in the experimental growth rates. Although we cannot yet attribute with certainty this deviation to a well identified cause, we note that the right angle at the edge of the crystal can easily give rise to numerical errors when computing the liquid parameters (velocity, concentration, ...) in its close vicinity. And once again, the trends are correctly depicted by the model. For all the investigated conditions, the growth front is perfectly stable and exhibits two-dimensional features (Figure

Figure 5. (a) Calculated and measured growth rates on the symmetry axis of the crucible, for three different temperatures: 1700, 1800, and 1900 °C. The crystal rotation rate is 40 rpm. (b) Same data plotted as Arrhenius laws.

geometries.9 They found an activation energy of 263 kJ/mol for the carbon face, with the horizontal substrate geometry. Our value is very close to theirs, though we cannot conclude that both values are related to the same mechanism. Actually, they obtained geometry dependent values together with the polarity dependent value. The interplay of the process (transport) and the surface kinetics is not clear. Furthermore, as they measured the activation energies in the 1500−1700 °C range, an interface kinetically limited process is more likely to occur than in the 1700−1900 °C range, studied in the present paper. The activation energy extracted from our simulations is a bit overestimated compared to the experimental value, about more than two times, but it keeps the same order of magnitude. The trend is thus correctly depicted by the model. A profile along the crystal radius of both the calculated and measured growth rates normal to the surface is given in Figure

Figure 7. Nomarski differential interference contrast (NDIC) microscopy of as grown surfaces for the samples obtained at (a) 1700 °C, (b) 1800 °C, and (c) 1900 °C after 8 h of growth. The rotation rate of the crystal is 40 rpm.

7). Surface steps are bunched, which is usual during solution growth. Samples grown at 1800 and 1900 °C exhibit a very similar surface morphology. Raman spectra collected on those two layers indicate that only the 6H-SiC polytype is formed. The layer grown at 1700 °C is different and presents straighter step edges (less wavy). In this case, Raman spectroscopy reveals that the 3C polytype has been grown. In addition, this layer observed by NDIC microscopy seems free of double positioning boundaries, which are usually easy to evidence by a simple observation of the growth surface morphology. This unexpected result is thus only related to the growth temperature, as all the other parameters are kept constant. We do not have a clear explanation for this result; it will be the topic of future investigations.

Figure 6. Calculated and measured growth rate profiles along the crystal radius. 0 corresponds to the symmetry axis of the crucible and 5 mm to the edge of the crystal. The experiment is done at 1800 °C with a crystal rotation rate of 40 rpm.

6 for an experiment achieved at 1800 °C. The growth rate is almost homogeneous along the radius, except at the edge of the 912

dx.doi.org/10.1021/cg201343w | Cryst. Growth Des. 2012, 12, 909−913

Crystal Growth & Design



Article

(7) Alexander.; Seki, K.; Kozawa, S.; Yamamoto, Y.; Ujihara, T.; Takeda, Y. Polytype stability of 4H-SiC seed crystal on solution growth. Mater. Sci. Forum 2011, 679−680, 24−27. (8) Durand, F.; Duby, J. C. Carbon solubility in solid and liquid silicona review with reference to eutectic equilibrium. J. Phase Equilib. 1999, 20 (1), 61−63. (9) Nasir Khan, M.; Nishizawa, S.-i.; Bahng, W.; Arai, K. Liquid-phase epitaxy on 6H-SiC Acheson seed crystals in closed vessel. J. Cryst. Growth 2000, 220, 75−81.

CONCLUSION As part of the investigation of top seeded solution growth (TSSG) as an alternative to the standard industrial PVT process for the growth of high quality single crystalline SiC ingots, a numerical model of the TSSG process has been developed. The main aim was to give quantitative outcomes in addition to qualitative information. We have first developed a specific crucible design with very low temperature gradients (1 °C/cm). Then, the carbon concentration and the supersaturation distributions have been computed to evidence the dissolution and growth areas. As expected, we found a direct correlation between the concentration fields and the convection patterns. This has emphasized the primary role of the convection loops on the carbon transport from the dissolution to the crystallization zones. Growth rates have been calculated for three different temperature ranges (1700, 1800, and 1900 °C) in pure silicon. In addition, the growth rate profile along the crystal radius has been computed. We obtained a good agreement between the computed and the experimental trends. Finally, we have evidenced that, under the conditions investigated, the limiting parameter is the carbon solubility in liquid silicon and not the carbon transport. The numerical tool is now ready to assist the future development of the process, with the next step being the increase of the growth rate. For that, higher growth rates will necessarily require the use of solvents able to increase the carbon solubility, such as Ti−Si or Cr−Si.



AUTHOR INFORMATION

Corresponding Author

*Telephone: +33 456 529 329. Fax: +33 456 529 301. E-mail: [email protected].



ACKNOWLEDGMENTS The authors warmly acknowledge Dr. F. Mercier from AIST Tsukuba (Japan) for fruitful discussions. This work has been supported by the French ANR project MINTEX Contract Nb. ANR-09-BLAN-0189-01 and by the French-Japanese SAKURA project (Contract Nb. 23544YE).



REFERENCES

(1) Hofmann, D. H.; Muller, M. H. Prospects of the use of liquid phase techniques for the growth of bulk silicon carbide crystals. Mater. Sci. Eng., B: Solid State Mater. Adv. Technol. 1999, 61−2, 29−39. (2) Chaussende, D.; Wellmann, P. J.; Pons, M. Status of SiC bulk growth processes. J. Phys. D: Appl. Phys. 2007, 40 (20), 6150−6158. (3) Kamei, K.; Kusunoki, K.; Yashiro, N.; Okada, N.; Tanaka, T.; Yauchi, A. Solution growth of single crystalline 6H, 4H-SiC using SiTi-C melt. J. Cryst. Growth 2009, 311 (3), 855−858. (4) Danno, K.; Saitoh, H.; Seki, A.; Daikoku, H.; Fujiwara, Y.; Ishii, T.; Sakamoto, H.; Kawai, Y. High-Speed Growth of High-Quality 4HSiC Bulk by Solution Growth using Si-Cr Based Melt. Mater. Sci. Forum 2009, 645−648 (645−648), 13−16. (5) Mercier, F.; Dedulle, J. M.; Chaussende, D.; Pons, M. Coupled heat transfer and fluid dynamics modeling of high-temperature SiC solution growth. J. Cryst. Growth 2010, 312 (2), 155−163. (6) (a) Mercier, F.; Nishizawa, S. i. Solution growth of SiC from silicon melts: Influence of the alternative magnetic field on fluid dynamics. J. Cryst. Growth 2011, 318 (1), 385−388. (b) Mercier, F.; Nishizawa, S. Numerical Investigation of the Growth Rate Enhancement of SiC Crystal Growth from Silicon Melts. Jpn. J. Appl. Phys. 2011, 50 (3). 913

dx.doi.org/10.1021/cg201343w | Cryst. Growth Des. 2012, 12, 909−913