Semiconductor Bridgman growth inside inertial flight mode orbiting systems of low orbital eccentricity and long orbital period

ARTICLE IN PRESS

Journal of Crystal Growth 273 (2004) 1–18 www.elsevier.com/locate/jcrysgro

Semiconductor Bridgman growth inside inertial flight mode orbiting systems of low orbital eccentricity and long orbital period X. Ruiza,,1, M. Ermakovb a

Applied Physics Laboratory, Universitat Rovira i Virgili, Tarragona, Spain b Institute for Problems in Mechanics, RAS, Moscow, Russia Received 3 June 2004; accepted 20 July 2004 Communicated by K.W. Benz Available online 13 October 2004

Abstract The present work numerically investigates the influence of a generic inertial flight mode orbiting system of low eccentricity and long orbital period on the final solid dopant distribution of a set of semiconductor crystals virtually grown by the Bridgman method. For the lowest 1 mg level it has been noted that, depending on the semiconductor, some disturbances could arise. For higher mg levels spot-like alternate dopant structures with a frequency equal to the orbital one have been computationally predicted in the solid phase. At higher levels, these structures spread trying to attain the opposite side of the crystal breaking any kind of symmetry. All these results are independent on the length of the sample and on the external thermal environment. Thus, to improve crystal homogeneity inside long period inertial flight mode spatial platforms some corrective strategies are strongly recommended. r 2004 Elsevier B.V. All rights reserved. PACS: 47.11.+j; 81.05.Cy; 81.10.Aj; 81.10.Mj Keywords: A1. Directional solidification; A2. Bridgman techniques; B2. Semiconducting materials

1. Introduction

Corresponding author.

E-mail addresses: [email protected] (X. Ruiz), [email protected] (M. Ermakov). 1 Also at Institut d’Estudis Espacials de Catalunya, IEEC, Barcelona, Spain.

A generic orbiting system could move around the Earth using two different basic strategies, the gravity-gradient stabilized and the inertial flight modes [1]. During the second one the orientation of the orbiting system remains fixed relative to some external celestial object, usually the sun. This

0022-0248/$ - see front matter r 2004 Elsevier B.V. All rights reserved. doi:10.1016/j.jcrysgro.2004.07.071

ARTICLE IN PRESS 2

X. Ruiz, M. Ermakov / Journal of Crystal Growth 273 (2004) 1–18

forces the components of the mg-vector to change periodically in a single frequency way. In the case of long orbital periods the frequency of this change is very low as well as the magnitude of the accelerations implied, consequently, for a long time the crystal growth community considered not significant the impact of this kind of environments on growth systems at high temperatures. However, recent reports concerning the growth of InP:S during the European Retrievable Carrier, EURECA-1, mission using the automatic mirror facility, AMF, concluded that, in that case, long crystal growth experiments should be carefully planned because the changes of orientation of the mg-vector relative to the temperature gradient seems to really affect the quality of the resulting solid phase [2–4]. A new question arose at that moment. It is interesting to mention that vibroconvection in doped melts or binary mixtures as well as their impact on the solid segregation has been, since the 1990s, a subject of constant attention. However, practically all works focuses on the consideration of oscillatory gravity signals of relatively high frequencies—of the order of 0.01–100 Hz—acting orthogonally to the thermal gradient. Quasi-steady [5,6], time-dependent [7,8] or time-averaged [9,10] formulations have been used to computationally solve the problem and the principal conclusions about are the unquestioned impact of the thermovibrational buoyancy-driven flows on the transfer of mass, heat and momentum in fluid systems as well as the irrefutable demonstration, by direct experimental checking, of the reliability of numerical models to analyze its g-sensitivity. To deep inside the former open problem we will present here a set of results concerning the impact of quasi-steady accelerations on the solid dopant segregation using relevant semiconductor materials and a 2D Bridgman configuration. Although our preliminary results on short samples and thermal conductive boundary conditions clearly stated a relation between quasi-steady accelerations and dopant distribution in the solid phase [11,12], the present paper goes further completing them in two different fronts, i.e., three different sample lengths and two different thermal environ-

ments, trying to locate the present computational predictions as near as possible to the real ESA facilities.

2. Model formulation The model used for the Bridgman growth system, shown schematically in Fig. 1, is a revised version of an old unsteady 2D model in which all the walls of the ampoule stay at rest except the solidification front which is defined as a flat moving interface [13,14]. Because the weak convection existing in specific mg environments, flat interfaces have been considered as appropriate hypothesis to conveniently approximate the real processes. The velocity of this plane front is supposed to be equal to the growth velocity from the beginning to the end. In this way, the Fu–Wilcox exponential law for initial transients is not tacked into account here [15]. Also, note that the consideration of a growth front moving at a constant velocity unfortunately eliminates also the possibility to analyze any link between the external environment and the variations of the growth rate. The transport of momentum, heat and dopant/ solute has been considered here time dependent, not quasi-steady, because of the continuous decreasing of the melt volume inside the finite length external ampoule. The time-dependent conservation equations of momentum, heat and species in the liquid phase have been made dimensionless using H 2 =n; H, n=H; DT=Th–Tm, c
0 as scales for time, length, velocity, temperature and concentration respectively. Th, Tm and H are the highest, the lowest, melting temperatures and the height of the domain, respectively. After that and because the variation of the domain of integration during the solidification process, a new transformation x ¼ x0 =SðtÞ being SðtÞ ¼ AR2Wt; AR ¼ L=H the aspect ratio of the cavity and W ¼ VH=n the nondimensional growth velocity, enabled to map the old variable x0 2 ½0; AR into a new one x 2 ½0; 1: In this way, the final dimensionless equations in 2D primitive
 variables with respect to a reference O; x; y in the transformed domain could be

ARTICLE IN PRESS X. Ruiz, M. Ermakov / Journal of Crystal Growth 273 (2004) 1–18

3

2D Bridgman model Th Tm + ∆T(t)

Tm One body furnace Perfectly conductive profile

Tl S(t)

ϑ

CRYSTAL

High Temperature Liquid y’

H

δ x’

O’ L

