Problem of attaining constant impurity concentration over ingot height
expand article infoMichael A. Gonik, Florin Baltaretu§
‡ Centre for Material Reseaches «Photon», Aleksandrov, Russia
§ Technical University of Civil Engineering of Bucharest, Bucharest, Romania
Open Access


The possibility of growing crystals with homogeneous impurity distribution over crystal height has been demonstrated in a study of segregation during silicon and germanium growth from thin melt layers using the submerged heater method. Numeric simulation of 200 mm diam. antimony-doped germanium crystallization has shown that, beginning from a 40 mm melt layer thickness, the exact problem solution with convection allowance is identical to the unidimensional heat exchange problem solution in the central ingot part. The conditions under which convection can be ignored in mass transport calculation are more rigorous: the melt layer height must be within 20 mm. In this case Tiller’s ratio can be used for calculating the longitudinal impurity distribution for predominantly diffusion-controlled mass transport pattern. Analysis of the existing attempts to describe the experimental crystal growth results using the simplified formulae shows that they only yield acceptable results if the actual growth rate or change in melt layer thickness during crystallization are taken into account, e.g. as in the formula suggested by Marchenko et al. One can therefore analytically describe the longitudinal impurity distribution in the ingot, e.g. B and P distribution in silicon, and recommend the degree of additional doping of the melt zone under the heater so that to provide a constant impurity concentration over the ingot height. Homogeneous material can be obtained after residual layer solidification in the end portion of the ingot if the growth rate is controlled through varying the cooling rate.


mathematical simulation, impurity diffusion, directional crystallization, submerged heater, germanium, multicrystalline silicon

1. Introduction

The directional crystallization method used for the production of multicrystalline silicon does not provide a homogeneous doping impurity distribution along the growing crystal. Detailed analysis [1] showed that the concentration of impurities such as boron and phosphorus which is described adequately well by Shell’s ratio grows slightly at an early growth stage, accumulates during the crystal growth and increases abruptly in the end portion of the ingot. This latter zone is typically 15–20 % of the ingot length which makes at least 30–50 mm of the material at the top of the ingot low-quality. Impurity concentration change along the ingot affects the resistivity of silicon and is therefore detrimental but this problem cannot be solved with standard methods. However, the task of producing a height-homogeneous ingot can be solved by using the melt-submerged heater technique for directional crystallization growth which is implemented either as the submerged heater method (SHM, [2]) or the axial heat flow method (AHP, [3]). The heater is placed at a short distance from the phase boundary and partitions the melt zone in two (the crystallization zone under the heater and the feeding zone above the heater), thus fundamentally changing the stream pattern. Analysis of the prerequisites for the technical implementation of this approach to producing multicrystalline silicon using the AHP method was reported earlier [4]. The abovementioned paper [4] reported results of testing a protective coating [5] of submerged heater casing for silicon crystallization from a standard G1 crucible.

Crystal growth from a thin melt layer significantly reduces convection intensity [6] and provides for a segregation pattern that is close to a diffusion one [7]. It was shown with assumptions for experimentally grown Te-doped GaSb [8] that even for crystals grown on Earth the impurity concentration distribution CS (x) can be described adequately well by an expression derived for the specific implementation conditions, i.e., the so-called diffusion controlled mode [9]:


where C0 is the initial impurity concentration in the melt, V is the growth rate, k is the equilibrium segregation ratio, D is the diffusion coefficient and x is the coordinate along the crystal axis.

Segregation in a thin melt layer reaches an equilibrium condition considerably faster [10], and under specific conditions the effective segregation ratio keff becomes unity [1]. This process was analytically described earlier [12] with the following formula:


It is assumed herein that different melt zones (under the heater and above the heater) had the impurity concentrations C1 and C2, respectively, CS is the initial impurity concentration in the melt, and h is the melt layer thickness. Obviously, if the feeding zone is doped to the concentration C2 = CS and the growth zone is doped so that the ratio C1 = C2/k is true, the crystal will have a constant impurity concentration equal to the impurity concentration in the feeding zone. Unfortunately, at the final stage of crystallization when there is no more fresh melt in the feeding zone and the submerged heater does not move relative to the crucible, Eq. (2) is no longer valid. However, nor true is the earlier suggested ratio [13] for the complete mixing case since for the SHM or the AHP methods the submerged heater is not removed from the melt even after crucible pulling is stopped. There is no free surface and hence no Marangoni convection, and the diffusion contribution to mass transport remains prevailing.

