Journal of Crystal Growth ∎ (∎∎∎∎) ∎∎∎–∎∎∎
Contents lists available at ScienceDirect
Journal of Crystal Growth journal homepage: www.elsevier.com/locate/jcrysgro
Numerical study of silicon crystal ridge growth G. Barinovs n, A. Sabanskis, A. Muiznieks † University of Latvia, Faculty of Physics and Mathematics, Zellu St 8, LV-1002 Riga, Latvia
art ic l e i nf o
Keywords: A1. Growth models A2. Czochralski method A2. Growth from melt B1. Semiconducting silicon
a b s t r a c t The size of the ridge-like protrusions appearing on the external surface of dislocation-free 〈100〉 silicon crystals grown from a melt was studied theoretically. According to existing models the growth of the ridges is caused by the presence of f111g crystal planes at the crystal–melt interface. They affect the height of triple phase line, free surface orientation and the crystal growth angle. A numerical 2-dimensional model was proposed for the calculation of the size of the crystal ridges. The model included the effect of the undercooling of the crystal–melt interface on the crystal growth angle. The numerical model estimated the effect of the ridge size on the free surface at the triple phase line. The size of the crystal ridges was calculated for the 〈100〉 silicon crystals grown by the Czochralski process. It was shown that the crystal ridge growth is sensitive to the physical conditions at the triple phase line. & 2013 Elsevier B.V. All rights reserved.
1. Introduction Nowadays industrial monocrystalline silicon is mainly produced by the Czochralski process (CZ) pulling a crystal from a melt in a crucible. Higher purity crystals are grown from a molten zone using Floating Zone (FZ) technique. The crystals grow cylindrical due to the crystal rotation with respect to the melt during the crystal growth. However, the cylindrical symmetry of crystals is broken in the direction where f111g planes “just touch” the cross-section of a crystal [1]. The facets, formed by f111g crystal planes at the crystal– melt interface (internal crystal surface), grow by different growth mechanisms, 2D nucleation, than the rest of the crystal, growing by rough growth. Depending on a crystal growth process and growth conditions, the crystal external surface in front of the f111g facets either bulges outwards forming ridges on the surface or it bulges inwards forming crystal constrictions. These deviations from the cylindrical symmetry are observed in dislocation-free crystals and they disappear if dislocations appear during crystal growth [1]. Usually, four bulges are observed for 〈100〉 crystal growth and three or six bulges are observed for 〈111〉 growth. For the crystal–melt–vapour triple phase line (TPL) to stay stationary during crystal pulling from a melt, all parts of the internal surface should grow with the crystal pulling speed. To achieve that, the f111g planes growing by 2D nucleation need larger crystal–melt interface undercooling than a surface growing by rough growth. Therefore the internal f111g planes lag in cooler regions of the system behind the rest of crystallization front
n
Corresponding author. Tel.: +371 67033769. E-mail address:
[email protected] (G. Barinovs). † A. Muiznieks died in April 2013, we inserted death dagger in the manuscript after his name, I believe that a note acknowledging his death would by appropriate.
growing by rough growth. Since for the CZ growth the crystal is pulled upwards, the TPL rises higher at places where it borders internal f111g planes. The TPL descends lower for the FZ process, where the crystal is pulled downwards. Different heights of the TPL change the melt orientation at the TPL and due to hydrostatic effects the free surface for FZ crystal bends outwards. The free surface bends inwards at crystal ridges for the CZ process. In addition to the free surface perturbation, the theory of Voronkov [2] introduces the deviation of the apparent growth angle εn from intrinsic growth angle ε at TPL by an angle χ
εn ¼ ε χ ðΔTÞ:
ð1Þ
This deviation is explained by the curvature driven surface diffusion of silicon atoms and its magnitude is mainly determined by the undercooling of crystal–melt interface at the TPL, ΔT. The angle χ can be calculated as
χ ðΔTÞ ¼ K ΔT=v1=3 p ;
ð2Þ
where vp is the crystal pulling speed and K is a constant that according to Voronkov could be found using experimental observation that χ ¼ 201 when vp ¼3 mm/min and ΔT ¼ 3:7 K [2]. The effect of a nonzero angle χ is a tendency for the crystal external surface to bulge outwards. The effect of the angle χ is significant when internal f111g facets are present since the f111g internal facets require considerable undercooling for growth at sufficient speed. Thus this mechanism competes with the bending of the melt surface inwards due to hydrostatic effects for the CZ process and may lead to the formation of crystal ridges. For the FZ process the χ angle amplifies the effect of the melt bulging outwards thus predicting the growth of higher ridges on the lateral surface of a silicon crystal. The ridge shape and size is controlled by the physical conditions (e.g. temperature gradients,
0022-0248/$ - see front matter & 2013 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.jcrysgro.2013.12.019
Please cite this article as: G. Barinovs, et al., Journal of Crystal Growth (2013), http://dx.doi.org/10.1016/j.jcrysgro.2013.12.019i
G. Barinovs et al. / Journal of Crystal Growth ∎ (∎∎∎∎) ∎∎∎–∎∎∎
2
crystallinity, free melt surface oscillations) at the TPL during the crystal growth.
In linear approximation the deviation of the melt angle from its value at δh ¼ δr ¼ 0, φm ðr 0 ; h0 Þ ¼ ε could be written as
Δφm ðδr; δhÞ ¼ φm ðr; hÞ φm ðr0 ; h0 Þ ¼
2. The numerical model We present a numerical 2D model to calculate the size of the crystal ridges for CZ growth. The position of the TPL at the crystal ridge with respect to the rest of crystallization front was determined by taking into account the perturbation to free melt surface by the crystal ridge and by the displacement of the TPL. The ridge growth was described by the system of first order differential equations with stationary solutions describing vertical crystal ridge growth. The effect of different physical conditions on the crystal ridge size was calculated. ! The velocity v of a point on the TPL line in laboratory frame of ! reference is a sum of crystal growth velocity, v TPL , and crystal ! pulling velocity, vp : ! ! ! v ¼ v TPL þ vp : ð3Þ We assumed that the point will stay in a vertical plane containing crystal axis and we introduced vertical coordinate h and radial coordinate r, see Fig. 1, where r is the crystal radius and h is the meniscus height. At some meniscus height, h0, the crystal will grow vertically with radius r0. At a crystal ridge the meniscus height and the radius will deviate from these values by δr ¼ r r 0 and δh ¼ h h0 , Fig. 2. The δr; δh values determine the position of the TPL line at the crystal ridge with respect to the position of the TPL for the rest of the line. We denoted the angular deviation of the growth direction from the vertical direction by αðr; h; ΔTÞ. Assuming that the crystal is ! pulled up, we could find vh and vr projections of v as ( vh ¼ dh=dt ¼ vp cos ½αðr; h; ΔTÞvTPL ðΔTÞ ð4Þ vr ¼ dr=dt ¼ sin ½αðr; h; ΔTÞvTPL ðΔTÞ:
ð7Þ
The derivative of the melt angle over p h ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi for small δh is estimated in Ref. [6] for CZ growth as ∂φm =∂h 2g ρm =γ m , where ρm is the density and γm is the surface tension of the molten silicon. The derivative of the melt angle over r depends on the height of the ridge as well as on its width. Let us assume that at a ridge the crystal–melt interface is formed entirely by f111g internal facet as shown in Fig. 2. This will be approximately true for 2D nucleation growth, since the curved part of the interface will be limited due to the considerable undercooling and the Gibbs–Thomson effect. The maximal width of such facet is qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi wmax ¼ 2 r 20 ðr 0 r f Þ2 ¼ 2 r 20 ðr 0 δh cotðφ111 Þ þ δrÞ2 ; ð8Þ where r0 is the crystal radius. For the estimate of the perturbation of the meniscus angle near the crystal ridge, we started with the perturbation treatment suggested in Ref. [7]. The meniscus surface can be described by the Young–Laplace equation stating that the fluid pressure is balanced by the surface curvature pressure
γ ð1=R1 þ 1=R2 Þ ¼ ρgh
ð9Þ
r [100]
The velocity of 2D nucleation in ½111 direction is given in Ref. [3]. Assuming that the rough growth is much faster than the 2D ! nucleation growth, the projection of the growth velocity, v TPL , in ½111 direction is equal to the speed of 2D nucleation vTPL ¼ 0:226e 140=ð3ΔTÞ ðΔTÞ2=3 = cos ðφ111 αÞ;
∂ φm ∂φ δr þ m δh: ∂r ∂h
crystal
ð5Þ
where φ111 is the angle between pffiffiffi h-axis and ½111 direction. For 〈100〉 growth φ111 ¼ arccosð1= 3Þ 54:74○ . The value of the angle α depends on the angle of the free surface at the TPL, φm, and, according to Voronkov, also on undercooling ΔT. Using the definition of the growth angle and Eq. (1), we can write
αðr; h; ΔTÞ ¼ φm ðr; hÞ εn ðΔTÞ ¼ φm ðr; hÞ ε þ χ ðΔTÞ;
h {111}
r0 r0 -rf
ð6Þ
where angles φm and α are measured counterclockwise from the h-axis. The growth angle ε for silicon is equal to 111 [4,5] and it was assumed to be independent of the crystallographic direction.
Fig. 2. The width of a crystal ridge is found from cross-section of internal f111g facet and the cylinder cross-section.
h vp
vp
r0
Crystal Melt
h0
vTPL
Crystal
r
Melt
m
! Fig. 1. a) The r; h coordinate system and main physical parameters for CZ growth. Dashed lines show the interfaces. v p is crystal pulling velocity, (b) zoom of the region ! ! around the triple phase line in front of internal facet. v TPL is velocity of the TPL. The angle α is the angle between the vertical direction and v TPL . The melt angle is φm.
Please cite this article as: G. Barinovs, et al., Journal of Crystal Growth (2013), http://dx.doi.org/10.1016/j.jcrysgro.2013.12.019i
G. Barinovs et al. / Journal of Crystal Growth ∎ (∎∎∎∎) ∎∎∎–∎∎∎
3
w TPL wmax 2
w
r0
rmax
r
r
wmax 2
Fig. 3. (a) Coordinate system chosen to describe the geometry of a ridge on a crystal surface. (b) Magnified view of the ridge showing the ridge height δr max and the ridge width wmax .
δ h/wmax
where R1 and R2 are the principal radii of the meniscus surface curvature. If the meniscus is oriented close to the vertical direction and the crystal ridge is small comparing to the capillary length, Eq. (9) can be written as [7] ð10Þ
The perturbation of the melt meniscus by a ridge in h ¼ h0 plane, Fig. 3, was analyzed. The meniscus was described by the function rðh; wÞ ¼ r 0 ðh; wÞ þ δrðh; wÞ;
ð11Þ
where r 0 ðh; wÞ describes the shape of the meniscus in the absence of the growth ridge and δrðh; wÞ describes meniscus perturbation by the ridge. Inserting Eq. (11) into the Young–Laplace equation (10), perturbation of meniscus satisfies the equation [7] ∂2 δrðh; wÞ 2
∂h
þ
∂2 δrðh; wÞ ¼ 0: ∂w2
ð12Þ
For boundary condition we requested that the perturbation vanishes at the ridge boundaries
δrðh; w ¼ 7 wmax =2Þ ¼ 0
ð13Þ
and the meniscus slope perturbation also vanishes at the ridge boundaries
Δφm ðh; w ¼ 7 wmax =2Þ
∂δrðh; w ¼ 7 wmax =2Þ ¼ 0: ∂h
ð14Þ
In Ref. [7] Voronkov presents approximate numerical solutions to Eq. (12) with the boundary conditions above. However, from the analysis of Laplace equation (12), the following solution satisfying the boundary conditions, Eqs. (13) and (14), can be found [8]
δrðh; wÞ ¼ r max enπδh=wmax cos ðnπ w=wmax Þ;
ð15Þ
where n is the positive odd integer. The validity of the solution for δh r 0 can be confirmed by substitution in Eq. (12) and boundary conditions. The solution for n ¼1 is plotted in Fig. 4. The perturbation to the meniscus diminishes moving down from the TPL line. The solutions corresponding to n 4 1 might be used to explain observation of additional crystal ridges forming besides the main ridge [9]. According to the solution, Eq. (15), the largest perturbation of the melt angle due to the ridge is in the middle of a ridge, w ¼0,
0.2
w/wmax
2
∂h
∂ rðh; wÞ ρgh þ ¼ : γ ∂w2 2
-0.
.4 -0 -00.6.8 -
∂ rðh; wÞ 2
0.4
2
0 -0.2 -0.4 0
0.2
0.4 0.6 δr/δ rmax
0.8
1
0 -0.1 -0.2 -0.3 -0.4 -0.5 -0.6 -0.7 -0.8 -0.9 -1
Fig. 4. The calculated meniscus shape perturbation by a crystal ridge for different height, δh, of the meniscus cross-section.
and is equal to
Δφm ðh ¼ h0 ; w ¼ 0Þ nπ
δrmax wmax
:
ð16Þ
For CZ growth this decrease in the melt angle will act oppositely to the effect of Voronkov0 s angle χ. We can observe that if a ridge growth outwards, then according to Eq. (8) the width of the ridge, which is equal to wmax , becomes smaller. If the ratio δr max =wmax increases, then according to Eq. (16) the meniscus bends inwards at the TPL and the outward growth stops. Thus, the narrowing of the crystal ridge by growing outwards restricts the height of the crystal ridge.
3. Calculations of the crystal ridge size in model systems The model presented above allowed us to calculate the size of the crystal ridges for a particular crystal growth experiment. The set of ordinary differential equations (Eq. (4)) with a guess for the initial conditions, rðt ¼ 0Þ, hðt ¼ 0Þ, was solved by using the GNU Octave software package [10]. The input parameters of the model were the crystal radius, pulling velocity and the temperature at the TPL. There have been several numerical models presented in the literature that allow the calculation of the temperature gradients in the crystal–melt system without analyzing the crystal ridge effect. The typical values of the are several thousands Kelvins per meter. On the contrary, we do not know any
Please cite this article as: G. Barinovs, et al., Journal of Crystal Growth (2013), http://dx.doi.org/10.1016/j.jcrysgro.2013.12.019i
G. Barinovs et al. / Journal of Crystal Growth ∎ (∎∎∎∎) ∎∎∎–∎∎∎
4 10
vertically. For different initial conditions the calculations converge to the same stationary values of r and h. Results of calculations of the effect of different temperature gradients on the ridge height and width are summarized in Table 1. According to the calculations, smaller temperature gradients lead to higher rise of the TPL in the front of an internal f111g facet and therefore to wider crystal ridge. Decreasing the temperature gradient even further, the bending of the free surface inwards is not able to compensate the effect of the χ angle and a wide constriction will form on the crystal surface.
H(mm)
8 6 4 2
0.355
0.356
0.357
0.358
0.359
0.36
0.361
r(mm)
4. Conclusions
Fig. 5. Calculated height of a CZ crystal ridge for 5000 K/m weighted average temperature gradient.
10
H(mm)
8 6 4 2 0
3.05
3.1
3.15
3.2
3.25
3.3
3.35
3.4
w(mm) Fig. 6. Calculated width of a CZ crystal ridge for 5000 K/m weighted average temperature gradient.
Table 1 Calculated 〈100〉 CZ silicon crystal ridge height and width
Based on perturbation treatment of the meniscus in the presence of a crystal ridge suggested by Voronkov, we presented an analytical solution to the Laplace equation describing the perturbation of the meniscus. A geometric model was suggested to relate the width and the height of a crystal ridge that allowed us to implement a numerical algorithm to calculate the size of the crystal ridges from known temperature gradients in the crystal–melt system. The calculations for 〈100〉 silicon growth show the decrease of the height and the width of the crystal ridges, increasing the temperature gradients in the crystal–melt system. We observed that the crystal ridge size is a sensitive probe of the temperature gradients at the crystallization front. Therefore the visually obtained information about the crystal ridge size during crystal growth process can be used to estimate the thermal gradients, which determines the point defect distribution [11] and can lead to the generation of dislocations in the single crystal. This might lead to a better control of the experimental crystal growth processes.
Acknowledgments
Temperature gradients (K/m)
δr (mm)
w (mm)
3000 5000 10000
0.55 0.36 0.18
13 3.3 1.4
The present work is carried out at the University of Latvia and has been supported by the European Regional Development Fund, project contract No. 2011/0002/2DP/2.1.1.1.0/10/APIA/VIAA/085. References
example of direct calculation of the temperature at the crystal– melt interface at the crystal ridge. For our calculations we will assume that the dependence of the crystal–melt interface undercooling on the r; h coordinates is linear:
ΔTðr; hÞ ¼
∂T ∂T δr þ δh: ∂r ∂h
ð17Þ
The nucleation will only take place close to the coolest place of an internal f111g facet at the outermost part of the ridge. A point on the external crystal surface with coordinates r; h, due to the crystal pulling, will rise to a height H ¼ h þ vp t. We can imagine that the external surface will be bounded by the TPL lines that have formed at different times. Thus the H as a function of r will show the final external surface. The calculated width and height of a crystal ridge is shown in Figs. 5 and 6 correspondingly. The calculation was done assuming temperature gradient components in r and h directions equal to 5000 K/m. The crystal radius was 4 in. and the crystal pulling velocity was assumed to be equal to 2 mm/min. The values of the height and width of the ridges start at some initial values and oscillating converge to the final values. After that, if the temperature gradients remain unchanged, the ridge grows vertically if the rest of the crystal also grows
[1] T. Ciszek, Non-cylindrical growth habit of float zoned dislocation-free [111] silicon crystals, J. Cryst. Growth 10 (1971) 263–268. [2] V.V. Voronkov, Mass transfer at the surface of a crystal near to its boundary with the melt, Sov. Phys., Crystallogr. 23 (1978) 137–141. [3] K.A. Jackson, Response to: Some remarks on the undercooling of the Si(111) facet and the “Monte Carlo modeling of silicon crystal growth” by Kirk M. Beatty (Kenneth A. Jackson, J. Crystal Growth 211 (2000) 13 by W. Miller, J. Cryst. Growth 325 (2011) 104). [4] T. Surek, B. Chalmers, The direction of growth of the surface of a crystal in contact with its melt, J. Cryst. Growth 29 (1975) 1–11. [5] T. Duffar, Crystal Growth Processes Based on Capillarity: Czochralski, Floating Zone, John Wiley & Sons, 2010. [6] V. Voronkov, B. Dai, M. Kulkarni, 303—Fundamentals and engineering of the Czochralski growth of semiconductor silicon crystals, in: Pallab Bhattacharya, R. Fornari, H. Kamimura (Eds.), Comprehensive Semiconductor Science and Technology, Elsevier, Amsterdam, 2011, pp. 81–169. [7] V.V. Voronkov, The effect of the faceting of the crystallization front on the external shape of crystals, Akad. Nauk SSSR Izv. Seriia Fiz. 49 (1985) 2467–2472. [8] P. de Gennes, F. Brochard-Wyart, D. Quere, Capillarity and Wetting Phenomena: Drops, Bubbles, Pearls, Waves, Springer, 2004. [9] K. Fritzler, E. Trukhanov, V. Kalinin, P. Smirnov, A. Kolesnikov, A. Vasilenko, In situ monitoring of floating-zone-grown Si(111) crystal structure using the behavior of ridgelike protrusions, Tech. Phys. Lett. 33 (2007) 521–523. [10] J.W. Eaton, D. Bateman, S. Hauberg, GNU Octave Manual Version 3, Network Theory Ltd, 2008. [11] W. von Ammon, E. Dornberger, H. Oelkrug, H. Weidner, The dependence of bulk defects on the axial temperature gradient of silicon crystals during Czochralski growth, J. Cryst. Growth 151 (1995) 273–277.
Please cite this article as: G. Barinovs, et al., Journal of Crystal Growth (2013), http://dx.doi.org/10.1016/j.jcrysgro.2013.12.019i