Tl Hyperbolic profile Three bodies furnace

Tm Tm + T(t)

Th

Fig. 1. Details of the 2D Bridgman growth arrangement considered here.

written as 1 qu qu þ ¼ 0; SðtÞ qx qy

(1a)

qu u xW qu qu 1 qp 1 q2 u q2 u þ þv ¼

þ 2 þ SðtÞ qx S ðtÞ qx2 qy2 qt SðtÞ qx qy  1
T Ra y þ Rac c ux ; þ Pr ð1bÞ qv u xW qv qv qp 1 q2 v q2 v þ þv ¼ þ 2 þ qt SðtÞ qx qy qy S ðtÞ qx2 qy2  1
T Ra y þ Rac c uy ; þ Pr ð1cÞ   qy u xW qy qy 1 1 q2 y q2 y þ þv ¼ þ ; qt SðtÞ qx qy Pr SðtÞ2 qx2 qy2 (1d)

  qc u xW qc qc 1 1 q2 c q2 c þ þv ¼ þ ; qt SðtÞ qx qy Sc SðtÞ2 qx2 qy2 (1e)   qy a0 1 q2 y q2 y ¼ þ qt Pr SðtÞ2 qx2 qy2

(1f)

being y ¼ T T m =DTð0Þ the nondimensional temperature, with DT ¼ T h T m and [ux, uy] the Cartesian components of the unitary vector u=g/ g. The definition of the different nondimensional groups, like the thermal and concentrational Rayleigh numbers, RaT and Rac, as well as the Prandtl number, Pr, is shown in Table 1 and is omitted here for the sake of brevity. Note that Eqs. (1a)–(1e) are directly referred to the transport of heat and matter in the liquid phase while that Eq. (1f) governs the heat transfer but in the ampoule walls. The relative diffusivity a0 ¼ aamp =a of this last equation is the ratio between the

ARTICLE IN PRESS X. Ruiz, M. Ermakov / Journal of Crystal Growth 273 (2004) 1–18

4

Table 1 Definition and numerical values of the different nondimensional parameters used Number

Microgravity Prandtl Thermal Rayleigh Solutal Rayleigh Schmidt Peclet Segregation coeff. Relative solid diff. Relative wall diff.

Pu Pr (  103) RaT Rac Sc Pe k Ds* a0 (  102)

Ge

GaSb

CdTe

GaAs

InP

Sn

Ga

Te

Zn

Se

S

Bi

4.227 mg 7.15 0.423 mg — 6.20 0.48 0.087 — 5.8

3.022 mg 44.33 0.007 mg — 76.6 6.4 0.34 — 12.2

23.26 mg 400.58 2.326 mg — 41.50 1.00 1.35 0.005 100.9

2.858 mg 68.06 0.336 mg — 108.45 4.67 0.3 — 14.6

5.796 mg 15.21 0.515 mg

2.244 mg 15.17 0.041 mg — 162.5 10.13 0.29 — 6.1

24.29 7.69 0.47 — 9.9

Bi–Sn

7.676 mg 17.05 0.192 mg

47.27 mg 50 17.05 0.03 — 12.1

(mg H3)/na nXa (mg bT DTH3)/na (mg bc c0 H3)/nD nXD VH/D — Ds/D aampXa

The factor mg should be replaced by the numbers 1, 5, 10, 50 for 1, 5, 10 and 50 mg, respectively.

thermal diffusivity of the ampoule to the thermal diffusivity of the liquid phase. No-slip conditions have been used to model a perfect contact between the internal part of the walls and the liquid (neither detaching effects nor bubbles have been considered here). Perfectly conducting and hyperbolic—yðx; 0; tÞ ¼ yðx; 1; tÞ ¼ tan hðpxÞ=tan hðpÞ; with p=0.8, by convenience— thermal profiles have been applied to the external upper and lower sides of the ampoule (neither radiation nor convection heat exchange have been considered) to mimic the thermal environment generated by one and three zones standard Bridgman furnaces. Across the inner solid walls and the liquid, heat flux continuity has been considered. In addition left and right sides of the domain have been defined as isotherms at different constant low and high temperatures, 0 and 1, respectively. Note that this kind of condition in relation with the growth interface is another approximation of the model—the incomplete Stefan problem. The ampoule material has been selected as a generic nonreactive low thermal conductivity one, consequently the thermal profile imposed on the external side of the ampoule walls could differ from that obtained inside, particularly in the hyperbolic case. In that case, the axial and radial gradients changed drastically. Fig. 2 illustrates these aspects, however, because of the low magnitude of the different Prandtl numbers used, the generic isotherm distribution presented there is

not too much dependent on the semiconductor considered. Concerning the dopant concentration, upper and lower liquid–wall contacts have been assumed to be impermeable and the amount of dopant being incorporated in the solid phase has been taken to be proportional to the amount of dopant in the melt at the interface. That is to say, cs=kcl, being k the segregation coefficient. This last distribution at the interface, cl, has been evaluated, as usual, through the interface dopant balance, dc/dx= Sc W(1 k)S(t)c, being Sc the Schmidt number. Concerning the geometric parameters the thickness of the wall has been fixed in 1 mm, the diameter of the ampoule in 10 mm and the sample length in 80-short samples—190-medium samples—or 300 mm-long samples—respectively. Taking into account the maximal translation length of the actual standard Bridgman mg furnaces— typically 300 mm for the ESA’s POLIZON or the DLR’s AGAT—computational results conveniently covers a wide spectrum of results. At this respect, it must be pointed out that in the case of an external linear thermal profile the only way to attain similar dopant distribution in the three cases is using a slightly modified thermal Rayleigh number based on the thermal gradient, not on the absolute difference of temperatures between the hot and cold sides. Such universality is due to the fact that the velocity of the liquid phase along the growing interface is a decisive factor and it is

ARTICLE IN PRESS X. Ruiz, M. Ermakov / Journal of Crystal Growth 273 (2004) 1–18

5

Fig. 2. Effect of the walls on the distribution of the thermal field. Left/right side shows the isotherms numerically calculated without/ with walls (the double arrows indicate the sense of evolution of the liquid domain).