Thus there is a need for a numeric simulation to study longitudinal segregation under the conditions set forth above and for determining the optimum crystallization heating mode providing for the growth of a homogeneous crystal over the maximum possible ingot distance. This study will determine the applicability limits of the impurity distribution analytical expressions to be further used for engineering calculations of silicon growth. Though only the “residual” melt portion is in question, the length of this zone may be as large as 100 mm or greater if large-size multicrystalline silicon ingots are grown because the melt layer height for the submerged heater method is chosen based on the h/d ~ 0.15 ratio. The model material for this study was antimony doped germanium since the antimony segregation ratio is 0.003, i.e., a priori lower than any of the segregation ratios of commercially important impurities in silicon, and furthermore the problem in hand becomes the most complex for small k. Moreover, the properties of germanium have been well studied and there is the possibility to conduct, if necessary, experiments at lower temperatures than would be required for silicon.

2. Mathematical model and numeric approaches

For submerged heater crystallization, by analogy with the conventional Bridgman method, the crucible is lowered into the cold zone of the growth chamber (Fig. 1a) and the crystal is grown until the melt in the zone above the heater is over. The unidimensional heat transfer problem can be described with the following equation:

Figure 1.

Crystal growth setup at (a) crucible pulling stage relative to submerged AHP heater and (b) at the end after stopping, and (c) domain schematic.


where λL,S is the liquid and solid state heat conductivity, ρL is the melt density, J is the crystallization heat, Tm is the melting point, Thot and Tcool are the temperatures at the hot and cool boundaries of the melt/crystal system at the AHP heater bottom and the crucible, respectively, and H is the height of the grown crystal. To ensure constant crystallization conditions at the phase boundary, including constant axial temperature gradient, one should maintain h and Thot constant. Meanwhile the cool boundary temperature is gradually reduced; for nontransparent melts, e.g. germanium and silicon, the reduction obeys a linear law in accordance with Eq. (3).

The heater is typically removed from the melt at the end of the process to avoid its freezing into the crystal, so the final portion crystallizes relatively rapidly and is later cut off from the ingot. Below we will consider a problem (Fig. 1b) for which crystal pulling is already stopped but the submerged heater is still in contact with the melt. This configuration allows controlling crystallization until its completion thus maintaining the required temperature at the top Ttop(t) and bottom Tbot(t) boundaries of the domain. It can be seen from Fig. 1 c that the domain for an axially symmetrical configuration is a half-disc with the diameter d and the height Δ which is cooled through simultaneous reduction of the temperatures at its ends, with the cooling rate being time-variable in the general case. The initial impurity distribution in the melt is accepted to be homogeneous over the entire domain and is set to unity.

We used the ANSYS Fluent software module for calculating crystallization and melting with a constant grid [14]. This module is designed for melts having a solidification zone between the solid and the liquid phases. To justify the use of this approach in this work, we considered a small temperature range (1 K) around the melting point, and therefore the solidus TS and liquidus TL temperatures were set as follows:

T S = Tm – 0.5; TL = Tm + 0.5.

The specific mixture enthalpy hmix expression can be written as follows:

where cS and cL are the specific heat capacities of the solid and the liquid phases, J is the latent melting heat and fL is the liquid phase fraction which is

The energy conservation equation will be written as follows:


Equation (4) is solved taking into account a linear temperature change at the bottom Tbot(t) = Tbot(0) – Rcoolt and the top Ttop(t) = Ttop(0) – Rcoolt surfaces, respectively. We assume that cooling has the constant rate Rcool. The convection heat exchange at the side boundary is

where α is heat convection coefficient.

