- Email: [email protected]
The effect of an axial magnetic field on periodic growth rate fluctuations in Czochralski crystal growth
466
Journal ofCrystal Growth 87 (1988) 466—480 North-Holland, Amsterdam
THE EFFECT OF AN AXIAL MAGNETIC FIELD ON PERIODIC GROWTH RATE FLUCFUATIONS IN CZOCHRALSKI CRYSTAL GROWTH R.A. CARTWRIGHT and A.A. WHEELER School of Mathematics, University of Bristol, University Walk, Bristol BS8 ITW, UK Received 15 June 1987; manuscript received in final form 25 September 1987
In this paper we develop a mathematical model to investigate the effect of growth rate fluctuations in magnetic Czochralski crystal growth. The pull rate is assumed constant and the growth rate fluctuates with frequency a rad/s, the crystal rotation rate. We explore the effect of the magnetic field and the crucible rotation on fluid flow and mass transfer in the vicinity of the crystal/melt interface. Numerical results are obtained for 0 N 32 and —0.1 < s 1.0, where N is the magnetic interaction parameter and s is the ratio of crucible to crystal rotation. Asymptotic results are obtained for the velocity field for large N and agreement is good for 2 —~oo where Sc is the N> 4. The analysis is extended to the concentration givethe asymptotic for NWe —* ~o, Sc —, c~,N/Sc Schmidt number. Again the asymptotic results agree field wellto with numericalresults solutions. conclude that the phase of the resulting concentration fluctuations appears to be independent of N and s. The spatially averaged effective distribution coefficient, k~ffand the amplitude of the concentration oscillations, A~,are linked. As keff -~ 1, A~increases and so improved macrosegregation is accompanied by a deterioration in microsegregation.
1. Introduction It is well known that crystals grown by the Czochralski technique may exhibit striations in which a spatially periodic variation of dopant concentration is observed in the direction of the axis of the crystal. This phenomenon has been attributed to periodic temporal fluctuations in the crystal growth rate, caused by thermal asymmetries in the melt (in which case they will have the same period as that of the rotating crystal) or unsteady thermal convection in the bulk of the melt, (Choe and Strudwick [1]). This paper investigates the effect of an axial magnetic field and crucible rotation on the flow near the crystal/melt interface and the concentration of impurity there. In particular we consider whether either of these two control mechanisms will help to damp out the fluctuations in concentration and so enable crystals with a more uniform distribution of impurity to be grown. In recent years semi-conductor chips have been progressively reduced in size, whilst at the same time device performance requirements have become
more stringent. These two facts result in the necessity for better crystal quality; a uniform distribution of impurities at prescribed levels. A magnetic field has been suggested (first by Chedzey and Hurle [2] in 1966) as a means of achieving this aim since it damps the flow and hence may reduce fluctuations in the system. Previous work has focussed on two areas of the problem. Wilson [3,4] and Wheeler [5] have considered the case of a periodic growth rate in the absence of a magnetic field and without crucible rotation. Cartwright, Hurle, Series and Szekely [6] studied the effect of magnetic fields and crucible rotation when the crystal grows at a constant rate. We consider here a more general model which includes a periodic growth rate, an axial magnetic field and crucible rotation. The work of the above investigators are special cases of our model. All of these workers have confined their attention to the boundary layer adjacent to the crystal. Other important work in this area is by Hjellming and Walker [7—9]who have conducted an analysis of the magnetohydrodynamic flow in the entire melt. Their model assumes that the imposed uniform
0022-0248/88/$03.50 © Elsevier Science Publishers B.V. (North-Holland Physics Publishing Division)
R.A. Cartwrighf, A.A. Wheeler / Axial magnetic field and periodic growth ratefluctuations
magnetic field is sufficiently large that inertia effects may be neglected. This has the mathematical advantage that the resulting governing equations are linear. The major insight of their work is to show that in the presence of such a high magnetic field the electrical conductivity of the crystal, although small, has a profound effect on the flow in the melt. It should be noted that our model does not incorporate the role played by natural convection or surface tension. The effects due to the memscus may be neglected because in the vicinity of the crystal/melt interface the flow is dominated by forced convection. A detailed argument based on the dimensionless groups of the system is given in Cartwright, El-Kaddah and Szekely [10]. With these assumptions the flow decouples from the temperature field and so we consider an isothermal melt. The resulting problem consists of two distinct parts: (i) the flow, (ii) the solute distribution. The first stage is to solve for the fluctuating velocity field and then in the second this solution is substituted into the mass transfer equations in order to determine the concentration field. In section 2 we derive our model and describe our solution method for each part of the problem and in section 3 we present our results. We find that the crucible rotation and the magnetic field strength have a profound effect on the effective distribution coefficient. In general terms, as the magnetic field strength increases the average effective distribution coefficient, kC~1, approaches unity. However, the amplitude of the concentration oscillations also increases and so although a magnetic field improves macrosegregation, microsegregation may not be beneficially affected. The effect of crucible rotation is more complicated but a similar correlation between increased keff and increased amplitude of concentration oscillations is observed. In practice, of course, the magnetic field also has an effect on the thermal asymmetries in the melt which cause the growth rate fluctuations. Therefore, in order to argue conclusively about the effect of magnetic fields on crystal quality it is necessary to consider the complete system. This paper demonstrates the effect of the field on a prescribed growth rate fluctuation.
467
2. The model
2.1. Introduction The result of an axial magnetic field will be to modify the flow in the melt. Here we are concerned with the consequent effect on the grown-in concentration in the crystal, which is determined by the concentration at the crystal/melt interface. Hence we shall consider the flow and solute in the boundary layer adjacent to the growing crystal. Our model consists of an infinite rotating disc representing the crystal/melt interface, which rotates with angular velocity w rad/s; the crystal rotation rate. We justify our model of the crystal as an infinite rotating disc by noting that a recent study by Brady and Durlofsky [11] of the flow adjacent to a finite rotating disc has revealed that the effect of the edges of a finite disk upon the flow in the momentum boundary layer, are, for moderate or large Reynolds number, confined to the outer half of the disc. The flow in the central half is qualitively described by the flow due to an infinite rotating disc over a wide range of conditions. It is the region close to the centre of the crystal to which our model is applicable. We locate our coordinate system coincident with the initial position of the crystal/melt interface which is given by the plane z = 0. The configuration of the system is shown in fig. 1.
e
I
‘i~y w
~
DISC ~
ICwt~)
~ w)
F u id
1.1 SW Fig. 1. Model configuration.
R.A. Carzwright, A.A. Wheeler / Axial magnetic field and periodic growth ratefluctuations
468
To model the rotation of the crucible the fluid is assumed to have a prescribed solid body rotation rate of ~ rad/s far from the crystal, that is, as z oo. The pulling rate of f rn/s is modelled by imposing “suction” in the axial direction at the rotating disc and the fluctuating growth rate, f * (t *) by an oscillatory motion of the disc about z = 0 in the axial direction. We assume that the density of the solid and liquid phases are the same. Hence, in this frame of reference, the fluctuating growth rate does not cause a fluctuating suction velocity at the interface (see Wilson [4]). We use a cylindrical polar coordinate system (r, 9, z) as shown in fig. 1, and we assume that the system is axisymmetric. A magnetic field is applied in the axial direction and because the fluid is electrically conducting the applied magnetic field and velocity field interact to induce electromagnetic effects. In general these effects are magnetic field line distortion due to the velocity field and an inhibition of motion across field lines. Although the applied magnetic field is steady, the velocity field is time dependent and so the induced electromagnetic field will in general also be
parameter defined in the next section and t~ is the radius of the crystal. We now consider the first part of the model; the flow.
—*
time dependent. However, we note that for a silicon melt the magnetic Prandtl number, Prm (= yap.0) where v is kinematic viscosity, ~ is electrical conductivity and p.~is magnetic permeability is 6 x iO~. For similar situations relevant to crystal growth Prm is also small. Thus we shall consider our model in the realistic limit Prm 0. In appendix A we show that if O(Prm) terms are neglected, the magnetic field is steady and is not distorted by the velocity field. Hjellming and Walker [7] showed that the electrical conductivity of the crystal, despite being small, can for large magnetic field be important. We shall ignore this effect in our model and assume the crystal to be a perfect electrical insulator. As noted by Cartwnght, Hurle, Series and Szekely [6] in a similar situation, we must therefore assume that —‘
2<1, N1~~’2o~(co/v)I~’
where v is the kinematic viscosity of the melt, is the ratio of the electrical conductivity of the crystal to that of the melt, N is the magnetic interaction
2.2. The flow We model the liquid phase by a viscous, incompressible, electrically conducting fluid. Thus the governing equations are the time dependent Navier—Stokes equations with an electromagnetic body force term (the Lorentz force) and the continuity equation. We seek an axisymmetric solution and assume that the axial velocity and axial magnetic field are independent of r, in which case the following similarity transformation is appropriate (von Karman [12]):
( ±~) 1/2
~
z
=
,
~(
=
~
_~-
=
cot ~,
~(~
=
wr ~
= 1/2’
‘
where u~,v9, v2 are the dimensional velocity cornponents in the radial, azimuthal and axial directions, v is the kinematic viscosity, t is dimensional time and overbar denotes a dimensionless quantity. The equations of motion then reduce to the following set of partial differential equations: ~ = ~ ‘ii~ —2(G~ 52’l NH~ ~ 2 ~ (la) = ~+—~—N(~—s), (ib) *
—
—
/
with boundary conditions —
-
-
H(~, t)
—
=
-
-
~
~
at —
~
-
-
t)
— —*
0,
-
—
t)
=
=
I(
0,
~)
-
-
G(~,t)
=
1, (2a)
-
G(~,t) —‘s, -
as~—’oo, (2b) where 1(1) is the position of the interface at time t, a = f/( i’w )1~~2is the dimensionless pulling rate, N = B~a/pwis the magnetic interaction parameter, -
R.A. Cartwright, A.A. Wheeler / Axial magneticfield and periodic growth ratefluctuations
B0 is the applied magnetic field strength, p is the density of melt, and terms O(Prm) have been neglected. It should be noted that the growth rate, a(~),is given by a(t)=a+I’(t). In all of the calculations presented in section 3, a sinusoidal growth rate was chosen but the model allows for any periodic function of time. Although F does not appear explicitly in the equations it is easily obtained using the continuity equation which gives
~
t)=
~).
~
tails of which are given in appendix B. The solution obtained is H(y, —
2+11 ($1x 2x+/13) 5), e~ (6a) +($4x+/15) e_2x1 + O(X G(y, t) s + (1 s) e_x ~ (1 s) I’(t)x e_x + O(X2),
0+
—
(6b)
(3)
al. [6,10]. As shown by Wilson [3], it is convenient to define a new spatial coordinate, y, which moves with the interface, but to consider velocities with respect to the fixed coordinate system. Thus we define 1(i),
G(y,
t)
=
t
=
G(~,t),
H(y,
H(~, t).
t)
Hence eqs. (1) and (2) become =
~ —
+ H~~(I’(t) H) + ~H2 —
2(G2
G~= ~
+
—
~2)
—
NH~,
(4a)
GH~+ G~(I’(t) H) —
—
N(G
—
s), (4b)
with H(y,
t)
=
—a,
—
—
as A
~—
t) 3(1—s)[2(1+2s+3sx)e~ —a+~A _(1_s)e_2x_(1+5s)]
_x4[s
The steady version of eqs. (1) and (2) is the same set of equations studied in detail by Cartwright et
=
469
H~(y, t)
0,
G(y,
t)
at y = 0, H~(y,t) 0, as y—*cc
G(y,
t)
=
—~
=
1, (5a)
~
~,
—~
cc,
where
a + (a~+ 4N )1/2 =
,
2
x=Ay,
and /3~flu..., $~ are coefficients which are functions of time given in appendix B. This asymptotic expansion shows clearly the structure of the solution in which solid body rotation (H( y, t) = a, G(y, t) = s) is combined with a boundary layer of thickness O(N”~’2). The second method we have used is a numerical integration of the equations. To do this we employed the method of lines and Gear’s method by using NAG subroutine DO3PBF. This routine —
solves the equations on a finite domain in y and so the boundary conditions at infinity were imposed at yf = 20. The values of G and H~at y = 7 were within 0.01% of the value required of them as y cc. Hence this was judged to be a suitably large interval. As a check, yf = 50 was also used to find the solution for randomly selected cases and the results agreed to within 1%. A mesh was used which clustered points close to y = 0 in the y direction and which was uniform in the t direction. The asymptotic solution (6) was used to provide an initial solution and the calculation was continued until the conditions —*
(Sb) jG(y, 2nir) G(y, 2(n + 1)IT) I 10~, I H~(y,2n~r) H~(y,2(n + 1)i~)I 10~, —
We have employed two methods to solve eqs. (4) and (5). First we have conducted an asymptotic analysis appropriate to the limit N cc, the de—‘
—
were both satisfied for all y
E
[0, yf I for some
R.A. Cartwright, A.A. Wheeler / Axial magneticfield and periodic growth ratefluctuations
470
positive integer n, in which case the solution was deemed to be temporally periodic and to have converged. 2.3. The solute distribution The transport of solute is governed by the convection diffusion equation which is vc, (7) ac/at * = D v2c —
~.
with boundary conditions f (t *) C + DC = f (t *) C at the interface, *
*
(g\ “
C(z,
C as z
t*)
—~
00
cc,
/
(9)
where D is the solute diffusivity, C is the concentration and v = (Vr, v 9, v~)is the velocity. Eq. (8) expresses conservation of solute at the interface and (9) describes the uniformity of concentration in the bulk. We employ the following non-dimensionalization
c(y, -
t*)/C
T)C(Z,
In the second approach we employed the NAG subroutine DO3PBF to integrate numerically eqs. (10)—(12). We used an interval 0 y 50 and the boundary condition at infinity (12) was replaced by the flux condition, 1 ëY(YF, T) = 1, YF = so, Sc h y~’ in order to ensure that solute was conserved. The mesh was clustered close to y = 0 in the y direction to ensure that there was adequate resolution in the boundary layer. It was assumed that if ~ 20, h(y, r) = h(20, T). Convergence was ë(y~, T)
—
(
judged to have occurred when E(y,2nir/Sc)—E(y,2(n+1)njSc)
e( y for all y where
E
2( n + 1) i’-/Sc)
tol,
[0,YF] for some positive integer n,
to! = iO~ for 4 for tol=10
N 8, N<8.
00, where cot */Sc, Sc = v/D is the Schmidt number,
3. Results
T =
3.1. Introduction in which case eqs. (7)—(9) become =
ë +
[K’(’r)
—
Sc h(y,
T)]
~
,
(10)
with (1—k0) Sc a(r) E(y, at y e(y,
=
7)
7)
+~~(y,7) =0
0,
—~
1 as y
—*
cc,
In this section we present results for some specific cases. The model permits the use of any periodic growth rate but we have chosen a(t)
=a(i +A sin
(11)
which gives
(12)
1(t)
=
aA(1
—
t),
cos
where k0 is the equilibrium distribution coefficient, and
where A is a parameter which governs the amplitude of the oscillations. When A > 1 backmelting
a(’r)
=
h(y,
~)
occurs and so we have restricted ourselves to 0 A <1. In particular we have investigated two values of A: A = 0 and A = 0.9. When A = 0, the
a(t), =
K(T)
H(y,
=
1(t),
t).
Again we have adopted two approaches to solve eqs. (10)—(12). The first is an asymptotic analysis2 valid the details limits N cc, solution Sc cc, are withgiven N/Scin cc. inThe of this appendix C and are discussed in the next section. —
—~
-~
steady case is recovered and we compare our resultslarge withamplitude Cartwrightgrowth et al. [6]. The case without A = 0.9 gives oscillations backmeltmg and enables a comparison with Wilson [3,4] to be made.
R.A. Cartwright, A.A. Wheeler
/ Axial magnetic field and periodic growth ratefluctuations
Table 1 Physical properties of liquid silicon at the melting point 1410°C 3 Density p = 2.533 x iO~kg/rn Kinematic viscosity v 0.35 x 10—6 m2/s Electrical conductivity a = 1.33 x 106/I? m Magnetic permeability = 4s~)< iO~ H/m Molecular diffusivity D = 1.25 X 10-8 m2/s
The physical properties of silicon are given in table 1 and other typical values for various experimental parameters in table 2. These two tables provide the seven dimensionless parameter values relevant to the Czochralski growth of silicon which are given in table 3. There are six parameters that may varied since Prm is very small for all realistic situations and in our analysis O(Prm) terms are neglected. We have chosen to primarily investigate the effects of varying N and s as these are the parameters which may most easily be altered experimentally. The remaining parameters are fixed at the values given. The exception to this is that, in order to make a comparison with Cartwright et al. [6], we used A = 0 and in this case varied a in addition to N and s. As explained above the calculation splits into two stages. First the velocity field is calculated and then this is used in the concentration model to produce effective distribution coefficients. Accordingly, the results for the velocity field are given first in section 3.2, and then in section 3.3 we analyse the concentration field results,
Table 3 Dimensionless numbers associated with a typical Czochralski system Rotation ratio s = D/~: —0.1 to 1.0 Magnetic interaction parameter N = B~a/p,~: 0 to 32 Pulling rate a = j’/( s~)~‘~= 0.03 Oscillation amplitude A = 0.9 Equilibrium distribution coefficient k 0 = 0.07, 0.35 Schmidt number Sc = v/D 8.75, 20 Magnetic Prandtl number Prm = pa~s0= 6 X iO~.
out as N increases. Because the fluctuations are represented by terms like N 2j ‘(t), for moderate magnetic fields (N> 4), the flow is perturbed only slightly as a result of the fluctuating growth rate. In table 4 we compare results for the velocity field using the asymptotic expression for H with those from the numerical integration of the equations. We find that at N = 10 the results are in excellent agreement differing by 0.1% and for N = 4 the agreement is still good being with 1%. For N> 4 we conclude that the asymptotic expressions (6) give a good approximation to eqs. (4).
Table 4 A comparison of numerical and asymptotic results obtained for H(oo) with A = 0
_____________________________________________ s
Table 2 Experimental parameters for a typical Czochralski system
Magnetic fields Pulling rate
Crystal L, = 3 X 10 m 2 m Crystal = 1.05 to 2.1 rad/s CruciblewL0=7.5X10 Crucible U = —0.63 to —1.05 rad/s B 0 = 0.1 to 0.23T m/s / = 0.025 x itr
Asymptotic results — 0.3960 —0.6249 —0.4004 —0.2598
0 0.5 0.8 —0.1
01
03 0.6
—01658
—011584
—01692 —0.1633
—0.12
—0.1166
—0.1173
0 0.4 0.9 —0 13 0 0.2
—0.1100 — 0.1181 —0.1055 —0 1037 — 0.0203 — 0.0266
—0.1102 — 0.1184 —0.1056 —0 1040 — 0.0205 — 0.0268
0.7 —0.15
— —
0.0240 0.0128
— 0.0242 —0.01 30
N=10 a = 0.1
—2
Length scales Rotation rates
Numerical results — 0.3291 —0.3767 —0.2270 —0.2311
N= 1 a=0.1
3.2. The flow From the asymptotic analysis presented in appendix B we find that H(y, t)— —a as N—~cc and the fluctuations in the velocity are damped
471
N = 10 a 0.01
472
R.A. Cartwright, A.A. Wheeler / Axial magneticfield and periodic growth rate fluctuations
It should be noted that when s = 1, the solution to (4) is Gm 1, H —a. Physically this corresponds to the solid body rotation when the disc and the far field fluid are corotating at the same rate. Thus, in this case the velocity field is independent of the applied magnetic field. 3.3. The solute distribution Although the two parameters N and s do not appear explicitly in the mass transfer equation (10) or boundary conditions (11) and (12), they do affect the axial velocity and, therefore, alter the advection of solute. The concentration of impurity in the solid crystal at the interface is given by k0ë(0, T) where k0 is the equilibrium distribution coefficient. One important feature of the solution is the concentration at the interface. We have found in all our calculations that ë(0, r) is penodic (with period 2~r/Sc).In fig. 2, we plot ë(0, T) over one period for the case s = 0.2, N = 1. This case is typical and shows that ë(0, T) closely resembles a sine wave with a phase shift. A discrete Fourier transform of ë(0, ‘r) confirms this. The extreme points of ~(0, r) are at i’ —
8
1.44
c(o,t)
1.39
a(t) .34
/
°
3
1 .28
I.23~
(0.37)2ir/Sc and T (0.87)2i7’/Sc and appear to be independent of N and s. The solution “lags” behind the growth rate (which has extreme points at (0.25)2w/Sc and (0.75)2ir/Sc). This is in agreement with the results of Wilson [3,4]. The other two quantities used to describe the solution are the amplitude which is defined as A~= e(o, ‘T)max — e(0, ~ (13) and the grown-in effective distribution coefficient given as
J
r2,r/Sc
keff = k0
ë(0,
‘r) a(i’)
2or/Sc
Fig.2.Aplotofc(0,t)forthecaseA=0.9,d=0.03,s=0.2, N = 1.0, k = 0.07, and Sc = 8.75. The growth rate a(t) is given for reference.
(14)
ft( r) d T The results of our calculations have been cornpared to other numerical work for the special
cases A = 0 (given by Cartwright et al. [6]) and N = s = 0 (given by Wilson [3,4]). Cartwright et al. [6] present their results graphically for J as a function of i~ and 9 for different values of s, where ~
— —
~
1
(~ —
k0 \
—
3~4, O=(N+0.427+0.59s) Sc1”2. ~=àSc In table 5 we compare the values of J obtained by Cartwright et al. [6] with the results of our calculations. The two sets of figures are in good agreement. A comparison with Wilson [3,4] is more straightforward as we have adopted the same formulation as she. For k 0=0.01, ~=0.03, A =0.9 and Sc = 20, Wilson [4] reports that E(0, i’)m~ = 1.65 and ë(0, i’)~ = 1.34 (with an error estimate of 1%) and our calculations give ë(0, T)m~ = 1.63, ~(0, i’)~ = 1.31 (with an error estimate of 0.1%) which is in good agreement. The comparison was made by estimating values from a graph (fig. 1 of Wilson [4]) and so this may account for part of the slight discrepancy. All other calculations were conducted for A = 0.9 over therangeON32, —0.1 s 1.0. We consider three different cases:
~II
t
di’
(a) (b) (c)
k=0.07, k=0.35, k = 0.07,
Sc=8.75, Sc=8.75, Sc = 20.
R.A. Cartwright, A.A. Wheeler / Axial magneticfield and pen odic growth ratefluctuations
Case (a) corresponds to carbon in silicon and case (b) to phosphorus in silicon. These were the two cases considered by Series, Hurle and Barraclough [13] for the steady situation. Case (c) is included
~
0.5 S 0.4
CONTOUR ERr O.44E*OO O.50E+OO O.57E+OO
2 3
~
~
S
O.75E~OO
7
O.SIE+O0
1
05
05
10
2~0
2/~
7
\
/
I
K
~
~
\ \
//
0~b
to show the effect of changing the Schmidt nurnber. Contour plots of keff and the amplitude A~, are given in figs. 3 and 4. The first point of note is that the contour plots of keff and A~are indis-
CONIQURT
-°•:O
I 2
~
//4
473
I
I I
I I
I
S 0.4
/////)///1/J
CONTOUR EEl 1 O.I5E*OO 2 O.25E*OO 3 O.34E*OO
2
1
~ S
O.62E+OO
7
0.725+00 1
~
~
0
20~
10
~
Fig. 3. Contour plots of k~8~ in (N, s) space for: (a) k0 = 0.07, Sc = 8.75; (b) k0 = 0.35, Sc = 8.75; (c) k0 = 0.07, Sc = 20.0. In each case A = 0.9 and a 0.03. The N axis is logarithmic with the exception that the data at N = 0.25 are replaced by those for N = 0.0.
R~A.Cartwright, A.A. Wheeler / Axial magnetic field andperiodic growth ratefluctuations
474
Table S A comparison between the results of Cartwright et aL [6] and numerical results with A =0
8=10 0.15 k0 = 0.5
o
1 = 0.15 k0 = 0.5 =
8=1 =
Cartwright etal.[6]
—0.19 0 0.4
0.995 0.935 0.891
0.995 0.94 0.90
0.8 0.95
0.934 0.980
0.94 0.98
0.19 0 0.4 0.8 0.95
0.782 0.411 0.337 0.433 0.671
0.78 0.43 0.35 0.45 0.70
0.030 0.249 0.882
0.03 0.25 0.89
0.01 0.1 1
0.5 0.4
=
s
Numerical results
The contours in (N, s) space for c(0, by the equation (1—s)(1+5s) =D, N
t)
are given
aconstantasN~cc. (16)
Hence for large N the dependence of both keff and A0 on N and s will be of this form. Thus the asymptotic analysis supports the numerical evidence that the contour plots are indistinguishable. Eq. (16) also represents the shape of the contours well: the line s = 1 is always a contour and for a large fixed value of N the contours will have a turning point at s = 0.4. We shall now consider in greater detail the asymptotic expressions for keff and A0. From appendix C we have an asymptotic expression for A given by (C.10): (17)
‘
tinguishable. In all three cases an increase in keti is accompanied by an increase in A 0 and vice versa, The second observation is that when S = 1, H —a (which is independent of N). It follows that keff and A0 are, therefore, also independent of N when s = 1. This causes a rapid change in the solution as s —~1 when N is small but the variation more keff gradual that in allis cases and for A large N. We found 0 are minimum for N = 0, s and 0.05 as ands —~ that increase N increases 1. they When N is both large,asthe
where y = 1 + M/k0. We see that aA0/as = 0 when ay/as = 0, that is, at s = 0.4. A0 has a maximum for a fixed value of N, at this value of s, if dA0/dy > 0 and a minimum2)1/2 if dA5/dy <0. dA0/dy = 0 when y = 1 and dA +2 )1/2, (1 — ASubstituting for0/dy <0s = for0.4,y>we1 have + (1 — y with A
minimum of keff or A0 for a fixed N is at s 0.4. The asymptotic analysis given in appendix C supports these findings. From (C.9) we have
With a = 0.03, k0 = 0.07, A = 0.9, N 32 the condition is satisfied and so dA0/dy <0. This
(1 k0)(1 +A sin t) + O(Sc’) 1 + k0(1 + A sin t) + M 2 cc, as Sc—’ cc,N cc, with N/Sc (15) c(0,
—
t)
—
-*
—
where M= —s)(1 c(y, (1 t)=ë(y,
3a, +T),5s)/3A
A=[a+(a2~l~4N)1~~’2]/2.
as
2—A2) 2(1—k0)(y—1)A +O(Sc’) k0(y Sc —p cc, N —~cc, N/Sc2 — cc,
2~___~l(1_A2)~~” 0 d’y Is.-0.4
suggests that A0 is minimum when s = 0.4 for a results.value fixed It is of interesting N. This to agrees note with that for thesufficiently numerical large N (N> 75 with the above data), dA 0/ d’y I ~=o.~> 0 and A0 will then have a maximum for s = 0.4. We can also consider the asymptotic forni’of kCff as given by (Cli). k~
2 88=k0+(1_k0)[2_~+
(y—1)
(~2 _A2)1/2j
1
(18)
R.A. Cartwright, A.A. Wheeler / Axial magnetic field and pen odic growth ratefluctuations
475
0.9 0.9 0.7 0.6
1I
S0.5
II II I I
1 I O.27E*OO KE~ IL_CONTOUR 2~ 0.~SE*OOI I 51 3 4 O.83E*OO I 0.105*01 0.54E*O0
0
I
I
0.3
I I
sI 7 I S I
I
0.125*01 0.14E*O1
I
O.ISE*01 O.18E+01 0.19E--O1
9 10
I
2
I
0.2
/
0.1
1I
0
a
0.5
.0
2.0
4.0
8.0
16.
32.
N
_____
S
I 04 0 35
I
O.12E-l-OO
2JI 3
O.14E*OO!\ 0.ISE*00 I
~
0.17E*00
6 O.2OE*0O 5 0.196-1-00 7CONTOUR 0.225-1-00 8 0.245*00 KEY
I
0.2
0.1 0
,) ))) )
I
1
I
I i I I Ii II
3
0 27E 00 O.25E*O0
b 0.5
2 1.0
I
I
5
2I 3 I
I
0.4
I
7
I
/;~~2~
2.0
N
4.0
8.0
16.
32.
0.2 0.1 0 Od~Q
9 10
I
0.5
I
I
O.*3E*OO
O.7OE*0OI O.97E*OO 0.125*01 I
I
II
II
II
0.23E*01
I
I
0.25E*O1 0.295*01
I
SI 0.15E*O1 sII 0.18E*01 7CONTOUR 0.20E*01 KEY
0.3
2
*
___ 9
5
0.5
I
I
/1
I
/ /1
1.0
2:0
N
1.0
2
/
8.0
32) 16.
32.
Fig. 4. Contour plots of A0 in (N, s) space for: (a) k0 = 0.07, Sc = 8.75; (b) k0 = 0.35, Sc = 8.75; (c) k0 = 0.07, Sc = 20.0. In each case A = 0.9 and a = 0.03. The N axis is logarithmic with the exception that the data at N = 0.25 are replaced by those for N = 0.0.
It is easy to show from (18) that k0
<
1 and
N—.oo
This is in agreement with the behaviour of the
numerical results. It should be noted that the changes in keff and A0 with N and s are due to changes in the axial flow. A strong axial flow toward the disc decreases k0ff and A0, whereas with a weak axial flow they increase.
R.A. Cartwright, A.A. Wheeler / Axial magneticfield and periodic growth ratefluctuations
476
It is difficult to substantiate the above analysis with experimental results due to the problems of isolating specific factors in experiments. In practice the effect of the magnetic field on striations is mixed. Series [14] reports that for a crystal grown in a 2000 G axial field the concentration fluctuations were 2% at the centre and 4% at the edge of the crystal. This compared with a control crystal grown in the absence of a magnetic field with corresponding fluctuations of 4% and 3%. The effect of the magnetic field is, thus, inconclusive, In reality, of course, the magnetic field also interacts with many other parts of the system. It has an effect on temperature fluctuations in the melt and hence may modify a ( t), the fluctuating growth rate which in our model is independent of N. We have assumed a constant pull rate, a, but this too may oscillate. Either of these modifications, once quantified could be incorporated into our model to provide a more sophisticated analysis. This paper provides insight into the effect of magnetic fields in fluctuating growth conditions and gives a basis for future developments.
4. Conclusions We have investigated the effect of a magnetic field and crucible rotation on the fluid flow and mass transfer in the vicinity of the crystal/melt interface. The axial flow weakens as N cc and the velocity boundary layer,For which is oflarge thickness 1”2) becomes thinner. a fixed value O(N of N the axial flow is maximum when s = 0.4. When s = 1 we have solid body rotation. Both the numerical and asymptotic results suggest that as N cc, keff i. For a fixed large value of N, keff is minimum at s = 0.4 and keff ~ 1 as s 1. The numerical results show that when 8 N 32, A 0 also has a minimum at s = 0.4 and increases as s 1, for a fixed value of N. An asymptotic analysis supports this and in addition indicates that for a large enough value of N (N> 75 for typical growth data), A0 will have a maximum at .s = 0.4 and will decrease when kCff increases. The phase of the concentration oscillations appears to be independent of N and s. —~
—‘
—‘
—~
~
AcknowLedgement One of us (R.A.C.) gratefully acknowledges the receipt of an SERC grant during the the period when this research was conducted.
Appendix A Due to the coupling between the magnetic and velocity fields, the fluctuations in the velocity field induce time dependence of the magnetic fields. However, in this appendix we show that, in the limit Prm 0, such time dependence is confined to terms of O(Prm) and hence negligible since Prm iO~ in realistic situations. Therefore, it is possible to replace the magnetic field by the steady applied field. The electromagnetic field is governed by Maxwell’s equations which are —‘
—
e0V’ E = Pc’
.
(A.i)
V’ B = 0, ~7xE= _8B/at*,
(A.2)
VXB=~t0(J+c0 8E/at*).
(A.4)
(A.3)
Ohm’s law gives J=
~
(E + ~ x
B),
(A.5)
and we wish to find the Lorentz force F= i.t>< B p
~A 6
‘
where E is the electric field, B is the magnetic field, J is the current density vector, p.0 is magnetic permeability, p is density, Pc is free charge density and e~ is permittivity of vacuum. Eqs. (A.3) and (A.5) give that ~ ~<
,
= ~
( \
+ ~
—
at
*
>
(v
x B )‘~ 1.
(A 7)
R.A. Cartwright, A.A. Wheeler / Axial magneticfield and periodic growth ratefluctuations
J from (A.4) to give, aE\
Substitute for fi
477
These equations can be non-dimensionalised as follows: -
1/2
~
Po
aB
(A.8)
j, \1/2
We also have the equation of continuity for the fluid, V. v = 0, (A.9) and substituting (A.2), (A.3) and (A.9) into the 0-component of (A.8) gives, in cylindrical polar coordinates with axial symmetry,
2
+
az
ar
-~—
—p-ar- (rBe))
r
a B9 r
r8r( —Tv—
—
)J +
a
0
0.
Similarly of (A.8) together with (A.2) and the (A.3)z-component gives, i
a
-,-
-~-,
2B 2\ ~ 2 ~~0a a ~ r-~-,-)+ -j—~-+ -~- (rv2 Br I
~B
—,-—
2B
aB 1at” 0a—
a
+
+ p.
=
—
rv~B~)
0.
(A.i1)
It should be noted that these equations are the same as (14) and (15) obtained by Cartwnght et al. [6] together with additional time dependent terms. The Lorentz force, F, is given by 1
~B
aB
F—— p.p ( ~ z
az
az
aE
—
=
~(-~---
0
strength, 02 is the angular velocity of the disc and p is kinematic viscosity. If we rewrite (A.10) in non-dimensional form it is seen that
a = 1.33 x
B B9
+ —fl— r +
~az — Braz aB,. —B9— — s~0p~
i0_6/~.
Bz
~B9
az
‘
m. 2/N m2,
o = 8.854 x 10—2 C and so e 5. 0W/a=1.398Xi0 So, as in the steady case we obtain L~ O(Prm) and to this order the additional time dependent terms provide no extra contribution. A similar argument applies to the time dependent terms in (A.11). —
~
The time dependent term in the body force is, >< B
at
*
/
From (A.3), I E I = O(B0w1) scale. So the magnitude of thiswhere term /is is a length B0c.,l
O~O—T~-—WBO) =
( at~ 8E xB aE
B~
For silicon, w=1.O5~2.1rad/s
I XB)),
1 I 8B9 — p.0p I\ Br — ar C0
aB ~
=
B0 is the externally applied magnetic field
~B
e0p.0~0a~
(~t)
=
~
where
L
—
L~+Prm(NGf—HLf)+Pr~4——Lii+Li)=0. I (002
—)
2.B~
(A.10)
~
~‘
8B
7/3Z v +p.0aIB~_ 2—+Br az9 rar
F
/ j~\1/2 r ~wi
o(NW2l(~~~)).
The steady terms in the body force are of magnitude, NW/alp.0,
+
3B Brar —i
—~(~x~)J.
and so the ratio of the time dependent terms to the steady terms is (~w/a) 12 (W/P) Prm, whichisnegligible.
478
R.A. Cartwright,
A.A.
Wheeler / Axial magnetic field and pen odic growth rate fluctuations
Appendix B
Expand g and h as power series in X1
In this appendix an asymptotic solution valid for large magnetic fields is developed. We have from (4)
g(x,
00
H~,H,,,,~+H,,,,(I’(t)_H)+~H~ —2(G2 — s2)
—
NH,,,
(B.1)
G~= G~,,+ GH~+ G~(I’(t)— H)
—
N(G
— s),
(B.2) with boundary conditions H(y,
t)
=
—a,
H,,(y, at
t)
0,
G(y,
t)
=
t)
t) —
y—’cc.
~, (B.3)
Let X=
(B.4)
[~+(~2+4N)h/’2]/2,
and put
t)
—à+h(x, t)/X3,
G(y, t)=g(x,
0
=
h0
=
1
=
(B.9) t).
(1 — s) ex + s, 2(1—s) ~ (1 + 2s + 3sx) e_x (1 —s)(1 + 5s) — 3
—
(1
2x
2
s) e (B.10)
=
—s(i
—
(B.ii)
s) I’(t),
—(i—s) ~(i+5s) 3
I’(t)—3sä],
—(1—s) r( 7+115) I’(t)—d(2+7s)], 9 (1 — s) I’(t)
$ 5_i~(1—s)[7I(t)—2a].
~
Eqs. (B.10) and (B.il) can be substituted into
g~~—(g—s) 1[g~(I’(t) = —A
(B.6)
(B.9) to give expressions for g(x, t) and h(x, t) and then using (B.5) we obtain G(y, t) and H(y, t).
(B.7)
AppendixC
+a(g—s)]
+X2g~—X4(g~h—gh~),
with boundary conditions t)
x’~hr(X,
2 I’(t) 2+f3 x e_x h1=/3o+($1x 2x+$3)e~ 2x, +(/34x+ 135) e_ where ~o= ~
/3 =
~
h (x,
~
=
~=
Eqs. (B.1)—(B.3) become 2—s2) h~~~—h~—2(g
=
—(1
(B.5)
t).
t)
g
~
=
t),
and substitute in (B.6) and (B.7) to obtain
x=Xy, H(y,
~ X~g,.(x, 00
h(x,
g 0, G(y,
=
r=0
1,
y=O
H,,(y, as
=
t)
=
0,
h~(x,
t)
=
In this appendix we conduct an asymptotic 0,
at x = 0, h~(x, t) —p0, as x
—b
cc.
g(
~,t)
g(x,
t)
=
1,
analysis of the concentration which supports the results produced numerically. We must solve eqs. (10)—(12),
—‘s, (B.8)
C~=
+
[I’(t)
—
H(y,
t)]
c,,,
(C.i)
R.A. Cartwnight, A.A. Wheeler
/ Axial magneticfield and periodic growth ratefluctuations
with boundary conditions,
with
(1—k0) Sc a(t) c(0, 1)+c,,(0, t)=0,
(C.2)
c(y, t)
(C.3)
—
1, as y
—
cc,
where c(y, t)
=
e(y, r).
t)
3(i
— —a + ~A —(1
—
—
s) e2AJ
s)[2(1 —
+
(1 + 5s)]
+
—
_~3(a +
3 + 5s) ) + O(A~), (1 —s)(i cc, with A/Sc — cc.
—
c”(3, t)=c
2c~’(.P, t)+O(Sc2), 0”(j, t)+Sc and find from eqs. (C.6)—(C.8) that at leading order ~
t)
1
=
+
(1
—
k
0)ä(1 + A sin t)
k0d(i +A sin t) + K
reveals that H(y, t)
sin t+â+K)j5], xexp[—A where K= (1 — s)(1 + 5s)/3. This gives 3(aA
—A3{d— -~(1— s)[2(1 + 2s+ 3sAy) e_AY 2XY_(1
—(1 —s) e asA—’cc,
3~.
+
5s)]}
where a =table A 3 we see that Sc From
+O(A4),
— t)
c”(O, t) —1
(C.4)
— 8.75
—‘
cc,
cc. We expand c(y, t) in powers of Sc’ 2c 3). (C.5) = c0 + Sc~c1+ Sc 2 + O(Sc
Then eq. (C.1) and boundary condition (C.3) give c0(y, t) = 1, c.(y, t) = 0 for i 1, which does not satisfyHence the boundary condition at the interface. (C.5) represents the(C.2) concentration in an outer region where y is 0(1) and we require an inner region where y = O(Sc~).We put j3 = Sc ~ in which case eqs. (C.1) and (C.2) give the following problem for c” which represents the concentration in this region
+
(1— k k 0)a(1 +A s~nt) 0a(1 +A sin t) +K
+O(Sc’). Hence
and so it is
appropriate to consider the realistic limits Sc
c(y,
Eq.
We put where A = [~+ (a2 + 4N)1~2]/2. For a typical crystal growth system a 0.03 and so the first two terms of the expansion are of comparable magnitude. If we assume a — A ~ then a similar analysis to that given in appendix B
A/Sc
t).
H*(j~,r) as A — cc, Sc
O(A4)
(C.8)
In this inner region we denote H as H”(j3, (C.4) gives
2s + 3sAy) e~
as A—’cc,
—
(1_k0) a(t) c*(j~,t)+c~’(y, t)=0,at ~3=0, (C.7) along with the matching condition c”(33, t)—i1,asy---cc.
In appendix B we derived an asymptotic form for H(y, t), H(y,
479
(C.9)
(1 —k 0)d(1 +A)
c, 0~= 1 + Cmth =
1
+
k0â(i —A) + K (1 — k k 0)d(1 — A)
0d(1 —A)+K
3), + O(Sc 1), +
O(Sc
which gives 1), (C.10) A= 2(1 k— k0)(y — i)A + O(Sc C 2 — A2) 0(’y where ‘y = 1 + K/k 0 a. Further, eqs. (C.9) and (14) give, [2 + (y — 1)2 kCff = k0 + (1 — k0) —
___________ (72_A2)1/2]’
~-c7 =c;~+(J’(t) -H(j~, t))c;,
(C.6)
(C.11)
480
RA. Cartwright, A.A. Wheeler / Axial magneticfield and periodic growth ratefluctuations
References [1] K.S. Choe and T.H. Strudwick, J. Crystal Growth 74 (1986) 338. [2] H.A. Chedzey and D.T.J. Hurle, Nature 210 (1966) 933. [3] L.O. Wilson, J. Crystal Growth 48 (1980) 435. [4] L.O. Wilson, J. Crystal Growth 56 (1982) 219. [5] A.A. Wheeler, Proc. Roy. Soc. (London) A379 (1982) 305. [6] R.A. Cartwright, D.T.J. Hurle, R.W. Series and J. Szekely, J. Crystal Growth 82 (1987) 327. [7] L. Hjellming and J. Walker, J. Fluid Mccli. 164 (1986) 237. [8] L. Hjellming and J. Walker, J. Crystal Growth 83 (1987) 51.
[9] L. Hjellming and J. Walker, J. Fluid Mech. 182 (1987) 335. [10] R.A. Cartwright, N. El-Kaddah and J. Szekely, IMA J. AppI. Math. 35 (1985) 175. [11] J.F. Brady and L. Durlofsky, J. Fluid Mech. 175 (1987) 363. [12] T. von Karman, Z. Angew. Math. Mech. 1 (1921) 233. [13] R.W. Series, D.T.J. Hurle and K.G. Barraclough, ImA J. Appl. Math. 35 (1985) 195. [14] R.W. Series, privative communication, 1987.
Journal ofCrystal Growth 87 (1988) 466—480 North-Holland, Amsterdam
THE EFFECT OF AN AXIAL MAGNETIC FIELD ON PERIODIC GROWTH RATE FLUCFUATIONS IN CZOCHRALSKI CRYSTAL GROWTH R.A. CARTWRIGHT and A.A. WHEELER School of Mathematics, University of Bristol, University Walk, Bristol BS8 ITW, UK Received 15 June 1987; manuscript received in final form 25 September 1987
In this paper we develop a mathematical model to investigate the effect of growth rate fluctuations in magnetic Czochralski crystal growth. The pull rate is assumed constant and the growth rate fluctuates with frequency a rad/s, the crystal rotation rate. We explore the effect of the magnetic field and the crucible rotation on fluid flow and mass transfer in the vicinity of the crystal/melt interface. Numerical results are obtained for 0 N 32 and —0.1 < s 1.0, where N is the magnetic interaction parameter and s is the ratio of crucible to crystal rotation. Asymptotic results are obtained for the velocity field for large N and agreement is good for 2 —~oo where Sc is the N> 4. The analysis is extended to the concentration givethe asymptotic for NWe —* ~o, Sc —, c~,N/Sc Schmidt number. Again the asymptotic results agree field wellto with numericalresults solutions. conclude that the phase of the resulting concentration fluctuations appears to be independent of N and s. The spatially averaged effective distribution coefficient, k~ffand the amplitude of the concentration oscillations, A~,are linked. As keff -~ 1, A~increases and so improved macrosegregation is accompanied by a deterioration in microsegregation.
1. Introduction It is well known that crystals grown by the Czochralski technique may exhibit striations in which a spatially periodic variation of dopant concentration is observed in the direction of the axis of the crystal. This phenomenon has been attributed to periodic temporal fluctuations in the crystal growth rate, caused by thermal asymmetries in the melt (in which case they will have the same period as that of the rotating crystal) or unsteady thermal convection in the bulk of the melt, (Choe and Strudwick [1]). This paper investigates the effect of an axial magnetic field and crucible rotation on the flow near the crystal/melt interface and the concentration of impurity there. In particular we consider whether either of these two control mechanisms will help to damp out the fluctuations in concentration and so enable crystals with a more uniform distribution of impurity to be grown. In recent years semi-conductor chips have been progressively reduced in size, whilst at the same time device performance requirements have become
more stringent. These two facts result in the necessity for better crystal quality; a uniform distribution of impurities at prescribed levels. A magnetic field has been suggested (first by Chedzey and Hurle [2] in 1966) as a means of achieving this aim since it damps the flow and hence may reduce fluctuations in the system. Previous work has focussed on two areas of the problem. Wilson [3,4] and Wheeler [5] have considered the case of a periodic growth rate in the absence of a magnetic field and without crucible rotation. Cartwright, Hurle, Series and Szekely [6] studied the effect of magnetic fields and crucible rotation when the crystal grows at a constant rate. We consider here a more general model which includes a periodic growth rate, an axial magnetic field and crucible rotation. The work of the above investigators are special cases of our model. All of these workers have confined their attention to the boundary layer adjacent to the crystal. Other important work in this area is by Hjellming and Walker [7—9]who have conducted an analysis of the magnetohydrodynamic flow in the entire melt. Their model assumes that the imposed uniform
0022-0248/88/$03.50 © Elsevier Science Publishers B.V. (North-Holland Physics Publishing Division)
R.A. Cartwrighf, A.A. Wheeler / Axial magnetic field and periodic growth ratefluctuations
magnetic field is sufficiently large that inertia effects may be neglected. This has the mathematical advantage that the resulting governing equations are linear. The major insight of their work is to show that in the presence of such a high magnetic field the electrical conductivity of the crystal, although small, has a profound effect on the flow in the melt. It should be noted that our model does not incorporate the role played by natural convection or surface tension. The effects due to the memscus may be neglected because in the vicinity of the crystal/melt interface the flow is dominated by forced convection. A detailed argument based on the dimensionless groups of the system is given in Cartwright, El-Kaddah and Szekely [10]. With these assumptions the flow decouples from the temperature field and so we consider an isothermal melt. The resulting problem consists of two distinct parts: (i) the flow, (ii) the solute distribution. The first stage is to solve for the fluctuating velocity field and then in the second this solution is substituted into the mass transfer equations in order to determine the concentration field. In section 2 we derive our model and describe our solution method for each part of the problem and in section 3 we present our results. We find that the crucible rotation and the magnetic field strength have a profound effect on the effective distribution coefficient. In general terms, as the magnetic field strength increases the average effective distribution coefficient, kC~1, approaches unity. However, the amplitude of the concentration oscillations also increases and so although a magnetic field improves macrosegregation, microsegregation may not be beneficially affected. The effect of crucible rotation is more complicated but a similar correlation between increased keff and increased amplitude of concentration oscillations is observed. In practice, of course, the magnetic field also has an effect on the thermal asymmetries in the melt which cause the growth rate fluctuations. Therefore, in order to argue conclusively about the effect of magnetic fields on crystal quality it is necessary to consider the complete system. This paper demonstrates the effect of the field on a prescribed growth rate fluctuation.
467
2. The model
2.1. Introduction The result of an axial magnetic field will be to modify the flow in the melt. Here we are concerned with the consequent effect on the grown-in concentration in the crystal, which is determined by the concentration at the crystal/melt interface. Hence we shall consider the flow and solute in the boundary layer adjacent to the growing crystal. Our model consists of an infinite rotating disc representing the crystal/melt interface, which rotates with angular velocity w rad/s; the crystal rotation rate. We justify our model of the crystal as an infinite rotating disc by noting that a recent study by Brady and Durlofsky [11] of the flow adjacent to a finite rotating disc has revealed that the effect of the edges of a finite disk upon the flow in the momentum boundary layer, are, for moderate or large Reynolds number, confined to the outer half of the disc. The flow in the central half is qualitively described by the flow due to an infinite rotating disc over a wide range of conditions. It is the region close to the centre of the crystal to which our model is applicable. We locate our coordinate system coincident with the initial position of the crystal/melt interface which is given by the plane z = 0. The configuration of the system is shown in fig. 1.
e
I
‘i~y w
~
DISC ~
ICwt~)
~ w)
F u id
1.1 SW Fig. 1. Model configuration.
R.A. Carzwright, A.A. Wheeler / Axial magnetic field and periodic growth ratefluctuations
468
To model the rotation of the crucible the fluid is assumed to have a prescribed solid body rotation rate of ~ rad/s far from the crystal, that is, as z oo. The pulling rate of f rn/s is modelled by imposing “suction” in the axial direction at the rotating disc and the fluctuating growth rate, f * (t *) by an oscillatory motion of the disc about z = 0 in the axial direction. We assume that the density of the solid and liquid phases are the same. Hence, in this frame of reference, the fluctuating growth rate does not cause a fluctuating suction velocity at the interface (see Wilson [4]). We use a cylindrical polar coordinate system (r, 9, z) as shown in fig. 1, and we assume that the system is axisymmetric. A magnetic field is applied in the axial direction and because the fluid is electrically conducting the applied magnetic field and velocity field interact to induce electromagnetic effects. In general these effects are magnetic field line distortion due to the velocity field and an inhibition of motion across field lines. Although the applied magnetic field is steady, the velocity field is time dependent and so the induced electromagnetic field will in general also be
parameter defined in the next section and t~ is the radius of the crystal. We now consider the first part of the model; the flow.
—*
time dependent. However, we note that for a silicon melt the magnetic Prandtl number, Prm (= yap.0) where v is kinematic viscosity, ~ is electrical conductivity and p.~is magnetic permeability is 6 x iO~. For similar situations relevant to crystal growth Prm is also small. Thus we shall consider our model in the realistic limit Prm 0. In appendix A we show that if O(Prm) terms are neglected, the magnetic field is steady and is not distorted by the velocity field. Hjellming and Walker [7] showed that the electrical conductivity of the crystal, despite being small, can for large magnetic field be important. We shall ignore this effect in our model and assume the crystal to be a perfect electrical insulator. As noted by Cartwnght, Hurle, Series and Szekely [6] in a similar situation, we must therefore assume that —‘
2<1, N1~~’2o~(co/v)I~’
where v is the kinematic viscosity of the melt, is the ratio of the electrical conductivity of the crystal to that of the melt, N is the magnetic interaction
2.2. The flow We model the liquid phase by a viscous, incompressible, electrically conducting fluid. Thus the governing equations are the time dependent Navier—Stokes equations with an electromagnetic body force term (the Lorentz force) and the continuity equation. We seek an axisymmetric solution and assume that the axial velocity and axial magnetic field are independent of r, in which case the following similarity transformation is appropriate (von Karman [12]):
( ±~) 1/2
~
z
=
,
~(
=
~
_~-
=
cot ~,
~(~
=
wr ~
= 1/2’
‘
where u~,v9, v2 are the dimensional velocity cornponents in the radial, azimuthal and axial directions, v is the kinematic viscosity, t is dimensional time and overbar denotes a dimensionless quantity. The equations of motion then reduce to the following set of partial differential equations: ~ = ~ ‘ii~ —2(G~ 52’l NH~ ~ 2 ~ (la) = ~+—~—N(~—s), (ib) *
—
—
/
with boundary conditions —
-
-
H(~, t)
—
=
-
-
~
~
at —
~
-
-
t)
— —*
0,
-
—
t)
=
=
I(
0,
~)
-
-
G(~,t)
=
1, (2a)
-
G(~,t) —‘s, -
as~—’oo, (2b) where 1(1) is the position of the interface at time t, a = f/( i’w )1~~2is the dimensionless pulling rate, N = B~a/pwis the magnetic interaction parameter, -
R.A. Cartwright, A.A. Wheeler / Axial magneticfield and periodic growth ratefluctuations
B0 is the applied magnetic field strength, p is the density of melt, and terms O(Prm) have been neglected. It should be noted that the growth rate, a(~),is given by a(t)=a+I’(t). In all of the calculations presented in section 3, a sinusoidal growth rate was chosen but the model allows for any periodic function of time. Although F does not appear explicitly in the equations it is easily obtained using the continuity equation which gives
~
t)=
~).
~
tails of which are given in appendix B. The solution obtained is H(y, —
2+11 ($1x 2x+/13) 5), e~ (6a) +($4x+/15) e_2x1 + O(X G(y, t) s + (1 s) e_x ~ (1 s) I’(t)x e_x + O(X2),
0+
—
(6b)
(3)
al. [6,10]. As shown by Wilson [3], it is convenient to define a new spatial coordinate, y, which moves with the interface, but to consider velocities with respect to the fixed coordinate system. Thus we define 1(i),
G(y,
t)
=
t
=
G(~,t),
H(y,
H(~, t).
t)
Hence eqs. (1) and (2) become =
~ —
+ H~~(I’(t) H) + ~H2 —
2(G2
G~= ~
+
—
~2)
—
NH~,
(4a)
GH~+ G~(I’(t) H) —
—
N(G
—
s), (4b)
with H(y,
t)
=
—a,
—
—
as A
~—
t) 3(1—s)[2(1+2s+3sx)e~ —a+~A _(1_s)e_2x_(1+5s)]
_x4[s
The steady version of eqs. (1) and (2) is the same set of equations studied in detail by Cartwright et
=
469
H~(y, t)
0,
G(y,
t)
at y = 0, H~(y,t) 0, as y—*cc
G(y,
t)
=
—~
=
1, (5a)
~
~,
—~
cc,
where
a + (a~+ 4N )1/2 =
,
2
x=Ay,
and /3~flu..., $~ are coefficients which are functions of time given in appendix B. This asymptotic expansion shows clearly the structure of the solution in which solid body rotation (H( y, t) = a, G(y, t) = s) is combined with a boundary layer of thickness O(N”~’2). The second method we have used is a numerical integration of the equations. To do this we employed the method of lines and Gear’s method by using NAG subroutine DO3PBF. This routine —
solves the equations on a finite domain in y and so the boundary conditions at infinity were imposed at yf = 20. The values of G and H~at y = 7 were within 0.01% of the value required of them as y cc. Hence this was judged to be a suitably large interval. As a check, yf = 50 was also used to find the solution for randomly selected cases and the results agreed to within 1%. A mesh was used which clustered points close to y = 0 in the y direction and which was uniform in the t direction. The asymptotic solution (6) was used to provide an initial solution and the calculation was continued until the conditions —*
(Sb) jG(y, 2nir) G(y, 2(n + 1)IT) I 10~, I H~(y,2n~r) H~(y,2(n + 1)i~)I 10~, —
We have employed two methods to solve eqs. (4) and (5). First we have conducted an asymptotic analysis appropriate to the limit N cc, the de—‘
—
were both satisfied for all y
E
[0, yf I for some
R.A. Cartwright, A.A. Wheeler / Axial magneticfield and periodic growth ratefluctuations
470
positive integer n, in which case the solution was deemed to be temporally periodic and to have converged. 2.3. The solute distribution The transport of solute is governed by the convection diffusion equation which is vc, (7) ac/at * = D v2c —
~.
with boundary conditions f (t *) C + DC = f (t *) C at the interface, *
*
(g\ “
C(z,
C as z
t*)
—~
00
cc,
/
(9)
where D is the solute diffusivity, C is the concentration and v = (Vr, v 9, v~)is the velocity. Eq. (8) expresses conservation of solute at the interface and (9) describes the uniformity of concentration in the bulk. We employ the following non-dimensionalization
c(y, -
t*)/C
T)C(Z,
In the second approach we employed the NAG subroutine DO3PBF to integrate numerically eqs. (10)—(12). We used an interval 0 y 50 and the boundary condition at infinity (12) was replaced by the flux condition, 1 ëY(YF, T) = 1, YF = so, Sc h y~’ in order to ensure that solute was conserved. The mesh was clustered close to y = 0 in the y direction to ensure that there was adequate resolution in the boundary layer. It was assumed that if ~ 20, h(y, r) = h(20, T). Convergence was ë(y~, T)
—
(
judged to have occurred when E(y,2nir/Sc)—E(y,2(n+1)njSc)
e( y for all y where
E
2( n + 1) i’-/Sc)
tol,
[0,YF] for some positive integer n,
to! = iO~ for 4 for tol=10
N 8, N<8.
00, where cot */Sc, Sc = v/D is the Schmidt number,
3. Results
T =
3.1. Introduction in which case eqs. (7)—(9) become =
ë +
[K’(’r)
—
Sc h(y,
T)]
~
,
(10)
with (1—k0) Sc a(r) E(y, at y e(y,
=
7)
7)
+~~(y,7) =0
0,
—~
1 as y
—*
cc,
In this section we present results for some specific cases. The model permits the use of any periodic growth rate but we have chosen a(t)
=a(i +A sin
(11)
which gives
(12)
1(t)
=
aA(1
—
t),
cos
where k0 is the equilibrium distribution coefficient, and
where A is a parameter which governs the amplitude of the oscillations. When A > 1 backmelting
a(’r)
=
h(y,
~)
occurs and so we have restricted ourselves to 0 A <1. In particular we have investigated two values of A: A = 0 and A = 0.9. When A = 0, the
a(t), =
K(T)
H(y,
=
1(t),
t).
Again we have adopted two approaches to solve eqs. (10)—(12). The first is an asymptotic analysis2 valid the details limits N cc, solution Sc cc, are withgiven N/Scin cc. inThe of this appendix C and are discussed in the next section. —
—~
-~
steady case is recovered and we compare our resultslarge withamplitude Cartwrightgrowth et al. [6]. The case without A = 0.9 gives oscillations backmeltmg and enables a comparison with Wilson [3,4] to be made.
R.A. Cartwright, A.A. Wheeler
/ Axial magnetic field and periodic growth ratefluctuations
Table 1 Physical properties of liquid silicon at the melting point 1410°C 3 Density p = 2.533 x iO~kg/rn Kinematic viscosity v 0.35 x 10—6 m2/s Electrical conductivity a = 1.33 x 106/I? m Magnetic permeability = 4s~)< iO~ H/m Molecular diffusivity D = 1.25 X 10-8 m2/s
The physical properties of silicon are given in table 1 and other typical values for various experimental parameters in table 2. These two tables provide the seven dimensionless parameter values relevant to the Czochralski growth of silicon which are given in table 3. There are six parameters that may varied since Prm is very small for all realistic situations and in our analysis O(Prm) terms are neglected. We have chosen to primarily investigate the effects of varying N and s as these are the parameters which may most easily be altered experimentally. The remaining parameters are fixed at the values given. The exception to this is that, in order to make a comparison with Cartwright et al. [6], we used A = 0 and in this case varied a in addition to N and s. As explained above the calculation splits into two stages. First the velocity field is calculated and then this is used in the concentration model to produce effective distribution coefficients. Accordingly, the results for the velocity field are given first in section 3.2, and then in section 3.3 we analyse the concentration field results,
Table 3 Dimensionless numbers associated with a typical Czochralski system Rotation ratio s = D/~: —0.1 to 1.0 Magnetic interaction parameter N = B~a/p,~: 0 to 32 Pulling rate a = j’/( s~)~‘~= 0.03 Oscillation amplitude A = 0.9 Equilibrium distribution coefficient k 0 = 0.07, 0.35 Schmidt number Sc = v/D 8.75, 20 Magnetic Prandtl number Prm = pa~s0= 6 X iO~.
out as N increases. Because the fluctuations are represented by terms like N 2j ‘(t), for moderate magnetic fields (N> 4), the flow is perturbed only slightly as a result of the fluctuating growth rate. In table 4 we compare results for the velocity field using the asymptotic expression for H with those from the numerical integration of the equations. We find that at N = 10 the results are in excellent agreement differing by 0.1% and for N = 4 the agreement is still good being with 1%. For N> 4 we conclude that the asymptotic expressions (6) give a good approximation to eqs. (4).
Table 4 A comparison of numerical and asymptotic results obtained for H(oo) with A = 0
_____________________________________________ s
Table 2 Experimental parameters for a typical Czochralski system
Magnetic fields Pulling rate
Crystal L, = 3 X 10 m 2 m Crystal = 1.05 to 2.1 rad/s CruciblewL0=7.5X10 Crucible U = —0.63 to —1.05 rad/s B 0 = 0.1 to 0.23T m/s / = 0.025 x itr
Asymptotic results — 0.3960 —0.6249 —0.4004 —0.2598
0 0.5 0.8 —0.1
01
03 0.6
—01658
—011584
—01692 —0.1633
—0.12
—0.1166
—0.1173
0 0.4 0.9 —0 13 0 0.2
—0.1100 — 0.1181 —0.1055 —0 1037 — 0.0203 — 0.0266
—0.1102 — 0.1184 —0.1056 —0 1040 — 0.0205 — 0.0268
0.7 —0.15
— —
0.0240 0.0128
— 0.0242 —0.01 30
N=10 a = 0.1
—2
Length scales Rotation rates
Numerical results — 0.3291 —0.3767 —0.2270 —0.2311
N= 1 a=0.1
3.2. The flow From the asymptotic analysis presented in appendix B we find that H(y, t)— —a as N—~cc and the fluctuations in the velocity are damped
471
N = 10 a 0.01
472
R.A. Cartwright, A.A. Wheeler / Axial magneticfield and periodic growth rate fluctuations
It should be noted that when s = 1, the solution to (4) is Gm 1, H —a. Physically this corresponds to the solid body rotation when the disc and the far field fluid are corotating at the same rate. Thus, in this case the velocity field is independent of the applied magnetic field. 3.3. The solute distribution Although the two parameters N and s do not appear explicitly in the mass transfer equation (10) or boundary conditions (11) and (12), they do affect the axial velocity and, therefore, alter the advection of solute. The concentration of impurity in the solid crystal at the interface is given by k0ë(0, T) where k0 is the equilibrium distribution coefficient. One important feature of the solution is the concentration at the interface. We have found in all our calculations that ë(0, r) is penodic (with period 2~r/Sc).In fig. 2, we plot ë(0, T) over one period for the case s = 0.2, N = 1. This case is typical and shows that ë(0, T) closely resembles a sine wave with a phase shift. A discrete Fourier transform of ë(0, ‘r) confirms this. The extreme points of ~(0, r) are at i’ —
8
1.44
c(o,t)
1.39
a(t) .34
/
°
3
1 .28
I.23~
(0.37)2ir/Sc and T (0.87)2i7’/Sc and appear to be independent of N and s. The solution “lags” behind the growth rate (which has extreme points at (0.25)2w/Sc and (0.75)2ir/Sc). This is in agreement with the results of Wilson [3,4]. The other two quantities used to describe the solution are the amplitude which is defined as A~= e(o, ‘T)max — e(0, ~ (13) and the grown-in effective distribution coefficient given as
J
r2,r/Sc
keff = k0
ë(0,
‘r) a(i’)
2or/Sc
Fig.2.Aplotofc(0,t)forthecaseA=0.9,d=0.03,s=0.2, N = 1.0, k = 0.07, and Sc = 8.75. The growth rate a(t) is given for reference.
(14)
ft( r) d T The results of our calculations have been cornpared to other numerical work for the special
cases A = 0 (given by Cartwright et al. [6]) and N = s = 0 (given by Wilson [3,4]). Cartwright et al. [6] present their results graphically for J as a function of i~ and 9 for different values of s, where ~
— —
~
1
(~ —
k0 \
—
3~4, O=(N+0.427+0.59s) Sc1”2. ~=àSc In table 5 we compare the values of J obtained by Cartwright et al. [6] with the results of our calculations. The two sets of figures are in good agreement. A comparison with Wilson [3,4] is more straightforward as we have adopted the same formulation as she. For k 0=0.01, ~=0.03, A =0.9 and Sc = 20, Wilson [4] reports that E(0, i’)m~ = 1.65 and ë(0, i’)~ = 1.34 (with an error estimate of 1%) and our calculations give ë(0, T)m~ = 1.63, ~(0, i’)~ = 1.31 (with an error estimate of 0.1%) which is in good agreement. The comparison was made by estimating values from a graph (fig. 1 of Wilson [4]) and so this may account for part of the slight discrepancy. All other calculations were conducted for A = 0.9 over therangeON32, —0.1 s 1.0. We consider three different cases:
~II
t
di’
(a) (b) (c)
k=0.07, k=0.35, k = 0.07,
Sc=8.75, Sc=8.75, Sc = 20.
R.A. Cartwright, A.A. Wheeler / Axial magneticfield and pen odic growth ratefluctuations
Case (a) corresponds to carbon in silicon and case (b) to phosphorus in silicon. These were the two cases considered by Series, Hurle and Barraclough [13] for the steady situation. Case (c) is included
~
0.5 S 0.4
CONTOUR ERr O.44E*OO O.50E+OO O.57E+OO
2 3
~
~
S
O.75E~OO
7
O.SIE+O0
1
05
05
10
2~0
2/~
7
\
/
I
K
~
~
\ \
//
0~b
to show the effect of changing the Schmidt nurnber. Contour plots of keff and the amplitude A~, are given in figs. 3 and 4. The first point of note is that the contour plots of keff and A~are indis-
CONIQURT
-°•:O
I 2
~
//4
473
I
I I
I I
I
S 0.4
/////)///1/J
CONTOUR EEl 1 O.I5E*OO 2 O.25E*OO 3 O.34E*OO
2
1
~ S
O.62E+OO
7
0.725+00 1
~
~
0
20~
10
~
Fig. 3. Contour plots of k~8~ in (N, s) space for: (a) k0 = 0.07, Sc = 8.75; (b) k0 = 0.35, Sc = 8.75; (c) k0 = 0.07, Sc = 20.0. In each case A = 0.9 and a 0.03. The N axis is logarithmic with the exception that the data at N = 0.25 are replaced by those for N = 0.0.
R~A.Cartwright, A.A. Wheeler / Axial magnetic field andperiodic growth ratefluctuations
474
Table S A comparison between the results of Cartwright et aL [6] and numerical results with A =0
8=10 0.15 k0 = 0.5
o
1 = 0.15 k0 = 0.5 =
8=1 =
Cartwright etal.[6]
—0.19 0 0.4
0.995 0.935 0.891
0.995 0.94 0.90
0.8 0.95
0.934 0.980
0.94 0.98
0.19 0 0.4 0.8 0.95
0.782 0.411 0.337 0.433 0.671
0.78 0.43 0.35 0.45 0.70
0.030 0.249 0.882
0.03 0.25 0.89
0.01 0.1 1
0.5 0.4
=
s
Numerical results
The contours in (N, s) space for c(0, by the equation (1—s)(1+5s) =D, N
t)
are given
aconstantasN~cc. (16)
Hence for large N the dependence of both keff and A0 on N and s will be of this form. Thus the asymptotic analysis supports the numerical evidence that the contour plots are indistinguishable. Eq. (16) also represents the shape of the contours well: the line s = 1 is always a contour and for a large fixed value of N the contours will have a turning point at s = 0.4. We shall now consider in greater detail the asymptotic expressions for keff and A0. From appendix C we have an asymptotic expression for A given by (C.10): (17)
‘
tinguishable. In all three cases an increase in keti is accompanied by an increase in A 0 and vice versa, The second observation is that when S = 1, H —a (which is independent of N). It follows that keff and A0 are, therefore, also independent of N when s = 1. This causes a rapid change in the solution as s —~1 when N is small but the variation more keff gradual that in allis cases and for A large N. We found 0 are minimum for N = 0, s and 0.05 as ands —~ that increase N increases 1. they When N is both large,asthe
where y = 1 + M/k0. We see that aA0/as = 0 when ay/as = 0, that is, at s = 0.4. A0 has a maximum for a fixed value of N, at this value of s, if dA0/dy > 0 and a minimum2)1/2 if dA5/dy <0. dA0/dy = 0 when y = 1 and dA +2 )1/2, (1 — ASubstituting for0/dy <0s = for0.4,y>we1 have + (1 — y with A
minimum of keff or A0 for a fixed N is at s 0.4. The asymptotic analysis given in appendix C supports these findings. From (C.9) we have
With a = 0.03, k0 = 0.07, A = 0.9, N 32 the condition is satisfied and so dA0/dy <0. This
(1 k0)(1 +A sin t) + O(Sc’) 1 + k0(1 + A sin t) + M 2 cc, as Sc—’ cc,N cc, with N/Sc (15) c(0,
—
t)
—
-*
—
where M= —s)(1 c(y, (1 t)=ë(y,
3a, +T),5s)/3A
A=[a+(a2~l~4N)1~~’2]/2.
as
2—A2) 2(1—k0)(y—1)A +O(Sc’) k0(y Sc —p cc, N —~cc, N/Sc2 — cc,
2~___~l
suggests that A0 is minimum when s = 0.4 for a results.value fixed It is of interesting N. This to agrees note with that for thesufficiently numerical large N (N> 75 with the above data), dA 0/ d’y I ~=o.~> 0 and A0 will then have a maximum for s = 0.4. We can also consider the asymptotic forni’of kCff as given by (Cli). k~
2 88=k0+(1_k0)[2_~+
(y—1)
(~2 _A2)1/2j
1
(18)
R.A. Cartwright, A.A. Wheeler / Axial magnetic field and pen odic growth ratefluctuations
475
0.9 0.9 0.7 0.6
1I
S0.5
II II I I
1 I O.27E*OO KE~ IL_CONTOUR 2~ 0.~SE*OOI I 51 3 4 O.83E*OO I 0.105*01 0.54E*O0
0
I
I
0.3
I I
sI 7 I S I
I
0.125*01 0.14E*O1
I
O.ISE*01 O.18E+01 0.19E--O1
9 10
I
2
I
0.2
/
0.1
1I
0
a
0.5
.0
2.0
4.0
8.0
16.
32.
N
_____
S
I 04 0 35
I
O.12E-l-OO
2JI 3
O.14E*OO!\ 0.ISE*00 I
~
0.17E*00
6 O.2OE*0O 5 0.196-1-00 7CONTOUR 0.225-1-00 8 0.245*00 KEY
I
0.2
0.1 0
,) ))) )
I
1
I
I i I I Ii II
3
0 27E 00 O.25E*O0
b 0.5
2 1.0
I
I
5
2I 3 I
I
0.4
I
7
I
/;~~2~
2.0
N
4.0
8.0
16.
32.
0.2 0.1 0 Od~Q
9 10
I
0.5
I
I
O.*3E*OO
O.7OE*0OI O.97E*OO 0.125*01 I
I
II
II
II
0.23E*01
I
I
0.25E*O1 0.295*01
I
SI 0.15E*O1 sII 0.18E*01 7CONTOUR 0.20E*01 KEY
0.3
2
*
___ 9
5
0.5
I
I
/1
I
/ /1
1.0
2:0
N
1.0
2
/
8.0
32) 16.
32.
Fig. 4. Contour plots of A0 in (N, s) space for: (a) k0 = 0.07, Sc = 8.75; (b) k0 = 0.35, Sc = 8.75; (c) k0 = 0.07, Sc = 20.0. In each case A = 0.9 and a = 0.03. The N axis is logarithmic with the exception that the data at N = 0.25 are replaced by those for N = 0.0.
It is easy to show from (18) that k0
<
1 and
N—.oo
This is in agreement with the behaviour of the
numerical results. It should be noted that the changes in keff and A0 with N and s are due to changes in the axial flow. A strong axial flow toward the disc decreases k0ff and A0, whereas with a weak axial flow they increase.
R.A. Cartwright, A.A. Wheeler / Axial magneticfield and periodic growth ratefluctuations
476
It is difficult to substantiate the above analysis with experimental results due to the problems of isolating specific factors in experiments. In practice the effect of the magnetic field on striations is mixed. Series [14] reports that for a crystal grown in a 2000 G axial field the concentration fluctuations were 2% at the centre and 4% at the edge of the crystal. This compared with a control crystal grown in the absence of a magnetic field with corresponding fluctuations of 4% and 3%. The effect of the magnetic field is, thus, inconclusive, In reality, of course, the magnetic field also interacts with many other parts of the system. It has an effect on temperature fluctuations in the melt and hence may modify a ( t), the fluctuating growth rate which in our model is independent of N. We have assumed a constant pull rate, a, but this too may oscillate. Either of these modifications, once quantified could be incorporated into our model to provide a more sophisticated analysis. This paper provides insight into the effect of magnetic fields in fluctuating growth conditions and gives a basis for future developments.
4. Conclusions We have investigated the effect of a magnetic field and crucible rotation on the fluid flow and mass transfer in the vicinity of the crystal/melt interface. The axial flow weakens as N cc and the velocity boundary layer,For which is oflarge thickness 1”2) becomes thinner. a fixed value O(N of N the axial flow is maximum when s = 0.4. When s = 1 we have solid body rotation. Both the numerical and asymptotic results suggest that as N cc, keff i. For a fixed large value of N, keff is minimum at s = 0.4 and keff ~ 1 as s 1. The numerical results show that when 8 N 32, A 0 also has a minimum at s = 0.4 and increases as s 1, for a fixed value of N. An asymptotic analysis supports this and in addition indicates that for a large enough value of N (N> 75 for typical growth data), A0 will have a maximum at .s = 0.4 and will decrease when kCff increases. The phase of the concentration oscillations appears to be independent of N and s. —~
—‘
—‘
—~
~
AcknowLedgement One of us (R.A.C.) gratefully acknowledges the receipt of an SERC grant during the the period when this research was conducted.
Appendix A Due to the coupling between the magnetic and velocity fields, the fluctuations in the velocity field induce time dependence of the magnetic fields. However, in this appendix we show that, in the limit Prm 0, such time dependence is confined to terms of O(Prm) and hence negligible since Prm iO~ in realistic situations. Therefore, it is possible to replace the magnetic field by the steady applied field. The electromagnetic field is governed by Maxwell’s equations which are —‘
—
e0V’ E = Pc’
.
(A.i)
V’ B = 0, ~7xE= _8B/at*,
(A.2)
VXB=~t0(J+c0 8E/at*).
(A.4)
(A.3)
Ohm’s law gives J=
~
(E + ~ x
B),
(A.5)
and we wish to find the Lorentz force F= i.t>< B p
~A 6
‘
where E is the electric field, B is the magnetic field, J is the current density vector, p.0 is magnetic permeability, p is density, Pc is free charge density and e~ is permittivity of vacuum. Eqs. (A.3) and (A.5) give that ~ ~<
,
= ~
( \
+ ~
—
at
*
>
(v
x B )‘~ 1.
(A 7)
R.A. Cartwright, A.A. Wheeler / Axial magneticfield and periodic growth ratefluctuations
J from (A.4) to give, aE\
Substitute for fi
477
These equations can be non-dimensionalised as follows: -
1/2
~
Po
aB
(A.8)
j, \1/2
We also have the equation of continuity for the fluid, V. v = 0, (A.9) and substituting (A.2), (A.3) and (A.9) into the 0-component of (A.8) gives, in cylindrical polar coordinates with axial symmetry,
2
+
az
ar
-~—
—p-ar- (rBe))
r
a B9 r
r8r( —Tv—
—
)J +
a
0
0.
Similarly of (A.8) together with (A.2) and the (A.3)z-component gives, i
a
-,-
-~-,
2B 2\ ~ 2 ~~0a a ~ r-~-,-)+ -j—~-+ -~- (rv2 Br I
~B
—,-—
2B
aB 1at” 0a—
a
+
+ p.
=
—
rv~B~)
0.
(A.i1)
It should be noted that these equations are the same as (14) and (15) obtained by Cartwnght et al. [6] together with additional time dependent terms. The Lorentz force, F, is given by 1
~B
aB
F—— p.p ( ~ z
az
az
aE
—
=
~(-~---
0
strength, 02 is the angular velocity of the disc and p is kinematic viscosity. If we rewrite (A.10) in non-dimensional form it is seen that
a = 1.33 x
B B9
+ —fl— r +
~az — Braz aB,. —B9— — s~0p~
i0_6/~.
Bz
~B9
az
‘
m. 2/N m2,
o = 8.854 x 10—2 C and so e 5. 0W/a=1.398Xi0 So, as in the steady case we obtain L~ O(Prm) and to this order the additional time dependent terms provide no extra contribution. A similar argument applies to the time dependent terms in (A.11). —
~
The time dependent term in the body force is, >< B
at
*
/
From (A.3), I E I = O(B0w1) scale. So the magnitude of thiswhere term /is is a length B0c.,l
O~O—T~-—WBO) =
( at~ 8E xB aE
B~
For silicon, w=1.O5~2.1rad/s
I XB)),
1 I 8B9 — p.0p I\ Br — ar C0
aB ~
=
B0 is the externally applied magnetic field
~B
e0p.0~0a~
(~t)
=
~
where
L
—
L~+Prm(NGf—HLf)+Pr~4——Lii+Li)=0. I (002
—)
2.B~
(A.10)
~
~‘
8B
7/3Z v +p.0aIB~_ 2—+Br az9 rar
F
/ j~\1/2 r ~wi
o(NW2l(~~~)).
The steady terms in the body force are of magnitude, NW/alp.0,
+
3B Brar —i
—~(~x~)J.
and so the ratio of the time dependent terms to the steady terms is (~w/a) 12 (W/P) Prm, whichisnegligible.
478
R.A. Cartwright,
A.A.
Wheeler / Axial magnetic field and pen odic growth rate fluctuations
Appendix B
Expand g and h as power series in X1
In this appendix an asymptotic solution valid for large magnetic fields is developed. We have from (4)
g(x,
00
H~,H,,,,~+H,,,,(I’(t)_H)+~H~ —2(G2 — s2)
—
NH,,,
(B.1)
G~= G~,,+ GH~+ G~(I’(t)— H)
—
N(G
— s),
(B.2) with boundary conditions H(y,
t)
=
—a,
H,,(y, at
t)
0,
G(y,
t)
=
t)
t) —
y—’cc.
~, (B.3)
Let X=
(B.4)
[~+(~2+4N)h/’2]/2,
and put
t)
—à+h(x, t)/X3,
G(y, t)=g(x,
0
=
h0
=
1
=
(B.9) t).
(1 — s) ex + s, 2(1—s) ~ (1 + 2s + 3sx) e_x (1 —s)(1 + 5s) — 3
—
(1
2x
2
s) e (B.10)
=
—s(i
—
(B.ii)
s) I’(t),
—(i—s) ~(i+5s) 3
I’(t)—3sä],
—(1—s) r( 7+115) I’(t)—d(2+7s)], 9 (1 — s) I’(t)
$ 5_i~(1—s)[7I(t)—2a].
~
Eqs. (B.10) and (B.il) can be substituted into
g~~—(g—s) 1[g~(I’(t) = —A
(B.6)
(B.9) to give expressions for g(x, t) and h(x, t) and then using (B.5) we obtain G(y, t) and H(y, t).
(B.7)
AppendixC
+a(g—s)]
+X2g~—X4(g~h—gh~),
with boundary conditions t)
x’~hr(X,
2 I’(t) 2+f3 x e_x h1=/3o+($1x 2x+$3)e~ 2x, +(/34x+ 135) e_ where ~o= ~
/3 =
~
h (x,
~
=
~=
Eqs. (B.1)—(B.3) become 2—s2) h~~~—h~—2(g
=
—(1
(B.5)
t).
t)
g
~
=
t),
and substitute in (B.6) and (B.7) to obtain
x=Xy, H(y,
~ X~g,.(x, 00
h(x,
g 0, G(y,
=
r=0
1,
y=O
H,,(y, as
=
t)
=
0,
h~(x,
t)
=
In this appendix we conduct an asymptotic 0,
at x = 0, h~(x, t) —p0, as x
—b
cc.
g(
~,t)
g(x,
t)
=
1,
analysis of the concentration which supports the results produced numerically. We must solve eqs. (10)—(12),
—‘s, (B.8)
C~=
+
[I’(t)
—
H(y,
t)]
c,,,
(C.i)
R.A. Cartwnight, A.A. Wheeler
/ Axial magneticfield and periodic growth ratefluctuations
with boundary conditions,
with
(1—k0) Sc a(t) c(0, 1)+c,,(0, t)=0,
(C.2)
c(y, t)
(C.3)
—
1, as y
—
cc,
where c(y, t)
=
e(y, r).
t)
3(i
— —a + ~A —(1
—
—
s) e2AJ
s)[2(1 —
+
(1 + 5s)]
+
—
_~3(a +
3 + 5s) ) + O(A~), (1 —s)(i cc, with A/Sc — cc.
—
c”(3, t)=c
2c~’(.P, t)+O(Sc2), 0”(j, t)+Sc and find from eqs. (C.6)—(C.8) that at leading order ~
t)
1
=
+
(1
—
k
0)ä(1 + A sin t)
k0d(i +A sin t) + K
reveals that H(y, t)
sin t+â+K)j5], xexp[—A where K= (1 — s)(1 + 5s)/3. This gives 3(aA
—A3{d— -~(1— s)[2(1 + 2s+ 3sAy) e_AY 2XY_(1
—(1 —s) e asA—’cc,
3~.
+
5s)]}
where a =table A 3 we see that Sc From
+O(A4),
— t)
c”(O, t) —1
(C.4)
— 8.75
—‘
cc,
cc. We expand c(y, t) in powers of Sc’ 2c 3). (C.5) = c0 + Sc~c1+ Sc 2 + O(Sc
Then eq. (C.1) and boundary condition (C.3) give c0(y, t) = 1, c.(y, t) = 0 for i 1, which does not satisfyHence the boundary condition at the interface. (C.5) represents the(C.2) concentration in an outer region where y is 0(1) and we require an inner region where y = O(Sc~).We put j3 = Sc ~ in which case eqs. (C.1) and (C.2) give the following problem for c” which represents the concentration in this region
+
(1— k k 0)a(1 +A s~nt) 0a(1 +A sin t) +K
+O(Sc’). Hence
and so it is
appropriate to consider the realistic limits Sc
c(y,
Eq.
We put where A = [~+ (a2 + 4N)1~2]/2. For a typical crystal growth system a 0.03 and so the first two terms of the expansion are of comparable magnitude. If we assume a — A ~ then a similar analysis to that given in appendix B
A/Sc
t).
H*(j~,r) as A — cc, Sc
O(A4)
(C.8)
In this inner region we denote H as H”(j3, (C.4) gives
2s + 3sAy) e~
as A—’cc,
—
(1_k0) a(t) c*(j~,t)+c~’(y, t)=0,at ~3=0, (C.7) along with the matching condition c”(33, t)—i1,asy---cc.
In appendix B we derived an asymptotic form for H(y, t), H(y,
479
(C.9)
(1 —k 0)d(1 +A)
c, 0~= 1 + Cmth =
1
+
k0â(i —A) + K (1 — k k 0)d(1 — A)
0d(1 —A)+K
3), + O(Sc 1), +
O(Sc
which gives 1), (C.10) A= 2(1 k— k0)(y — i)A + O(Sc C 2 — A2) 0(’y where ‘y = 1 + K/k 0 a. Further, eqs. (C.9) and (14) give, [2 + (y — 1)2 kCff = k0 + (1 — k0) —
___________ (72_A2)1/2]’
~-c7 =c;~+(J’(t) -H(j~, t))c;,
(C.6)
(C.11)
480
RA. Cartwright, A.A. Wheeler / Axial magneticfield and periodic growth ratefluctuations
References [1] K.S. Choe and T.H. Strudwick, J. Crystal Growth 74 (1986) 338. [2] H.A. Chedzey and D.T.J. Hurle, Nature 210 (1966) 933. [3] L.O. Wilson, J. Crystal Growth 48 (1980) 435. [4] L.O. Wilson, J. Crystal Growth 56 (1982) 219. [5] A.A. Wheeler, Proc. Roy. Soc. (London) A379 (1982) 305. [6] R.A. Cartwright, D.T.J. Hurle, R.W. Series and J. Szekely, J. Crystal Growth 82 (1987) 327. [7] L. Hjellming and J. Walker, J. Fluid Mccli. 164 (1986) 237. [8] L. Hjellming and J. Walker, J. Crystal Growth 83 (1987) 51.
[9] L. Hjellming and J. Walker, J. Fluid Mech. 182 (1987) 335. [10] R.A. Cartwright, N. El-Kaddah and J. Szekely, IMA J. AppI. Math. 35 (1985) 175. [11] J.F. Brady and L. Durlofsky, J. Fluid Mech. 175 (1987) 363. [12] T. von Karman, Z. Angew. Math. Mech. 1 (1921) 233. [13] R.W. Series, D.T.J. Hurle and K.G. Barraclough, ImA J. Appl. Math. 35 (1985) 195. [14] R.W. Series, privative communication, 1987.