determined by the gradient of temperature in this region. No change of the thermal Rayleigh number is necessary in the case of external hyperbolic profiles or with the solutal Rayleigh number. These adjustments are not surprising; in fact similar considerations about have also been effected a few years ago in the literature [16]. The variation of the mg-vector inside real inertial flight mode orbiting systems has been modeled using a circularly polarized vibration. The eccentricity of the real orbit should be, thus, as low as possible to accomplish modeling conditions. Under these conditions both gravity components have been defined as harmonic signals with the same amplitude and frequency but shifted a phase of p=2: To be near to long orbit flight parameters— note, for instance, that the orbital period in the case of EURECA-1 satellite is near 90 min—the present work will focus on very low orbital frequencies between 100 and 250 mHz—equivalently, between 67 and 167 min of orbital period. Microgravity amplitudes cover also an enough

wide range of values typical of these kinds of environments, i.e., between 1 and 50 mg. Because of the aim of this choice is to clarify the impact on the segregation mainly generated by vibroconvective flows we simplify the reality considering only sinusoidal signals of zero average, equivalently, considering not relevant for the problem the action of a realistic 0.1 mg constant background along the growth process—in fact, only a drift to be systematically substracted to each signal. In relation with the different thermophysical parameters used we will like to make here a few considerations. To analyze the obtaining of the crystalline ZnyCd1 yTe we will approach the Zn–Cd–Te system inside the Te-rich region as a Zn-doped CdTe stoichiometric melt. It seems that the consideration of the much-simplified situation of zinc as a dilute specie in a CdTe melt properly supports experimental checking. But, because of experimental measurements indicate that the diffusion coefficient of zinc in solid CdTe at elevated temperatures is high, solid diffusion is

ARTICLE IN PRESS X. Ruiz, M. Ermakov / Journal of Crystal Growth 273 (2004) 1–18 2.0 In concentration

also accounted for in this case [17–19]. On the other hand for Bi-1 at% Sn alloy, thermosolutal convection has been considered to obtain more accurate results. This is because the low value of the segregation coefficient implies an important enhancing of tin concentration near the growing interface with the subsequent increasing of solutal convective effects in this small region near the interface [7,8]. However, the melting temperature has been regarded here as a constant supposing that tin concentration values at the interface would be sufficient reduced to not change it substantially. Obviously, the authors know that using a concentration-dependent melting temperature the solute increasing during the solidification causes a decrease in its value and consequently a reduction in the value of the growth velocity, but all these effects have not been considered here. So, in this sense, the present simplified thermosolutal calculations could be considered only as a first approach to the problem. Finally, note that no thermal dependencies of the different thermophysical parameters have been considered and that the quantitative values of the growth velocities have been chosen in order to totally ensure morphological stability of the different growing interfaces. Table 1 presents the quantitative values of the different nondimensional groups used in the calculations Taking into account the low values of the thermal Rayleigh numbers involved it has been decided to use a standard 2D finite volume formulation based on the semi-implicit SIMPLE scheme. The accuracy has been of the second order and we used 9600 (161  61) rectangular nonregular elements strongly clustered near the solid–liquid boundaries to well resolve flow gradients in these areas. In addition, solid walls of this conjugate problem have been recovered with a 3200 (161  11) rectangular regular elements. All cases started at rest, that is to say, both longitudinal and transversal velocity components equals initially zero. A linear distribution for the temperature, the basic conductive solution, and a uniform distribution of the dopant in the liquid phase have also been supposed in all cases as starting values. By simplicity, exhaustive thermohydrodynamic checking of the present algorithm

1.5 1.0 Experimental results [22] Present calculations (diffusive case)

0.5 0.0

0

1

(a)

2

4

3

Crystal length -cm1.0 Sn : Bi

Normalized axial macrosegregation

6

0.8

Experimental results [23] Low frequency vibration { g x = 0, g y = g 0sin2π f } - Numerical -normalized- conditions: Hyperbolic thermal profile g0 = 10 µg f = 30 µHz Vg = 5.8 mm/h

0.6

0.4

µg y

0.2

0

(b)

50

100 Time

150

200

Fig. 3. Experimental and computational solid segregation under: (a) diffusive, (b) vibroconvective conditions.

are omitted here, but additional information could be consulted in the proposed references [20,21]. Only a detail is presented here, Figs. 3a and b, to clearly show the potentialities of fully timedependent schemes, like the present one, in order to reproduce the growth results previously reported [22,23].

3. Results Because of the control of the segregation in the solid phase is the most challenging problem in relation with the growth of crystalline semiconductors we will concentrate results and discussions on this topic. In addition, to be as exhaustive as possible, the solid phase has been evolved systematically until roughly a 90% of the cavity becomes solidified. As an example of the results obtained, Figs 4a and b show several segregation patterns for short

ARTICLE IN PRESS X. Ruiz, M. Ermakov / Journal of Crystal Growth 273 (2004) 1–18

7

Sn : Bi

cmax : 1.001, cmin : k, deltac = 0.01 f = 0.25 mHz, T = 67 min (growth : 5.15 periods)

cmax : 1.0059, cmin : k, deltac : 0.01 f = 0.25 mHz, T = 67 min (growth : 5.15 periods)

cmax : 1.0013, cmin : k, deltac = 0.01 f = 0.20 mHz, T = 83 min (growth : 4.12 periods)

cmax : 1.0068, cmin : k, deltac : 0.01 f = 0.20 mHz, T = 83 min (growth : 4.12 periods)

cmax : 1.0212, cmin : k, deltac : 0.01 f = 0.15 mHz, T = 111 min (growth : 3.09 periods)

cmax : 1.0104, cmin : k, deltac : 0.01 f = 0.15 mHz, T = 111 min. (growth : 3.09 periods)

cmax : 1.0031, cmin : k, deltac : 0.01 f = 0.10 mHz, T = 167 min (growth : 2.06 periods) 1 microg [Ra(T) = 0.04; Ra(c) = 0]

cmax : 1.0164, cmin : k, deltac : 0.01 f = 0.10 mHz, T = 167 min (growth : 2.06 periods) 5 microg [Ra(T) = 0.205; Ra(c) = 0]