The momentum conservation equation is solved using the porosity enthalpy method which considers the solidification zone as a porous medium. The momentum absorption term added to Eq. (4) as a flow for bringing the rates in the solid phase zones to zero becomes as follows:

where ɛ = 0.001 is the small number used for avoiding division by zero and Amush is the solidification zone constant.

The properties of germanium (density, heat capacity and heat conductivity) used for the calculations are shown in Table 1 with allowance for their stepwise change at the phase boundary:

Physical properties of germanium and silicon used in calculations.

Properties Notation Germanium [15–17] Silicon [18,19]
Melting point, K Tm 1210 1683
Density, kg/m3:
Crystal ρS 5260 2330
Melt pL 5550 2530
Latent melting heat, kJ/kg J 460 164
Heat capacity, J/(kg × K):
Crystal C S 400 700
Melt C L 1140 1100
Heat conductivity, W/(m × K):
Mass specific heat capacity, J/(kg × K)
Crystal λS 17 22
Melt λL 39 67
Kinematic viscosity, m2/s ν 1.35 × 10-5 1.5 × 10-5
Heat expansion coefficient, K-1 β 0.97 × 10-4 1.32 × 10-4
Segregation ratio k
Sb 0.003
P 0.35
B 0.8
Diffusion coefficient, cm2/s D
Sb 1 × 10-4
P 3.3 × 10-4
B 2.7 × 10-4


The impurity distribution is set in a similar manner using mixture concentration description:

where the generalized scalar transfer equation


allows calculating the impurity distribution by introducing the UDS user defined function for the concentration. The transfer equation for the solidification zone imitating the crystallization front allows for a segregation-related impurity source which was calculated as follows:


The parameters for which the numeric calculations were carried out are summarized in Table 2.

Numeric calculation parameters.

Parameter Notation Germanium Silicon
Heat convection coefficient, W/(m2 × K) α 25 25
Ambient temperature, K T 800 1300
Cooling rate, K/h R cool 10–60 15–40
Initial impurity concentration in eq. (1), cm-3 С 0
Sb 2.7 × 1017
P 6.0 × 1017
B 2.62 × 1017
Initial impurity concentration in Eq. (2) above AHP Heater, cm-3 С 2
P 2.1 × 1017
B 2.1 × 1017
Initial impurity concentration in Eq. (2) under AHP heater, cm-3 С 1
P 2.62 × 1017
B 6.0 × 1017

3. Applicability limits of simplified expressions

Figure 2 shows typical established temperature field patterns in the domain, stream function in the melt and velocity vectors certain time after crystallization onset. Since the top and bottom boundaries of a thin layer (Δ/d < 0.15) are maintained under isothermal conditions, no convection occurs in the middle of the domain. Convection driving force only emerges at the edge, due to the difference between the side boundary temperature and the ambient temperature, and a vertex with a clockwise stream forms at the wall. The bulk temperature field is unidimensional and only undergoes distortion at the edge. The swirl only exists at the side wall and can be considered to be localized in the domain having Δ = 10 mm. With an increase in Δ (Fig. 3) the stream penetrates deeper and can even affect the middle part. This however does not influence the temperature profile in most of the domain. Therefore convection can be ignored for heat transfer calculations with such Δ/d.

Figure 2.

(a) Temperature distribution. (b) stream function and (c) velocity vectors at periphery of domain in 500 s after Ge solidification onset.

Figure 3.

Stream velocity (m/s) in the domain at an early crystallization stage for different melt thicknesses: (a) 10 mm, (b) 20 mm and (c) 40 mm.

