Effects of ampoule rotation on flows and dopant segregation in
vertical Bridgman crystal growth
C.W. Lan*
Chemical Engineering Department, National Taiwan University, Taipei 10617, Taiwan, ROC
Received 30 July 1998; accepted 11 September 1998
Abstract
Computer simulation is conducted to study the effects of ampoule rotation on the flows and the dopant segregation in vertical Bridgman (VB) crystal growth. A benchmark problem for the growth of gallium-doped germanium in a graphite ampoule is used as an example. Both pseudo-steady state and fully transient calculations are conducted. It is found that the convection and dopant segregation can be affected significantly by ampoule rotation, and the effect is similar to that by an axial magnetic field. Ampoule rotation decreases the intensity of convection and stretches the flow cell axially. As a result, both radial and axial segregations are affected. When the convection is weak, the flow can be suppressed almost completely by moderate ampoule rotation and the dopant segregation becomes diffusion-controlled. For stronger convection, the elongated flow cell by ampoule rotation may bring dopant mixing into the bulk melt reducing axial segregation at the early stage of the growth. However, if the cellular flow cannot be suppressed completely, ampoule rotation may induce larger radial segregation due to poor mixing. 1999 Elsevier Science B.V. All rights reserved.
PACS: 44.25.#f; 47.27.Te; 81.10.Fq; 02.60.c6; 02.70.Fj
Keywords: Ampoule rotation; Dopant segregation; Bridgman; Convection; Mixing
1. Introduction
The control of melt flow is important in crystal growth. Besides eliminating unstable flow that causes growth striations, it also affects dopant seg-regation. A well-known approach to control the melt flow is the use of magnetic fields (e.g. Refs. [1,2]). For vertical Bridgman (VB) crystal growth,
* Fax: #886 2 2363 3917; e-mail: lan @ruby.che.ntu.edu.tw.
the effects of magnetic fields, especially the axial ones, to suppress the convection and thus to con-trol dopant segregation have been studied exten-sively (e.g. Refs. [3—5]). However, the control of melt flow by ampoule rotation has not been investi-gated. Ampoule rotation is often used in VB crystal growth to achieve an axisymmetric thermal field. The effect of ampoule rotation on the flow and dopant segregation is often believed to be small in practice, especially for slow rotation. However, from fluid mechanics point of view, the importance of the ampoule rotation (say several tens of rpm)
0022-0248/99/$ — see front matter 1999 Elsevier Science B.V. All rights reserved. PII: S 0 0 2 2 - 0 2 4 8 ( 9 8 ) 0 0 8 6 8 - 9
depends on the intensity of buoyancy flow caused by radial thermal gradients. For the VB configura-tion, because of the stabilized thermal fields, as compared with that in other crystal growth pro-cesses, the thermal convection is usually very weak. Therefore, a typical ampoule rotation speed, say 20 rpm, could have a profound effect in the flow and the dopant transport. In fact, the contribution of rotation to flow stability is long recognized [6]. For example, the critical Rayleigh number for the Reyleigh—Benard problem in a cylinder increases with the rotation speed (Taylor number). Accord-ing to the Taylor—Proudman theorem, the inviscid fluid (Prandtl number Pr"0) in rotation is always thermally stable for all adverse temperature gradi-ents. Therefore, the suppression of the buoyancy convection by ampoule rotation is feasible. Al-though the flow and dopant segregation, both radial and axial, have been studied extensively (e.g. Refs. [3—5,7—11]) with and without a magnetic field, the effect of ampoule rotation for VB crystal growth, to the best of our knowledge, has not been investigated.
In this study, the effect of steady ampoule rota-tion on VB crystal growth is investigated numer-ically. The effect of ampoule rotation on the buoyancy flow and dopant segregation, both radial and axial, at different Rayleigh numbers is dis-cussed. For comparison purposes, we consider the gallium-doped germanium (GaGe) system used by Adornato and Brown [7]. Extensive benchmark comparison for this system was also conducted by Lan and Chen [9] and Liang and Lan [10] using different numerical approaches including the finite volume method (FVM) and the finite element method (FEM). Therefore, its numerical analysis provides a good candidate for the theoretical study here. In the next section, the mathematical model and its numerical solution are briefly described. Section 3 is devoted to the results and discussion, followed by conclusions in Section 4.
2. Model and numerical solution
The schematic of the VB crystal growth used in this study is depicted in Fig. 1. This system is the same as the one used by Lan and Chen [9];
Fig. 1. Schematic sketch of vertical Bridgman (VB) crystal growth.
however, the growth dynamics is also considered here. The furnace is described by an effective heat-ing profile ¹ (z, t), which is assumed linear in this study. To start crystal growth from a stationary state, this profile is moved upward at speed º. For pseudo-steady-state calculations, the thermal pro-file is kept stationary and the ampoule is moved downward at a speed º . The system is assumed to be axisymmetric, and the initial dopant distribution in the melt is assumed uniform at C. The flow and temperature fields, as well as the growth front (the melt/crystal interface, h(r, t)), are represented in a cylindrical coordinate system (r, z).
The melt is assumed incompressible and Newto-nian, while the flow is laminar. The Boussinesq approximation is also adopted. If the stream func-tiont and vorticity u are defined in terms of radial (u) and axial (v) velocities as
u"! 1 o r jt jz, v" 1 o r jt jr, (1)
and
u"jujz!jv
jr, (2)
the unsteady-state governing equations for stream function, vorticity, azimuthal velocity (vF), temper-ature (¹), and dopant concentration (C) in the conservative-law form can be written as follows:
Equation of motion: !o ju jt#jrj
ur jt jz!jzjur jt jr #j jr 1 r j jr(k ru)# j jz 1 r j jz(k ru) !o b2gj¹ jr#o b1g jC jr#o j jz vF r"0. (3) Circulation equation: !o jvF jt# 1 r j jrrvF jt jz! j jzrvF jt jr #j jrk r j jr(rvF)# j jzk r j jz(rvF)"0. (4) Stream equation: j jz 1 o r jt jz# j jr 1 o r jt jr#u"0. (5) Energy equation: !oGCNrj¹ jt! j jr(oGCNGru¹)! j jz(oGCNGrv¹) #j jzrkG j¹ jz# j jrrkG j¹ jr"0, i"(m,c,a), (6) Dopant equation: !rjC jt! j jr C o jt jz! j jz C o jt jr #j jzrD jC jz# j jrrD jC jr"0, (7) wherek is the melt viscosity. Also, oG is the density,CNG the specific heat, and kG the thermal
conductiv-ity of phase i; i represents the melt (m), the
crystal (c), or the ampoule (a). Also, g is the gravi-tational acceleration,b2 and b1 are the thermal and solutal expansion coefficients, and D is the dopant diffusivity in the melt. The associated thermal (Ra2) and solutal (Ra1) Rayleigh numbers are defined as follows:
Ra2,b2¸g¹ /(l a ), Ra1,b1¸gC/(l a ),
where ¹ is the melting temperature of the pure material, l "k /o , and a "k /(o C ). The characteristic length ¸ is the domain length as shown in Fig. 1. All the physical properties used in Ra2 and Ra1 are obtained at ¹ .
To solve the previous equations, boundary con-ditions need to be specified. Most boundary condi-tions can be found elsewhere [9,12]. In brief, the no-slip boundary condition for velocities is used at solid boundaries. For example, at the melt/solid interface,
vF"rX, (8)
whereX is the ampoule rotation speed. The upper melt interface is assumed stress free. For pseudo-state calculations, the solute boundary conditions can also be found elsewhere [9]. For time-depen-dent calculations, boundary conditions for the dopant at the upper and melt/crystal interfaces are obtained through the interface dopant balance: !(n ) DC)"X*"0 (9) and
!(n ) DC)"!(1!K)Cjh
jt(n ) eX)
"0, (10) where K is the segregation coefficient obtained from the phase diagram; K,C/C at the growth interface, where C is the dopant concentration in the crystal. In the above equations, the dopant diffusion in the solid is neglected. The effectiveness of mixing can be measured by the effective segrega-tion coefficient K as defined before [12]:
where CM is the laterally averaged dopant concen-tration in the growth interface in the crystal side, while1C2 the volumetric average of dopant con-centration in the melt. The zero flux condition is used at the centerline and the melt/ampoule inter-face. Furthermore, in this paper, the solutal effect and the non dilute effect on the growth temperature and melt flow are not considered; Ra1"0. These effects can be included easily if necessary.
Finally, heat transfer from the system to the ambient is governed by the simple Newton’s cool-ing law:
!(nk )¹)"h[¹!¹ (z, t)], (12) where h is the heat transfer coefficient and ¹ (z, t) the effective ambient temperature, which is a linear distribution. The top and bottom temperatures of the system are equal to the effective ambient temperature. Again, more realistic boundary condi-tions can be used if necessary. However, for comparison, we have chosen the same thermal boundary conditions used before by Lan and Chen [9].
The above governing equations and their asso-ciated boundary conditions can only be solved nu-merically. We have developed efficient FVM schemes using Newton’s method and the DASPK solver [9,12,13] for solving these equations. Excel-lent global dopant conservation can be easily achieved even for a coarse mesh. Detailed descrip-tion of the numerical method can be found else-where [12].
3. Results and discussion
The benchmark system for gallium-doped ger-manium (GaGe) is used here to illustrate the effect of ampoule rotation on the flow and dopant segre-gation. The physical properties and some input parameters are listed in Table 1 for reference. The mesh used here is 31;21 in the crystal, 31;101 in the melt, and 11;121 in the ampoule. Such a mesh is found satisfactory, and further mesh refinement does not improve accuracy much. We have per-formed extensive mesh refinements as before [9], and the detailed benchmark comparison of FVM
Table 1
Physical properties and some input parameters [9]
GaGe o"5.5 g cm\ o "5.5 g cm\ ¹ "937°C *H"460 J g\ h"46.571 W cm\ °C\ k"0.17 W cm\ °C\ k "0.39 W cm\ °C\ Cp"Cp "0.39 J g\ °C\ b2"5;10\ K\ b1"0 (mol% Ga)\ D"2.1;10\ cm/s K"0.087 Graphite o "1.8 g cm\ k "3.26 W cm\ °C\ Cp "1.814 J g\ °C\ Other input parameters
¸"7 cm
R"0.5 cm R "0.7 cm
¹&"1112°C ¹!"762°C
º "!4;10\ cm/s (for pseudo-steady-state calculations) º"#4;10\ cm/s (for fully transient calculations)
and FEM calculations can be found elsewhere. An additional comparison for 10 rpm ampoule rota-tion at Ra2"2.489;10 (equivalent to the normal gravity condition) and º "º"0 is conducted here. Fig. 2 shows the results obtained by the pres-ent (solid-line) and the FEM (dashed-line) ap-proaches; the FEM approach is based on the primitive variable (º»P) formulation [9]. As shown, the calculated results for both flow pattern (t on the left hand side) and rvF, as well as the growth front, are in excellent agreement. Now, we are ready to present our results. In the following discussion, we will start from pseudo-steady-state (PSS) calculations first. The PSS results can be obtained quickly and the physical significance of ampoule rotation on the flow and dopant mixing can be grasped easily. However, detailed global behavior and axial dopant segregation still require fully transient calculations, which will be presented later.
Fig. 2. Comparison of the calculated results by FVM (t/u) and FEM (º»P) for 10 rpm ampoule rotation and Ra2"2.489;10.
For the cases without ampoule rotation, the PSS results at different convection levels are shown in Fig. 3. In each plot, the left-hand side represents the stream function, while the right-hand side the iso-concentration lines. As shown, at Ra2"0, the straight streamlines are due to ampoule translation. Also, the concave interface is mainly caused by the difference of the thermal conductivity in the melt (k ) and the crystal (k); k is twice smaller than k . As the buoyancy force is considered, at Ra2"2.489;10 (equivalent to 0.01 g), a cellular flow upon the growth front is induced due to radial thermal gradients. Although the flow intensity is low, the local dopant field near the growing inter-face is affected significantly, which can be seen from the distortion of the iso-concentration lines. As a result, the radial dopant segregation increases (C /C"13.876); the mean value of C/C is 11.494, i.e., 1/K, and K"0.2393. For the con-venience of representation, hereafter the nor-malized dopant concentration C/C or C/C is represented by C or C. Although the convection in Fig. 3b is quite weak, in the early study by Kim et al. [3] and Adornato and Brown [7] using a
Fig. 3. Calculated flow patterns (LHS) and dopant fields (RHS) for 0 rpm ampoule rotation: (a) Ra2"0; (b) Ra2"2.489;10; (c) Ra2"2.489;10; *t"t /10 for negative t and *t"t /10 for positive t, and *C"(C !C )/10.
FEM scheme, the convection level considered was in this order of magnitude. The poor global conservation of their FEM scheme for dopant cal-culations made their calcal-culations for higher con-vection levels difficult unless an extremely fine mesh was used [9]. Under normal gravity condi-tions in Fig. 3c, Ra2"2.489;10, because of the much stronger convection, the dopant near the growth front is well mixed, and the radial dopant segregation decreases (C "13.047); K"0.2204. Away from the growth interface, the diffusion mechanism is still dominant. Apparently, the larger radial segregation in Fig. 3b is due to poor dopant mixing. In the case of no mixing in Fig. 3a, the radial dopant segregation is due to the deformed growth interface.
With 20 rpm ampoule rotation, the flow patterns are significantly changed, as shown from Fig. 4a—Fig. 4c. At Ra2"0, the melt is rigidly rotating, and the melt flow due to ampoule transla-tion is not affected at all. The most significant change is found in Fig. 4b at Ra2"2.489;10. As shown the flow cell disappears due to ampoule rotation. As a result, the dopant transport upon the growth front becomes diffusion controlled. In other words, the suppression of the thermal convection
Fig. 4. Calculated flow and dopant fields for 20 and 100 rpm ampoule rotation, respectively: (a) and (d) Ra2"0; (b) and (e) Ra2"2.489;10; (c) and (f) Ra2"2.489;10.
by 20 rpm ampoule rotation is effective when the convection level is low. The radial dopant segrega-tion is also greatly reduced (C is reduced).
At a much higher convection level shown in Fig. 4c (Ra2"2.489;10), although the flow intensity is greatly reduced by rotation (from 5.567;10\ g/s in Fig. 3c to 8.954;10\ g/s in Fig. 4c), the effect of convection on the dopant field is still significant. More importantly, the rotation pushes the flow cell toward the ampoule wall and stretches it axially. As shown, the main flow cell fills almost the whole melt domain. As a result, the elongated cell brings dopant mixing deeply into the
bulk melt. Therefore, even though the radial dopant segregation increases due to poorer local mixing near the growth front, as will be shown shortly, the axial segregation is reduced due to the deeper dopant mixing into the bulk melt; K" 0.1153.
As the rotation speed is further increased to 100 rpm, the results at Ra2"0 are still not affected (Fig. 4d). The thermal convection at lower Ra2, as shown in Fig. 4e, is suppressed almost completely. Therefore, the flow and dopant fields are nearly the same as that in Fig. 4d; the dopant transport now is in the diffusion-controlled limit. The flow cell at higher Ra2 may not be suppressed completely un-der 100 rpm ampoule rotation. As shown in Fig. 4f, the buoyancy cell is further stretched and its inten-sity is greatly reduced. Interestingly, because of the poor dopant mixing by the weakening flow, the radial dopant segregation increases significantly (C "15.844); K"0.217.
The effect of ampoule rotation can be further illustrated by the axial dopant distribution in the melt. Taking the dopant distribution at the centerline for Ra2"0 and Ra2"2.489;10 as an example, the results from Figs. 3 and 4 are present-ed in Fig. 5. As shown, 20 rpm ampoule rotation has a profound effect in the dopant distribution. The dopant penetrates deeply into the bulk melt for Ra"2.489;10 due to 20 rpm ampoule rotation. As a result, the axial segregation may be reduced; we will use dynamic calculations to support this shortly. Furthermore, increasing ampoule rotation speed pushes the profile toward the diffusion-con-trolled limit (Ra2"0).
The PSS calculations presented so far have illus-trated the physical significance of ampoule rota-tion. However, during crystal growth the melt height is changing. Therefore, to understand the detailed global segregation behavior, transient cal-culations are necessary, even though they are quite time consuming. The effect of ampoule rotation for crystal growth over a long period of time for Ra2"2.489;10 is illustrated in Fig. 6; the heater moving speed º"4;10\ cm/s. Fig. 6a and Fig. 6d are steady-state solutions at º "º"0 for 0 and 20 rpm, respectively. To start the growth, the heating profile is suddenly moved upward at
Fig. 5. Axial dopant distribution at the centerline from Figs. 3 and 4.
shown in Fig. 6b and Fig. 6e, respectively. Similar to the PSS results, the dopant mixing for 20 rpm ampoule rotation penetrates deeply into the bulk melt, while the dopant mixing away from the growth front for 0 rpm is poor. Nevertheless, at the end of growth, t"7500 s, because the melt zone is short, the dopant is well mixed in the whole melt for 0 rpm. On the contrary, the dopant mixing is less for 20 rpm due to the weaker convection.
The axial dopant segregations in the grown crys-tal for Fig. 6a—Fig. 6f are further illustrated in Fig. 7 by taking the laterally averaged dopant con-centration in the grown crystal (C/C) versus the solidified fraction. For comparison purposes, we also include the analytical result by Smith et al. [14] for pure diffusion. As shown, good agreement is found for Ra2"0 except at the end of growth. The discrepancy is due to the breakdown of the Smith model near the end boundary, which is based on an infinite domain assumption. The segregation for complete mixing (the Scheil equation [15]) is included as well. As shown, as compared with the case of no rotation, the axial dopant segregation for 20 rpm ampoule rotation, as expected from Fig. 6d, is closer to the complete-mixing result at the first half of the growth. Near the end of the growth, due
Fig. 6. Calculated results for 0 and 20 rpm ampoule rotation, respectively: (a) and (d) t"0; (b) and (e) t"4000 s; (c) and (f)
t"7500 s. Ra2"2.489;10 and º"5;10\ cm/s.
to the stronger melt convection, the result of no rotation becomes closer to the complete-mixing one; also see Fig. 6c.
The effect of 100 rpm ampoule rotation on the axial dopant segregation is also illustrated in Fig. 8. Obviously, increasing ampoule rotation speed moves the segregation curves toward the diffusion-controlled limit, i.e., Ra2"0. In other words, faster ampoule rotation decreases flow intensity and thus increases axial dopant segrega-tion. As a result, the effectiveness of dopant mixing decreases. Fig. 9 shows the radial dopant segrega-tion, measured by C"P0!C"P as a function of
Fig. 7. Effect of 20 rpm ampoule rotation on axial dopant segre-gation for different Rayleigh numbers.
Fig. 8. Effect of 100 rpm ampoule rotation on axial dopant segregation for different Rayleigh numbers.
radial dopant segregation increases with the in-creasing ampoule rotation speed. Again, this is due to the poor dopant mixing caused by the weaken-ing flow. Such an effect will be more significant for larger crystal diameter. Nevertheless, we also ex-pect that the radial dopant segregation can be re-duced if the diffusion-controlled limit can be reached as illustrated in Fig. 4e. However, for the
Fig. 9. Effect of ampoule rotation on radial dopant segregation as a function of growth distance.
present case (Ra2"2.489;10), this requires faster ampoule rotation to suppress the thermal convec-tion sufficiently.
4. Conclusions
The effect of ampoule rotation on the flows and dopant segregation in vertical Bridgman crystal growth is investigated numerically. From the cal-culated results, it is clear that a moderate ampoule rotation speed can significantly affect the flows and further the dopant mixing. Similar to axially mag-netic damping, the Coriolis force due to ampoule rotation can suppress thermal convection. The flow inhibition by ampoule rotation can be regarded as an extension of the Taylor—Proudman theorem for inviscid melt. According to the theory, an addi-tional advantage is that the three-dimensional motion can be suppressed as well [6]. Nevertheless, due to the poorer local mixing in front of the growth interface, in the present study, ampoule rotation may induce larger radial dopant segrega-tion, unless the thermal convection is suppressed sufficiently. In summary, ampoule rotation is be-lieved to have profound effects on vertical Bridg-man crystal growth. However, how to better control the melt convection and dopant transport
by ampoule rotation still requires further study, especially, through crystal growth experiments. References
[1] H.A. Chedzey, D.T.J. Hurle, Nature 210 (1966) 933. [2] H.P. Utech, M.C. Flemming, J. Appl. Phys. 37 (1966) 2021. [3] D.H. Kim, P.M. Adornato, R.A. Brown, J. Crystal Growth
89 (1988) 339.
[4] S. Motakef, J. Crystal Growth 104 (1990) 833.
[5] M. Yao, A. Chait, A.L. Fripp, W.J. Debnam, J. Crystal Growth 173 (1997) 467.
[6] S. Chandrasekhar, Hydrodynamic and Hydromagnetic Stability, Dover, New York, 1961, p. 76.
[7] P.M. Adornato, R.A. Brown, J. Crystal Growth 80 (1987) 155.
[8] R.A. Brown, D.H. Kim, J. Crystal Growth 109 (1991) 50.
[9] C.W. Lan, F.C. Chen, Comput. Meth. Appl. Mech. Eng. 131 (1996) 191.
[10] M.C. Liang, C.W. Lan, J. Crystal Growth 167 (1996) 320. [11] S. Kuppurao, S. Brandon, J.J. Derby, J. Crystal Growth
155 (1995) 103.
[12] C.W. Lan, M.C. Liang, J. Crystal Growth 186 (1998) 203.
[13] P.N. Brown, A.C. Hindmarsh, L.R. Petzold, SIAM J. Sci. Comput. 15 (1994) 1467.
[14] V.G. Smith, W.A. Tiller, J.W. Rutter, Can. J. Phys. 33 (1955) 723.