cmax : 1.0121, cmin : k, deltac : 0.01 f = 0.25 mHz, T = 67 min (growth : 5.15 periods)

cmax : 1.0657, cmin : k, deltac : 0.03 f = 0.25 mHz, T = 67 min (growth : 5.15 periods)

cmax : 1.0140, cmin : k, deltac : 0.01 f = 0.20 mHz, T = 83 min (growth : 4.12 periods)

cmax : 1.0771, cmin : k, deltac : 0.03 f = 0.20 mHz, T = 83 min (growth : 4.12 periods)

cmax : 1.0212, cmin : k, deltac : 0.01 f = 0.15 mHz, T = 111 min (growth : 3.09 periods)

cmax : 1.1152, cmin : k, deltac : 0.03 f = 0.15 mHz, T = 111 min. (growth : 3.09 periods)

0

1

2

3

4

5

cmax : 1.0334, cmin : k, deltac : 0.01 f = 0.10 mHz, T = 167 min (growth : 2.06 periods) 10 microg [Ra(T) = 0.41; Ra(c) = 0]

(a)

6

7 0

1

2

3

4

5

6

7

cmax : 1.1781, cmin : k, deltac : 0.03 f = 0.10 mHz, T = 167 min (growth : 2.06 periods) 50 microg [Ra(T) = 2.05; Ra(c) = 0]

(xy polarized rotational vibrations; linear profile)

Fig. 4. (a) Sn:Bi, (b) Bi–Sn solid segregation patterns obtained in the case of short samples after the virtual growth process. The different rectangles show the different duration of one orbital period, and the arrows, the direction of the growth evolution—by simplicity, no walls have been plotted here.

ARTICLE IN PRESS X. Ruiz, M. Ermakov / Journal of Crystal Growth 273 (2004) 1–18

8

Bi - Sn

cmax : 0.9380, cmin : k, deltac : 0.05 f = 0.25 mHz, T = 67 min (growth : 5.15 periods)

cmax : 1.1221, cmin : k, deltac ; 0.1 f = 0.25 mHz, T = 67 min (growth : 5.15 periods)

cmax : 0.9490, cmin : k, deltac ; 0.05 f = 0.20 mHz, T = 83 min (growth : 4.12 periods)

cmax : 1.2053, cmin : k, deltac ; 0.1 f = 0.20 mHz, T = 83 min (growth : 4.12 periods)

cmax : 0.9676, cmin : k, deltac ; 0.05 f = 0.15 mHz, T = 111 min (growth : 3.09 periods)

cmax : 1.3303, cmin : k, deltac ; 0.1 f = 0.15 mHz, T = 111 min. (growth : 3.09 periods)

0

1

2

3

4

5

6

7

0

1

2

3

4

5

cmax : 0.9968, cmin : k, deltac ; 0.05 f = 0.10 mHz, T = 167 min (growth : 2.06 periods) 1 microg [Ra(T) = 0.192; Ra(c) = -47.278]