Interestingly, the stream which forms in quite a short time even for Δ = 40 mm as can be seen from Fig. 4 changes but slightly thereafter, until certain portion of the domain crystallizes and the melt layer starts depleting noticeably. This is well illustrated by Fig. 5 where it can be seen how the vertex affected zone decreases in size, the convection intensity drops and the stream structure degenerates with the growth of the crystal. A completely different pattern is observed if one considers the effect of convection on mass transport. It can be seen from Fig. 4b that the initially homogeneous impurity distribution undergoes distortion in a short time due to convection and continues changing, including in the middle of the domain. It is therefore difficult to conclude about a predominantly diffusion pattern of mass transport, at least for Δ = 40 mm. The situation is completely different if the domain height is 10 mm (Fig. 6). The convection emerging at the side wall causes transfer impurity such that its excess quantity is accumulated immediately behind the vortex where an antimony-enriched crystal grows. Meanwhile closer to the middle of the domain the inhomogeneity is completely eliminated and the concentration isocontours are parallel to the flat phase boundary in the middle part. Thus a diffusion crystallization mode is established in almost the entire zone and the axial antimony distribution in germanium calculated using Eq. (1) will hold in almost the whole ingot. This is confirmed by the axial impurity profile pattern shown in Fig. 7. It can well be seen that the antimony distribution profiles over the height of the domain (in the crystal and in the melt) is almost the same for different crystallization stages until a 75 mm radius is achieved. Only the extreme edge portions exhibit differences in the absolute values and the pattern of the profile.

Figure 4.

Formation dynamics for (a) velocity field and (b) antimony concentration in germanium during crystal growth from 40 mm thick melt layer in 500, 750 and 1000 s after solidification onset. Initial concentration in the domain is unity.

Figure 5.

Stream pattern in germanium melt for (a) ∆ = 10 mm and (b) ∆ = 20 mm and its degeneration with decreasing crystal growth melt layer thickness.

Figure 6.

Change in antimony relative concentration in melt and growing germanium crystal in (a) 500, (b) 1250 and (c) 1750 s after solidification onset for domain with ∆ = 10 mm.

Figure 7.

Axial change in antimony relative concentration in (a, c, e, g) 1000 and (b, d, f, h) 2000 s after solidification onset for different crystal radii: (a, b) 25 mm, (c, d) 60 mm, (e, f) 75 mm and (g, h) 95 mm).

In practice only A.G. Ostrogorsky et al. used Tiller’s ratio to describe the SHM growth experimental results. The best fit with theory was achieved by analyzing the Te distribution in seed grown small-diameter GaSb single crystals [7] and later [8] where a disc-shaped seed with a diameter equal to the crucible diameter was used by analogy with the AHP method. However to achieve a closer matching the authors had to use the diffusion coefficient D = 1 × 10-5 cm2/s for the calculations as per Eq. (1) which is several times smaller than the conventional value for tellurium [20, 21]. The authors justified their choice by the fact that the diffusion coefficient in the vicinity of the growing crystal within the boundary diffusion layer is far lower compared with that in the liquid bulk.

We however would like to offer a simpler explanation related with uncontrolled melt layer thickness change and hence growth rate in the experiments. This explanation can well be accepted taking into account that for the submerged heater design in question A.G. Ostrogorsky et al. never used a heating element or reduced the furnace temperature but only pulled the crucible to the cool zone within the preset temperature gradient. The latter technical solution was dictated by the fact that there is not thermocouple at the crucible bottom in the SHM method which would provide Tcool temperature readings and thus allow temperature control during crystal growth in accordance with Eq. (3). Therefore just like in many other methods there is no ground to assume that the growth rate (phase boundary velocity) V is equal to the pulling rate P. Although the crucible bottom temperature decreases, crystallization delays by several minutes or a longer time for the following reasons. First, there is a transient process of establishing during which the melt/crystal system has to be sufficiently cooled to a certain new condition. To start crystal growth, it is insufficient to reduce the temperature solely by lowering the crucible within the existing temperature gradient. The temperature should be somewhat lower in order to allow crystallization heat removal. Otherwise, as follows from Eq. (3), the Tcool temperature will decrease while the crystal height H will remain the same. Furthermore one can expect that the temperature gradient established before pulling will be closer to the temperature gradient in the melt since in the initial position [8] the melt zone was ~60 mm while the crystal was only 6 mm thick. In this case, as follows from Eq. (3), the melt layer h will grow with increasing distance between the heater bottom and the crucible bottom. This is illustrated by Table 3 demonstrating for the example of germanium crystallization that even if the furnace temperature gradient is the same as the one in the crystal before growth onset (20 K/cm) the growth will be accompanied by an increase in h to 10.87 mm after pulling for 1 h. Only setting a slightly higher temperature gradient (23.5 K/cm) will keep the layer thickness the same (h = 10.0, see Table 3) which should be ensured by choosing correct crystallization mode. If the furnace temperature gradient is close to the melt temperature gradient, then even with allowance for the crystallization heat the layer thickness will be 11.81 mm after pulling for 1 h (Table 3). An additional factor that is specific for submerged heater growth is the presence of a small but much overheated (by at least 30–50 K) melt portion moving from the zone above the heater to the zone under the heater immediately after the start of pulling. This also contributes to the increase in h. Therefore a general technical recommendation is to achieve the target puling rate gradually in order for the crystal to have sufficient time for cooling and thus to grow in accordance with the preset cooling rate.