cmax : 1.5248, cmin : k, deltac ; 0.1 f = 0.10 mHz, T = 167 min (growth : 2.06 periods) 5 microg [Ra(T) = 0.96; Ra(c) = -236.39

cmax : 1.4228, cmin : k, deltac : 0.1 f = 0.25 mHz, T = 67 min (growth : 5.15 periods)

cmax : 2.0673, cmin : k, deltac ; 0.15 f = 0.25 mHz, T = 67 min (growth : 5.15 periods)

cmax : 1.6058, cmin : k, deltac ; 0.1 f = 0.20 mHz, T = 83 min (growth : 4.12 periods)

cmax : 2.1228, cmin : k, deltac ; 0.15 f = 0.20 mHz, T = 83 min (growth : 4.12 periods)

cmax : 1.8260, cmin : k, deltac ; 0.1 f = 0.15 mHz, T = 111 min (growth : 3.09 periods)

cmax : 2.1797, cmin : k, deltac ; 0.15 f = 0.15 mHz, T = 111 min. (growth : 3.09 periods)

0

1

2

3

4

5

6

7 0

cmax : 2.1611, cmin : k, deltac ; 0.1 f = 0.10 mHz, T = 167 min (growth : 2.06 periods) 10 microg [Ra(T) = 1.92; Ra(c) = -472.779]

(b)

1

2

3

4

5

cmax : 2.2594, cmin : k, deltac ; 0.15 f = 0.10 mHz, T = 167 min (growth : 2.06 periods) 50 microg [Ra(T) = 9,60; Ra(c) = -2363.895

(xy polarized rotational vibrations; hyperbolic profile) Fig. 4. (Continued)

6

7

6

7

ARTICLE IN PRESS X. Ruiz, M. Ermakov / Journal of Crystal Growth 273 (2004) 1–18

9

Fig. 5. Snapshots of the streamfunction evolution at four different instants. Note the existence of a secondary vortical structure near the growing interface in the case of Bi–Sn alloy (the double arrows indicates the sense of the evolution of the liquid domain bounded by two upper and lower walls, not plotted by simplicity).

samples of Sn: Bi and Bi–Sn. The dopant distribution in the different virtually grown crystals clearly depends on the particular vibroconvective liquid reality driving dopant towards the growing interface. In effect, in the case of Sn:Bi and as it could be observed in Fig. 5a, big thermally driven convective cell Hadley type is the responsible of driving dopant toward the growing interface. On the contrary, in the case of Bi–Sn alloy, solutally driven flow near the growing interface redirects the solute towards this interface. This convective cell circulates solute-rich material near the solid–liquid interface as a consequence of the solute rejection during the growth process. As the time proceeds this secondary vortex increases its size and strength until it becomes more important than the generated by natural convection away the growing area. This fact, recently

reported for a constant gravity level, is also accomplished in the present case, that is to say, when a harmonic time-dependent gravity shakes the liquid melt. All these aspects could be easily observed in Fig. 5, but it must be emphasized that the flow patterns evolved in a harmonic virtual reality; consequently the sense of rotation of the different vortical structures depends on the value of the gravity at that instant. Although Figs. 4 and 5 show the computed behavior inside short samples, in the case of medium or large ones the situation is similar, so details are omitted here for the sake of brevity. A first qualitative attempt to organize the results generated by this multiparametric study has been based here on the consideration of four different basic scenarios. Strictly speaking, all these scenarios are not independent, i.e., they are different

ARTICLE IN PRESS 10

X. Ruiz, M. Ermakov / Journal of Crystal Growth 273 (2004) 1–18

stages of evolution depending on the flow intensity, but for simplicity in the description, the authors have maintained the independence. The first scenario, labeled type I, shows a situation in which no appreciable global impact on the solid phase is finally obtained. In the second scenario, labeled type II, dopant structures begin to appear in the solid phase, but without periodicity indicating the existence of a weak irregular interaction between dopant distribution and transport mechanisms inside the liquid phase at high temperature. The third scenario, labeled type III, is characterized by the existence of periodically alternated dopant structures in the solid phase as a combined result of the convective mechanisms— vertical extent—and the growth velocity—horizontal extents. That unfortunately implies dopant stalagmithic-like structures along the crystal, or equivalently radial and longitudinal pernicious segregation. Finally the fourth scenario, labeled type IV, is similar to the third one but because of the high convective liquid level, the dopant structures spreads trying to attain the opposite side of the domain and loosing the center line symmetry. As an example Table 2 synthesizes the virtual growth campaign for the case of short samples and for the two types of above-defined thermal profiles, linear and hyperbolic. The left side corresponds to the linear results while that the right side corresponds to the hyperbolic ones. Looking at that table and in the case of the Ge:Ga, it could be noted that a gravity level of 1 mg finally generates a scenario type I for all frequencies as a consequence of the effect of the quasi diffusive flow regime established in the melt coupled with the very low value of the distribution coefficient. If the gravity level increases to 10 or 50 mg, the increasing of the convective level in the liquid phase is not enough to generate a minimal regular impact and only type II scenario appears for the lowest frequency cases. In this way, the theoretical prediction of the gallium distribution in the solid phase under inertial flight mode conditions amazingly gives a relative homogeneity in the vertical— radial—direction but not in the horizontal—longitudinal—one. These results are practically independent on the thermal boundary conditions used.

In the case of GaSb:Te and despite of the low values of the thermal Rayleigh numbers, type III scenarios have been computationally obtained for all frequencies using both kinds of thermal boundary conditions. For the case of selenium inside the gallium arsenide matrix the things are totally different. In effect, even in the lowest 1mg case, scenarios type II and III appeared at low frequencies in a way practically independent of the thermal conditions applied. More strongly effects could be observed in the cases of Indium phosphate doped with sulfur and in the tin bismuth diluted alloy. In the case of cadmium telluride—distribution coefficient greater than one—there is no interface rejection and the global tendency of the dopant level is obviously a decreasing one. Zinc segregation appears insensitive at 250 mHz—only type I scenarios arose, but as the frequency diminishes and the gravity level increases only structures type II appears, practically independent on the thermal boundary conditions. The absence of type III structures in clearly convective regimes is a direct consequence of the intense solid diffusion strongly favored with the long time needed to grow the crystals—see Fig. 6. Finally, it must be pointed out that in the case of the nondiluted tin bismuth alloy the nonnegligible value of the concentrational Rayleigh numbers change the pattern of segregation in the solid phase and despite the low 1 mg value in which type II scenarios are reported, the rest of cases presents a clear type IV scenarios in the solid phase. As before all the above-discussed characteristics are practically independent of the length of the sample. With the aim to evaluate the reliability of a simplified analytical approach using scaling arguments we plan to work with, as usual, the longitudinal segregation defined cavg ðxÞ ¼ R1 c ðx; yÞ dy [24–26]. However, as it could be 0 s noted in Fig. 7, the averaging procedure introduces in the present cases problematic masking effects if comparing, for instance, with the solid concentration profiles recorded at different points of the growing interface. Longitudinal segregation is of lower magnitude than the one corresponding to punctual concentrations and its shape is more regular and harmonic than the shape of whatever

ARTICLE IN PRESS X. Ruiz, M. Ermakov / Journal of Crystal Growth 273 (2004) 1–18

11

Table 2 Computationally predicted solid dopant scenarios in the case of short samples Linear—hyperbolic thermal profiles

Ge:Ga

GaSb:Te

CdTe:Zn

GaAs:Se

InP:S

Sn:Bi

Bi–Sn

T/F-

1 mg

5 mg

10 mg

50 mg

min/mHz

RaT/Rac

RaT/Rac

RaT/Rac

RaT/Rac

67/250 83/200 111/150 167/100

0.4231/0 I–I I–I I–I I–I

2.1160/0 I–I I–I I–II II–II

4.2312/0 I–I I–I II–II II–II

21.1560/0 I–I I–I II–II II–II

67/250 83/200 111/150 167/100

0.0072/0 I–III I–III I–III I–III

0.0362/0 III–III III–III III–III III–III

0.0724/0 III–III III–III III–III III–III

0.3619/0 III–III III–III III–III III–III

67/250 83/200 111/150 167/100

2.3259/0 I–I I–I I–I I–I

11.6259/0 I–I I–I I–I II–II

23.2591/0 I–I I–I I–I II–II

116.2953/0 I–I I–I II–II III–III

67/250 83/200 111/150 167/100

0.3355/0 I–II II–II II–III II–III

1.6776/0 II–III II–III III–III III–III

3.3553/0 II—II III–III III–III III–IV

16.7764/0 II–III III–III III–III IV–IV

67/250 83/200 111/150 167/100

0.5147/0 III–III III–III III–III III—III

2.5734/0 III–III III–III III–III III–III

5.1469/0 III—IV III–IV III–IV III–IV

25.7346/0 IV–IV IV–IV IV–IV IV–IV

67/250 83/200 111/150 167/100

0.0406/0 I–III III–III III–III III–III

0.2028/0 III–III III–III III–III III–III

0.4055/0 III–III III–III III–III III–III

2.0276/0 III–III III–III III–III III–IV

67/250 83/200 111/150 167/100

0.192/ 47.278 II–II II–II II–II II–II

0.96/ 236.39 III–III III–III IV–IV IV–IV

1.92/ 472.779 IV–IV IV–IV IV–IV IV–IV

9.60/ 2363.895 IV–IV IV–IV IV–IV IV–IV

of the above-mentioned concentration profiles, even the frequency could also be masked in some cases. General quantitative considerations have thus been effected using two longitudinal profiles of concentration at two different heights—see Fig. 8. The first immediate comment is that both longitudinal profiles are not symmetrical and that

in the cases in which a strong incorporation/ rejection of dopant away of the interface is generated during the continuous solidification process—Ge:Ga, GaAs:Se, CdTe:Zn and Bi–Sn, the dopant level progressively decreases/increases in the liquid phase and consequently in the solid one as a consequence of the closedness of the

ARTICLE IN PRESS 12

X. Ruiz, M. Ermakov / Journal of Crystal Growth 273 (2004) 1–18

Fig. 6. Effect of the solid diffusion in the distribution of Zn inside the CdTe matrix. The arrows indicate the direction of the growth evolution and, also by simplicity, no walls have been plotted here.

domain itself. Then the absence of a well defined ‘plateau’ indicates that, at least from a segregation point of view, no steady state is attained in short samples. A second remark is that the longitudinal profiles of concentration corresponding to type III and IV scenarios exhibit an oscillatory behavior which follows the orbital frequency with amplitudes strongly depending on the magnitude of the external gravity. If the gravity level increases, the amplitude of the oscillation increases too, but simultaneously the degree of nonlinearity of the profile, i.e., several harmonics appears trying to resume the combined effect of the two vibrations mutually orthogonal. This fact is independent on the semiconductor considered and on the geometric and thermal boundary conditions applied. As

it could be observed in Figs. 9a–c the behavior of the concentration profiles in medium and long samples is similar to the short ones, so no new relevant considerations could be added at this respect. Only to another time emphasize here the smoothing effect of a strong solid diffusion in the CdTe:Zn medium and long cases.

4. Discussion First of all it must be pointed out that, from a qualitative point of view, the presented results completely explains why long crystal growth experiments should be carefully planned, but, except in the case of Fig. 3, no other detailed

ARTICLE IN PRESS X. Ruiz, M. Ermakov / Journal of Crystal Growth 273 (2004) 1–18

GaAs : Se (50 µg, 100 µHz)

GaSb : Te (50 µg, 100 µHz)

cs(x, 0.15) cs(x, 0.75) cavg(x)

1.00

1.2

0.75 0.8

cs

cs

13

0.50 0.4 0

20

40

60

80

0

50

100

150

Time

200

250

Time

Fig. 7. Possible problems appearing with the consideration of the longitudinal segregation in vibroconvective problems.

Scenario type III

Scenario type I

0.6

0.4

Ge: Ga (10 µg, 200 µHz, hyperbolic)

0.5

GaSb: Te (50 µg, 250 µHz, linear)

csol(x, 0.75)

csol(x, 0.15)

1.0

InP: S (10 µg, 150 µHz, linear)

0.2

Sn : Bi (50 µg, 150 µHz, hyperbolic)

Ge: Ga (50 µg, 200 µHz, linear)

0

1

2

3

4

5

6

7

0

Scenario type II

1.0

1

2

Scenario type IV

3

4

5

6

7

0.0

InP: S (50 µg, 100 µHz, hyperbolic)

1.5

InP: S (50 µg, 150 µHz, linear)

GaAs: Se (1 µg, 100 µHz, linear)

csol(x, 0.75)

csol(x, 0.15)

1.0 GaAs: Se (1 µg, 250 µHz, hyperbolic)

0.5 Ge: Ga (10 µg, 150 µHz, hyperbolic)

0.5 Bi: Sn (10 µg, 150 µHz, hyperbolic) Bi- Sn (5 µg, 100 µHz, linear)

0.0

0.0 0

1

2

3

x

4

5

6

7

0

1

2

3

x

4

5

6

7

Fig. 8. Normalized longitudinal concentration profiles for several different semiconductors at two different heights, y ¼ 0:15 and 0.75.

quantitative comparison is presented here to validate calculations because the lack of reliable quantitative information concerning this, till now, forgotten effect. The computational work presented here is, thus, intrinsically predictive in nature.

In connection with the generic problem of the tolerability limits for the different materials and gravity levels, it could be very useful for experimentalists to have an easy quick tool to predict the maximum level of dopant to be obtained under similar conditions. To answer

ARTICLE IN PRESS X. Ruiz, M. Ermakov / Journal of Crystal Growth 273 (2004) 1–18

14

Growth length (cm)

0 1.0

4

7

11

14

18

22

25

29

T

AR = 30, S(t f ) = 28

Ra =79.33 T

Ra =50.25 AR = 19, S(t f ) = 18

0.8

T

Ra =21.16

csol (x, 0.15) / c0

AR = 8, S(t f ) = 7

0.6

0.4

0.2

Ge : Ga ( linear ; 50µg , 250µHz) 0

10

0

4

20

30

(a)

8

40 50 Growth time (hours)

Growth length (cm) 12 16

60

20

70

80

24

28

InP : S ( linear; 50 µg, 250µHz)

csol (x, 0.15) / c0

1.2

0.8

T

Ra =25.73

AR = 8, S(t f) = 7

0.4

T

Ra =61.12

T

AR = 19, S(t f) = 18

Ra =96.50

AR = 30, S(t f) = 28

0.0

(b)

0

2

4

7

9 11 Growth time (hours)

13

16

Fig. 9. Short, medium and long concentration patterns of (a) Ge:Ga, (b) InP:S and (c) Bi–Sn under different conditions (the arrows indicate the direction of the growth evolution and no walls have been plotted by simplicity).

ARTICLE IN PRESS X. Ruiz, M. Ermakov / Journal of Crystal Growth 273 (2004) 1–18

15

Growth lenght (cm) 0

5

9

14

18

23 27 T c Ra =9.6; Ra =-2364

3.0 AR = 30, S(t f ) = 28

2.5

T

Rac =9.6 Ra =-2364

AR = 19, S(t f ) = 18 T

csol (x, 0.15) / c0

AR = 8, S(t f ) = 7

Rac =9.6 Ra =-2364

2.0

1.5

1.0

Bi - Sn ( hyperbolic ; 50µg , 100µHz)

0.5

0.0

(c)

0

4

7

11

15

19

22

Growth time (hours) Fig. 9. (Continued)

this pragmatic—as well as simplistic—question a two-step strategy is applied here. In the first step, upper and lower boundaries for each one of the three lengths considered will be determined. Using this information the second step will fit finally the variation between a minimum and maximum absolute values. The raw prediction is then finished. To illustrate this strategy Fig. 10 shows results for the first step in the case of short samples. The lower boundaries are, as it is also indicated in the picture, the values of the corresponding distribution coefficients. The methodology used to trace the upper boundaries is very simple and it has been based on the consideration of the maximum values of the solid dopant concentration for each one of the different virtual samples. After that, and taking into account only the lowest and the highest values of all data— points 1 and 2 in Fig. 9—a region has been bounded for each material using the simplest function, i.e., a straight line, for all cases in which thermal convection is the main convective mechanism in the melt. In the cases in which the solutal

convection is relevant a parabola fits better the range of conditions analyzed here. This difference is not unusual when the solutal convection is dominant, in fact this behavior justify that due to the well-known peculiarities of the radial segregation, with a maximum in the poor mixing region, solutal convection could be particularly pernicious in microgravity environments [27]. Finally, note that because the distribution coefficient of the Zn in CdTe is greater than the unity, what is computed by this method is the lowest boundary. The upper one is, as in the other cases, the constant value of the distribution coefficient. After finishing this first set of calculations, each semiconductor has associated three straight lines/ parabolas, which should be reduced to only one simply considering the one, which delimitates the highest values of the concentration. Table 3 presents the final results of this process for the seven semiconductors used, enabling predictive raw fitting between 1 and 50 mg and 100 and 250 mHz.

ARTICLE IN PRESS X. Ruiz, M. Ermakov / Journal of Crystal Growth 273 (2004) 1–18

16

1.0

Theoretical Ga dopant threshold (short samples): csol, max(x,y) = 0.025 log10(µg) + 0.619

2

k

0.9

1

k µ Hz 100 150 200 250 Linear

0.60

0.0

0.4

µ Hz 100 150 200 250 Hyperbolic

0.8

Ge : Ga 1.2

2

GaAs : Te

1.6

csol, max(x,y)

Theoretical Te dopant threshold (short samples): csol, max(x,y) = 0.353 log10(µg) + 1.037

0.8

Theoretical Zn dopant threshold (short samples): csol, max(x,y) = - 0.106 log10(µg) + 0.946 0.0

1.6

0.4

0.8

1.2

csol, max(x,y)

1

2 1.6

Theoretical Sn dopant threshold (short samples): 2 csol, max(x,y) = - 0.575 log10(µg) + 1.692 log10(µg) + 1.044 3

2.4

2

1.4

1.6

k 1.2

csol, max(x,y)

csol, max(x,y)

0.65

CdTe :Zn

k 1

1

Bi - Sn

1.0

0.8

0.0

0.4

0.8

log10(µg)

1.2

1.6

0.0

0.4

0.8

1.2

log10(µg)

1.6

Fig. 10. Computationally predicted dopant thresholds for several semiconductors in the case of short samples (mg=1, 5, 10 and 50 in this plot). Table 3 Final quantitative predictions of the different upper thresholds Material

a

csol, max(mg)=a log10 (mg)2+b log10 (mg)+c(with mg=1, 5, 10, 50) Ge:Ga u 0 GaSb:Te u 0 CdTe:Zn l 0 GaAs:Te u 0 InP:S u 0 Sn:Bi u 0 Bi–Sn u

0.575

b

c

0.025 0.029

0.106 0.353 0.268 0.253 1.692

0.979 1.001 1.246 1.037 1.038 1.015 1.464

u/l—upper/lower dopant threshold (the value of the segregation coefficient is the other boundary).

5. Conclusions The present multiparametric computational study based on seven relevant semiconductor

materials, clearly conclude the existence of a strong link between the orbital characteristics of a generic inertial flight mode orbiting system and the final dopant distribution in virtually grown

ARTICLE IN PRESS X. Ruiz, M. Ermakov / Journal of Crystal Growth 273 (2004) 1–18

crystals using a Bridgman furnace. This relation has been defined with the help of four scenarios, which appears sequentially as the vibroconvective level increases in the liquid phase at high temperature. Thus, to improve crystal homogeneity in long crystal growth experiments inside inertial flight mode spatial platforms some strategies should be seriously considered in order to efficiently mitigate the pernicious effects of these very low orbital frequencies. Perhaps these strategies could range between passive isolation mounts with cutoff frequencies even lower than 0.01 Hz [28,29] or smart mechanical devices constant and slowly orienting the axis of the growth system in order to align it with the direction of the residual gravity [30–33].

Acknowledgements This work was supported by DGCYT Grant BFM2003-00657 and MCYT Grant BQU200305042-C02-02.

References [1] E.S. Nelson, An examination of anticipated g-Jitter on Space Station and its effects on materials processes, NASA Technical Memorandum 103775, 1991. [2] A.N. Danilewsky, K.W. Benz, InP growth from In solutions under reduced gravity, J. Crystal Growth 97 (1989) 571. [3] A.N. Danilewsky, St. Boschert, K.W. Benz, The effect of the Orbiters attitude on the mg-growth of InP crystals, Microgravity Sci. Technol. X/2 (1997) 106. [4] St. Boschert, A.N. Danilewsky, K.W. Benz, Numerical simulation of the influence of the orbiters attitude on the mg-growth of InP:S crystals from an In solution during the EURECA-1 flight, J. Crystal Growth 205 (1999) 92. [5] J.I.D. Alexander, S. Amiroudine, J. Ouazzani, F. Rosenberger, Analysis of the low gravity tolerance of Bridgman–Stockbarger crystal growth. II. Transient and periodic accelerations, J. Crystal Growth 113 (1991) 21. [6] A.I. Fedoseyev, J.I.D. Alexander, Investigation of vibrational control of convective flows in Bridgman melt growth configurations, J. Crystal Growth 211 (2000) 34. [7] V. Timchenko, P.Y.P. Chen, E. Leonardi, G. de Vahl Davis, R. Ababbaschian, A computational study of transient plane front solidification of alloys in a Bridgman apparatus under microgravity conditions, Int. J. Heat Mass Transfer 43 (2000) 963.

17

[8] V. Timchenko, P.Y.P. Chen, E. Leonardi, G. de Vahl Davis, R. Ababbaschian, A computational study of binary alloy solidification in the MEPHISTO experiment, Int. J. Heat Mass Transfer 23 (2002) 258. [9] R. Monti, R. Savino, Microgravity experiment acceleration tolerability on space orbiting laboratories, J. Spacecrafts Rockets 33 (1996) 707. [10] Y. Zhao, I. Alexander, Effects of g-jitter on experiments conducted in low-earth orbit: a review, Proceeding of the 41th AIAA Aerospace Sciences Meeting and Exhibit, Reno, Nevada, 2003. [11] X. Ruiz, M. Ermakov, Segregation analyses of semiconductor single crystals grown inside modeled inertial mode orbiting systems, European Low Gravity Society, ELGRA, Mu¨nchen, Germany, 2003. [12] X. Ruiz, M. Ermakov, Inertial flight mode and semiconductor segregation patterns, 54th International Astronautic Congress, Bremen, Germany, 2003. [13] S.A. Nikitin, V.I. Polezhaev, A.I. Fedyushkin, Mathematical simulation of impurity distribution in crystals prepared under microgravity conditions, J. Crystal Growth 52 (1981) 471. [14] S.A. Nikitin, V.I. Polezhaev, A.I. Fedyushkin, Convection and impurity distribution in crystal growth in low-gravity environments, Adv. Space Res. 3 (1983) 65–78. [15] Ta Wei Fu, W.R. Wilcox, Rate change transients in Bridgman Stockbarger growth, J. Crystal Growth 51 (1981) 557. [16] A.I. Feonychev, G.A. Dolgikh, Effects of constant and variable accelerations on crystal growth onboard Spacecraft by the method of directional crystallization, Cosmic Res. 39 (2001) 365. [17] S. Kuppurao, S. Brandon, J.J. Derby, Modeling the vertical Bridgman growth of cadmium zinc telluride. I. Quasi-steady analysis of heat transfer and convection, J. Crystal Growth 155 (1995) 93. [18] S. Kuppurao, S. Brandon, J.J. Derby, Modeling the vertical Bridgman growth of cadmium zinc telluride. II. Transient analysis of zinc segregation, J. Crystal Growth 155 (1995) 103. [19] A. Yeckel, J.J. Derby, Buoyancy and rotation in smallscale vertical Bridgman growth of cadmium zinc telluride using accelerated crucible rotation, J. Crystal Growth 233 (2001) 599. [20] M. Ermakov, X. Ruiz, New possibilities of computer laboratory COMGA for modeling of convective processes, Grid Generation: Theory and Applications—Proceedings, Computing Center RAS, Moscow, Russia, June 2002, p. 207. [21] M. Ermakov, M.S. Ermakova, X. Ruiz, Numerical Modeling in Crystal Growth Processing, Fifth International Conference on Single Crystal Growth and Heat and Mass Transfer, Obninsk, Russia, September 2003, p. 537. [22] Th. Duffar, M.D. Serrano, C.D. Moore, J. Camasel, S. Contreras, P. Dusserre, A. Rivoallant, B.K. Tanner, Bridgman solidification of GaSb in space, J. Crystal Growth 192 (1998) 63.

ARTICLE IN PRESS 18

X. Ruiz, M. Ermakov / Journal of Crystal Growth 273 (2004) 1–18

[23] J.P. Garandet, J.I.D. Alexander, S. Corre, J.J. Favier, Composition variations induced by g-jitter in Bridgman growth of Sn–Bi alloys in microgravity, J. Crystal Growth 236 (2001) 543. [24] J.P. Garandet, J.J. Favier, D. Camel, Segregation Phenomena in Crystal Growth from the Melt, in: D.T. Hurle (Ed.), Handbook of Crystal Growth, vol. 2, 1994. [25] J.P. Garandet, S. Corre, S. Gavoille, J.J. Favier, J.I.D. Alexander, On the effect of gravity perturbations on composition profiles during Bridgman crystal growth in space, J. Crystal Growth 165 (1996) 471. [26] J.P. Garandet, S. Corre, J.J. Favier, Acceptable microgravity levels in Bridgman solidification—scaling analyses and experiments, Microgravity Sci. Technol. XI/2 (1998) 59. [27] R.A. Brown, Theory of transport processes in single crystal growth from the melt, A.I.Ch.E. J. 34 (1988) 881.

[28] J. de Carufel, Performance of the microgravity vibration isolation mount on the Space Shuttle and the Mir, Proceedings of Spacebound, Canada, 2000. [29] M.S. Whorton, g-Limit: a microgravity vibration isolation for the International Space Station, Proceedings of Spacebound, 2000, Canada. [30] R. Fortezza, R. Monti, R. Savino, Residual g orientation—as feasibility study on a concept aimed at mitigating the effects of the quasi-steady accelerations on board ISS, Microgravity Sci. Technol. XI/3 (1998) 119. [31] R. Monti, R. Savino, g-Sensitivity of microgravity experimentation—fundamentals of disturbance response, Microgravity Sci. Technol. XI/2 (1998) 53. [32] R. Savino, Residual-g and g-jitter effects on the measurement of thermophysical properties in microgravity, Adv. Space Res. 29 (2002) 559. [33] R. Monti, D. Paterna, R. Savino, Counter-measures to mitigate residual g-effects on microgravity experiments on the Space Station, Acta Astronaut. 50 (2003) 209.