Change in germanium melt layer thickness h as a function of furnace temperature gradient for crucible pulling to cool zone*.

No. Crystallization Stage Temperature Reduction under Condition: gradTc, K/cm Cooling Zone Gradient, K/cm T cool, °C h, mm
1 Before pulling Crystallization heat is constant 20.0 20.0 927.0 10.00
2 After lowering by 5 mm 20.0 20.0 917.0 10.00
3 20.0 10.0 922.0 11.43
4 Crystallization heat is released but ignored in cooling rate calculations 23.5 20.0 917.0 10.87
5 23.5 10.0 922.0 12.38
6 With account of crystallization heat release 23.5 23.5 913.5 10.00
7 23.5 10.0 920.2 11.81

Once the melt is finally cooled to the temperature required for the onset of growth the crystallization front velocity V starts increasing but reaches the preset pulling rate R nonmonotonically. Calculations of V based on preset Tcool and Thot temperatures in a CsI (Tl) single crystal growth experiment [22] showed that during the time of establishing (until the establishment of a new growth rate) the instantaneous growth rate was several times greater than the growth rates at the target growth rate change from 2 to 4 mm/h due to melt/crystal system cooling. This is also corroborated by our calculation results which will be analyzed below. Thus one can quite reasonably assume that the melt layer in the cited work [8] was not constant but quite probably started decreasing after a certain increase so the crystallization front almost approached the submerged heater bottom and followed it thereafter. Thus the growth rate finally became equal to the pulling rate after passing its peak value. Otherwise the experiment would be stopped due to heater freezing into the crystal. This seems to be the case in earlier work [23] which included unfinished InSb crystal growth experiments under microgravity conditions on board the ISS where a noticeable jump (by approx. 3–4 times) of the heat conductivity coefficient was observed during the solid to liquid state transition [24]. This leads to the formation of significantly different axial temperature gradients in the melt and in the crystal. Crystallization from a very thin layer [24] is suggested by the short time in which a steady-state Te distribution was achieved over the entire 15–20 mm crystal length which evidently testifies to a diffusion mass transport pattern. On the contrary, during growth of similar crystals on Earth [25] the axial temperature gradient in a similar furnace was insufficient for the crystallization front to follow the submerged heater. In all the experiments the melt seemed to occupy the entire ampoule in some time, reaching at least a 30–40 mm height. This melt portion most likely solidified rapidly during furnace cooling (after the end of the experiment) and was polycrystalline as the presented photos suggest. The single crystal regions were within 2–12 mm in length which is comparable with the size of the seed.

With the above assumptions about the variation patterns of the actual crystallization front velocity and the thickness of the melt layer from which the GaSb crystal grew, Curves 3 and 4 (Fig. 8) calculated in accordance with Eqs. (1) and (2) using the standard data on the tellurium diffusion coefficient D = 3 × 10-5 cm2/s are quite close to Curve 2 calculated for D = 1 × 10-5 cm2/s in the assumption that V = R. Thus Curves 3 and 4 (Fig. 8) describe earlier experimental data [8] adequately well. In comparison with Curve 2 they better describe the experimental data for the initial crystal part because Eq. (1) yields a convex curve unless the change of V in time is taken into account whereas the experimental curves have a concave bend, just like Curves 3 and 4. Thus, approximations of Eqs. (1) and (2) can be used in impurity distribution calculations for specific shapes of the crystallizing layer, provided however that the actual crystal growth rate and melt layer thickness patterns are taken into account.

Figure 8.

Calculation of Te distribution along GaSb axis using (13) Tiller’s ratio (Eq. (1)) and (4) Eq. (2): (1, 2) V = 5 mm/h, D = 3 × 10-5 and 1 × 10-5 cm2/s, respectively; (3) D = 3 × 10-5 cm2/s for V growth from 0 to 15 mm/h and further drop to 5 mm/h; (4) h first reaches 12 mm, then drops to 2 mm when growth rate exceeds pulling rate.

4. Growth of length homogeneous crystals

Even if the melt zone under the submerged heater in the AHP method is not intentionally doped, the impurity distribution over the silicon ingot height as shown in Fig. 9 (Curves 1) will differ from that for standard directional melt crystallization because the thin ~ 10–20 mm melt layer achieves a steady impurity concentration relatively rapidly, this occurring much faster for boron with a close to unity segregation ratio than for phosphorus. However if the melt zones under and above the submerged AHP heater are doped as described above, the boron and phosphorus concentrations will be constant over the entire growth process in accordance with Eq. (2) regardless of the melt layer thickness which can be chosen as large as 50–100 mm depending on the crucible size. Certainly these element concentrations will be at their respective levels depending on the initial B and P concentrations in silicon.

Figure 9.

(a, b) boron and (c, d) phosphorus distributions over silicon ingot length: (b, d) initial establishing stage: (1) without special melt zone doping (C1 = C2 = CS/k); (2) with additional doping of zone under AHP heater (C2 = CS, C1 = C2/k).

It can be seen from Fig. 10 that it is no problem to ensure constant boron concentration in the final portion of crystallizing silicon. Whichever crystallization rate in the 1 to 60 mm/h range is chosen the B concentration will change by 30 at.% worst. The situation is completely different for phosphorus distribution. To achieve at least the same result as for boron one should set a very low crystallization rate, ~1 mm/h. The process will take 50–100 which is almost one order of magnitude longer than the crystallization time of the main ingot part. Increasing the growth rate to 60 mm/h will lead to at least a threefold increase in the impurity concentration, to say nothing about possible structural quality problems of the material. Obviously to satisfy both the conditions, i.e., maintain the P concentration within the desired range and complete crystallization within a reasonable time, one should change the growth rate at a certain process stage. Figure 1 illustrates this possibility: the total crystallization time will be ~15 h if the growth rate is abruptly reduced almost to zero after the onset of crystallization. This can be achieved by changing the cooling rate in time, though this is not such a simple task and cannot be solved simply by replicating the same time function Rcool(t) as that for V. This is illustrated by Fig. 11a which suggests that the growth rate is a nonlinear function of the cooling rate.

Figure 10.

(a) boron and (b) phosphorus distributions in the final silicon portion for (a) 1, (b) 10 and (c) 60 mm/h crystallization rate.

Figure 11.

Typical phase boundary rate V change pattern for (a) germanium melt layer solidification and (b) phosphorus distribution for silicon crystallization at variable rate V = V0(1 – 0,5x0.3)

5. Summary

The ANSYS Fluent software package along with the in-house mathematical tools were used for the study of 2D segregation in thin melt layers. We chose antimony doped germanium as simulation object due to the low antimony segregation ratio of 0.003. Natural convection during AHP semiconductor growth can be efficiently suppressed and its effect on diffusion can be made negligible for some flat layer thicknesses, i.e. the so-called diffusion-controlled crystal growth mode can be established. This mode allows describing segregation using simple expressions which yield good results if the changes in the actual growth rate and melt layer thickness during crystallization are taken into account. Additional melt doping under the AHP heater provides for a constant impurity concentration in silicon during crucible pulling, and upon the solidification of the residual layer this can be achieved by varying the growth rate (from high to low) through changing the cooling rate in time. By and large the required cooling rate vs time pattern can be determined by solving the inverse task; this can be the objective of the next work stage.


  • 3 Golyshev V. D., Gonik M. A. A temperature field investigation in case of crystal growth from the melt with a plane interface on exact determination thermal conditions. Crystal Properties and Preparation. 1991;36–38: 623–630.
  • 4 Gonik M. A. Directional crystallization of multicrystalline silicon in a weak melt convection and gas exchange. Izvestiya Vysshikh Uchebnykh Zavedenii. Materialy Elektronnoi Tekhniki = Materials of Electronics Engineering. 2015; 18(2): 95–102. (In Russ.).
  • 5 Filonov K. N., Kurlov V. N., Klassen N. V., Kudrenko E. A., Shteinman E. A. Peculiarities of nanostructured silicon carbide films and coatings obtained by novel technique. Izvestiya Rossiiskoi Akademii Nauk. Seriya Fizicheskaya = Bulletin of the Russian Academy of Sciences: Physics. 2009; 73(10): 1457–1459. (In Russ.)
  • 6 Burago N. G., Golyshev V. D., Gonik M. A., Polezhaev V. I., Tsvetovsky V. B. The nature of forced and natural convection and its influence on the distribution of impurities in a crystal during growth using the OTF 1a method. Trudy III Mezhdunarodnoi konferentsii «Kristally: rost, svoistva, real’naya struktura i primenenie» = Proceedings of the III International Conference «Crystals: growth, properties, real structure and application». Alexandrov, 1997; 1: 239–259. (In Russ.)
  • 10 Ostrogorsky A. G., Mosel F., Schmidt M. T. Diffusion-controlled distribution of solute in Sn-1% Bi specimens solidified by the submerged heater method. J. Crystal Growth. 1991; 110(4): 950–954.
  • 12 Golyshev V. D., Gonik M. A., Tsvetovsky V. B., Frjazinov I. V., Marchenko M. P. Proceedings of the 3rd International Conference on Single Crystal Growth, Strengthen Problems, Heat and Mass Transfer. Obninsk (Russia), 2000: 125–134.
  • 14 Voller V. R., Prakash C. A fixed grid numerical modelling methodology for convection-diffusion mushy region phase-change problems. Int. J. Heat Mass Transfer. 1987; 30(8): 1709–1719.
  • 15 Smirnova O. V., Kalaev V. V., Makarov Yu. N., Abrosimov N. V., Riemann H., Kurlov V. N. Three-dimensional unsteady modeling analysis of silicon transport in melt during Cz growth of Ge1-xSix bulk crystals. J. Crystal Growth. 2007; 303(1): 141–145.
  • 16 Golyshev V. D., Gonik M. A., Tsvetovsky V. B. Study of thermal conductivity close to the melting point. High Temperatures-High Pressures. 2003/2004; 35–36(2): 139–148.
  • 17 Sato Y., Nishizuka T., Tachikawa T., Hoshi M., T. Yamamura, Y. Waseda. Viscosity and density of molten germanium. High Temperatures-High Pressures. 2000; 32: 253–260.
  • 18 Ratnieks G., Muižnieks A., Mühlbauer A. Modelling of phase boundaries for large industrial FZ silicon crystal growth with the needle-eye technique. J. Crystal Growth. 2003; 255(3–4): 227–240.
  • 20 Raffy C., Duffar T. Internal Report, CEA-Grenoble (France), SES No. 15/95, 1995.
  • 21 Müller G. Convection and Inhomogeneities in Crystal Growth from the Melt. In: Crystals: Growth, Properties and Applications; 12. Berlin; Heidelberg: Springer, 1988: 1–136.
  • 23 Ostrogorsky A. G., Marin C., Churilov A., Volz M. P., Bonner W. A., Duffar T. Reproducible Te-doped InSb experiments in Microgravity Science Glovebox at the International Space Station. J. Crystal Growth. 2008; 310(2): 364–971.
  • 24 Nakamura S., Hibiya T., Yamamoto F. Thermal conductivity of GaSb and InSb in solid and liquid states. J. Appl. Phys. 1990; 68(10): 5125–5127.