Document Type : Full Lenght Research Article

**Authors**

Physics Department, Bu-Ali Sina University

**Abstract**

For a Czochralski growth of Ge crystal, thermal fields have been analysed numerically using the three-dimensional finite volume method (FLUENT package). The arrangement used in a real Czochralski crystal growth lab included a graphite crucible, heat shield, heating device, thermal insulation and chamber including two gas outlets. We have considered two cases for calculations, which are configuration containing (a) only gas and (b) melt and gas, related to initial stages of the growth process (seeding process). It has been assumed that the growth system is in steady state, fluids are incompressible Newtonian fluids and the flow is laminar. It was shown that the thermal field in the growth setup is completely three-dimensional. Especially, the temperature field at the melt free surface has not a uniform radial distribution due to the three-dimensional orientation of Argon flow above it.

**Keywords**

**Main Subjects**

1. Introduction

The Czochralski crystal growth technique is one of the most essential methods for producing large size germanium single crystals. Because of the multidisciplinary nature of the crystal growth technology, understanding of the growth phenomena as well as the interaction between the many factors influencing the final crystal quality, not only is not so easy but also is a very complicated task. However, because of high temperature, opaque materials and high purity requirements, improvement of any extensive experiments in the crystal growth environment is proved to be time consuming and extremely expensive. Therefore, numerical modelling and computer simulation has become an essential and indispensable tool for understanding, developing and optimizing the crystal growth processes and related equipment [1-41], that provides a meaningful and cost-effective alternative tool to experimental trials. The ultimate goal of numerical modelling in crystal growth technology is to control and optimize the growth process.

The transport processes (heat and mass transfer including melt and gas convection) problem in Czochralski crystal growth technique represents one of the biggest challenges of crystal growth modelling [6,9,24,34]. Crowly [42] developed a mathematical model of the heat transfer in the holm region in the germanium Czochralski furnace including two moving boundaries, the phase change surface and the air-liquid meniscus applying the enthalpy method. Dupret etal. [43] implemented a simulation of the temperature field and the position of the crystallization front during growth of germanium crystal using radiative heat transfer. Bogaert and Dupret [44,45] computed the time-dependent Czochralski growth of germanium crystal including heat exchange by conduction and diffuse grey radiation, in and between all the setup constituents, together with all the transient effects induced by the growth of crystal without melt convection. Bykova [46] simulated a two-dimensional time-dependent Ge crystal growth by the AHP method under microgravity conditions and showed that axial microaccelerations will have no impact on the forced melt flow along the crystallization front. Abbasoglu [47] performed a transient 3D numerical simulation to examine the roll of crystal and crucible rotations on the fluid flow and the radial segregation of silicon during the growth of Ge_{x} Si_{1−x }crystals by the Czochralski technique using microgravity conditions. Recently, Honarmandnia etal. [48,49] have carried out a 2D global simulation of an RF Czochralski apparatus for different stages of germanium crystal growth in order to analyse of thermal field, the convexity of the crystal–melt interface and thermo-elastic stresses of the grown crystal. Besides, an asymmetrical configuration leads to asymmetric and 3D thermal fields in the growth setup, which in turn, changes the crystal quality, Figure 1. For this reason, a 3D calculation is required and has to be performed.

It is usual to provide a heat shield assembly disposed above the molten source material (Si and Ge) and surrounding the ingot as it is pulled upward from the melt to shield the grown crystal against heat radiated from the crucible and the heater surrounding the crucible. Heat shield assemblies are typically constructed of graphite. Because of a relatively high emissivity of graphite, a conventional heat shield has a high ability to emit radiant energy from its surface, and so it radiates a substantial amount of heat toward the grown crystal, thereby inhibiting cooling of it. Consequently, the crystal pulling rate from the melt can be increased in crystal puller setups. On the other hand, the heat shield is refereed sometimes to a gas flow guide shield. Because it guides the inert inlet gas (usually Argon) toward the ingot to aid the cooling of the single crystal (reducing the possibility of microdefect nucleation). Another important effect of the guided gas flow is that the melt free surface is blown by it in such a manner that the gas flow takes the impurity away from the free surface regularly. Therefore, the obtained grown crystal will contain a reduced concentration of impurities and so a crystal with improved quality is achieved.

The goal of this work is to apply an appropriate 3D numerical approach of a resistance heated CZ furnace for a Ge crystal before the seeding process (i.e., including only melt and gas), and to analyze the structure of fluid flow, temperature field and global heat transfer using finite volume method (CFD FLUENT Package [50]). Our attention is specially focused on the region where the crystal is grown. It is worth to note that the obtained numerical results are interesting and quite important for the growers, which present detail information about the nonsymmetrical thermal field in the growth setup. During the growth process, from seeding to cooling crystal, a real-time control and monitoring of the thermal field is crucial task because it directly effects on the shape and quality of the grown crystal. For example, a nonsymmetrical thermal field can change the cylindrical shape (circular cross section) of the crystal to a non-cylindrical shape [51-53]. Therefore, improvement and reduction of the non-cylindrical symmetric conditions of the growth setup, and precisely controlled the growth procedure are possible using the simulation results presented here.

**2. Model description**

*2.1. Governing equations*

Figure 1 represents the schematic diagram of the CZ growth, which is used in this study. Our model covers the following key assumptions: (1) The system is not time dependent (i.e., steady state). (2) Melt and gas are incompressible Newtonian fluids keeping the Boussinesq approximation. (3) The flow is laminar. (4) Viscous dissipation is not included. (5) The surface of melt is flat, i.e., there is not any meniscus at the crucible wall. (5) The thermo-physical properties of the fluids are constant except for the density difference in the buoyancy force term. (6) Fixed temperature is set on the chamber walls.

Fig 1. Schematic model of the resisting heated Czochralski system for Ge growth.

Thus, the basic equations are;

(a) Fluid flow in the melt and gas:

(1)

which is the Navier-Stokes equation with Boussinesq approximation.

(b) Continuity equation:

(2)

c) Energy equation:

(convection and conduction in the melt and gas) (3)

(conduction in the solid parts) (4)

where the subscripts f and s denote fluids (gas and melt), solids (crucible, heat shield, insulation and chamber), respectively. The boundary conditions are:

a) at the solid-fluid interfaces:

(5)

b) at the melt-gas interface:

(6)

(7)

(8)

where the subscripts m and g denote melt and gas, respectively. The Equation (7) represents the thermal Marangoni phenomena at this interface.

c) at the outer surfaces of the chamber (the water cooled wall):

(9)

d) at the solid-gas interfaces:

(10)

which denotes the gas cooling combined with surface to surface heat radiation exchange from those surfaces.

e) at the gas inlet:

(11)

(12)

f) at the gas outlet:

(13)

In the above equations, = (u, v, w) is the fluid velocity vector in the cylindrical coordinate system (*r*,*φ* *z*), p the pressure, T the temperature, the z-directional unit vector, g the acceleration due to gravity, β the thermal expansion coefficient, µ the dynamic viscosity, α the thermal diffusivity, γ the surface tension, ε the emissivity, k the thermal conductivity, ρ the density, the unit normal vector, the unit tangential vector along the meridional direction and σ the Stefan-Boltzmann constant.

*2.2. Numerical method*

In the next step of the calculation procedure, the current problem of great interest is to obtain a very high-order accurate and efficient numerical solution for the system of governing equations in the computation domain. It is as follows:

1. The governing equations were solved employing the finite-volume based CFD FLUENT Package with three-dimensional double precision second-order discretization. The finite volume method is a common way used in several CFD codes and has some advantages in memory usage and solution speed, especially for large problems (like our problem).

2. Since the domain of interest has a complex geometrical shape, an unstructured grid of more than 5×10^{5} finite control volume was employed, Figure 3.

3. The convergence is checked by monitoring residuals of momentum, continuity and energy equations which were set at 10^{-5}, 10^{-5} and 10^{-6}, respectively.

*2.3. The calculation conditions*

All parts of the considered CZ setup have cylindrical symmetry except for gas ducts and outlets. Placing of two gas ducts in the bottom insulation and the same outlets in the chamber breaks the 2D rotational symmetry of the system. This means that the CZ setup has still a symmetry axis (z) and two perpendicular planes of symmetry (zx and zy), Figure 2. The gas ducts and outlets are located in the zx plane.

The temperature at the central point of the melt-gas interface (i.e., the position of the seed crystal) is 1210 K, i.e., 20 K above the melting point of Ge (real condition before the seeding). It is not a boundary condition but it must be a calculation result by adaption of the heater power. The inlet gas flow rate is set to be 1 lit/min. The thermophysical properties used for our simulations are selected from the FLUENT materials database.

In order to analyse and obtain detailed information about the effects of the asymmetrical geometry and the gas flow on the thermal conditions during the Ge growth, we have considered two computational cases, the configuration contains (a) only gas (i.e., without any melt and crystal), and (b) gas and melt, confirming to the real condition just before the Ge seed contacts the melt. We will discuss the flow and temperature fields in these two casesspecifically.

Fig 2. Two perpendicular planes of symmetry zx and zy in the growth setup.

Fig 3. The applied unstructured grid used for the calculation.

Fig 4. Temperature field (right hand side) and the vector velocity (left hand side) in the symmetry plane zx for Case a - configuration containing only gas.

**3. Results and discussion**

*3.1. Case a: Configuration contains only gas*

Figures 4 and 5 show the temperature field (right hand side) and flow arrows (left hand side) in the setup for both perpendicular planes of symmetry zx and zy. In both planes, a counter clockwise gas vortex caused by free convection (buoyancy driven) with average speed exists and occupies the upper entire volume of the configuration (i.e., independent of the asymmetrical conditions). This flow carries the heat from the upper side of the heat shield to the middle space (location of the seed holder) and so prevents the inlet gas flow in that part. For this reason, the inlet gas flow is completely deflected from the central part toward the chamber side wall. Then it is directed toward the central part along the upper side of the heat shield, downward into the crucible, passing the area between the heat shield and the crucible, and then again downward to the space below the crucible. In that part of the system, the gas flow is guiding into two ducts placed in the bottom insulation (zx plane) and speeding it up (average speed ), and then exits the setup via the chamber outlets with average speed . For this reason, the behaviour of gas flow is completely three-dimensional and non-cylindrical symmetric in that region.

Fig 5. Temperature field (right hand side) and the flow arrows (left hand side) in the symmetry plane zy for Case a - configuration containing only gas.

Fig 6. A three-dimensional view of the temperature field in the symmetry planes zx and zy for Case a - configuration containing only gas.

The temperature gradient in the setup is considerably affected by the orientation of gas convection, and so it is completely three-dimensional too, Figure 6. In the growth setup, the role of gas flow is quite important as other mechanisms of the heat transport phenomena (such as conduction and radiation). Therefore, the isotherms have a curvature commensurate to the structure and intensity of gas flow. The isotherms are pushed along the gas flow (proportional to its velocity) in the upper part of the system and deflected toward the bottom ducts and outlets. The temperature maximum in the setup is and is located at the upper part of the heating elements.

Fig 7. Temperature field (right hand side) and the flow arrows (left hand side) in the symmetry plane zx for Case b - configuration containing melt and gas.

Fig 8. Temperature field (right hand side) and the flow arrows (left hand side) in the symmetry plane zy for Case b - configuration containing melt and gas.

Fig 9. Temperature field (right hand side) and the flow arrows (left hand side) in the symmetry planes zx and zy in the melt for Case b.

Fig 10. A three-dimensional view of the temperature field in the symmetry planes zx and zy for Case b - configuration containing melt and gas.

** **

*3.2 Case b: Configuration contains melt and gas*

Figures 7 and 8 represent the temperature field (right hand side) and flow arrows (left hand side) in the growth system for both perpendicular planes of symmetry (zx and zy). The orientation of argon flow is completely similar to Case a, which is totally three-dimensional as has been already described. The gas average speed is which is a little higher than Case a.

In the molten material, a counterclockwise eddy that is produced by the free convection overlapped with the thermal Marangoni flow happens and occupies its entire volume with an average speed , Figure 9. This vortex has expanded from the crucible wall to the centerline. The melt flow is started along the hot crucible sidewall (buoyant flow) and the melt-gas interface where the temperature gradient via the surface tension produces the thermal Marangoni convection [40,47].

The temperature field in the growth setup is distinctly influenced by the gas flow similar to Case a. The influence of gas flow is more essential than conduction and radiation in the gas part. Therefore, the isotherm shapes have a deviation commensurate to the gas convection approach and flow arrangement. The gas isotherms are moved downward along the inner surface of the heat shield, upward adjacent to inner crucible sides, deflected downward in the space between the crucible and heater and deformed toward the gas ducts and outlets, Figure 10. The temperature maximum in the furnace is and is located at the lower part of the heating zone contrary to Case a. In the Ge melt which its location is in the lower portion of the crucible sidewall and close to its bottom.

Figure 11 illustrates the temperature variation of the melt surface and its sensitivity to the gas flow above the melt surface. It shows a non-cylindrical symmetric structure, which arises from the different gas flow in the zx and zy symmetry planes. In the zx plane, the argon flow rate more accelerates and becomes stronger than zy plane and so it can cool the melt surface more efficiently in that direction. Therefore, the temperature maximum of the melt surface is located in the zy symmetry plane with more compared to the zx plane. It should be noted here that the non-symmetrical conditions of thermal field in the growth setup is quite similar in both cases considered here which is arises directly from different argon flow in the zx and zy symmetry planes.

It has been mentioned that any non-symmetrical condition for setup thermal field is not desirable and can lead to crucial problems related to the crystal quality. Consequently, we have to prevent or reduce these non-symmetrical thermal conditions of the growth furnace by special considerations during the system design, for example, using four symmetrical gas ducts and outlets in the presented growth setup.

**4. Conclusions and outlook**

We have presented and exhibited the results of a three-dimensional global numerical calculation of the flow structure and temperature field for a non-cylindrical Ge crystal puller containing only gas (Case a), and melt and gas (Case b). From these simulations, we can conclude:

1- The presented numerical solution is only the first result of the considered non-cylindrical Ge Czochralski configuration including seed and crystal. Therefore, the obtained results are useful for the seeding process because the detailed information about the thermal field at the melt surface is a crucial factor for a successful seeding process (before, during after touching the seed to the melt).

Fig 11. Temperature variation of the melt free surface, which shows a non-rotationally symmetric distribution.

2- Although the temperature gradient is symmetrical around the center of melt surface but it is completely non-cylindrical symmetric close to the crucible wall which arises from the three-dimensional orientation of argon flow above it. The temperature difference along the y-direction is more than x-direction, Figure 11. Therefore, this condition can influence on the control of crystal diameter and its uniformity. For this reason, further calculations including the seed and grown crystal have to be performed.

] N. Kobayashi, Computational simulation of the melt flow during Czochralski growth, J. Crystal Growth, 43, 357-363 (1978).

[2] N. Kobayashi, Hydrodynamics in Czochralski growth-computer analysis and experiments, J. Crystal Growth, 52, 425-434 (1981).

[3] T. Tsukada, N. Imaishi and M. Hotawa, Theoretical study of the flow and temperature fields in CZ single crystal growth, J. Chem. Eng. Jpn, 21, 184-191 (1988).

[4] J.J. Derby, L. J. Atherton and P. M. Gresho, An integrated process model for the growth of oxide crystals by the Czochralski method, J. Crystal Growth, 97, 792-826 (1989).

[5] J.J. Derby and Q. Xiao, The role of internal radiation and melt convection in Czochralski oxide growth: deep interfaces, interface inversion, and spiraling, J. Crystal Growth, 113, 575-586 (1991).

[6] J.J. Derby and Q. Xiao, Heat transfer and interface inversion during the Czochralski growth of yttrium aluminium garnet and gadolinium gallium garnet, J. Crystal Growth, 139, 147-157 (1994).

[7] T. Tsukada, N. Imaishi and M. Hozawa, Global analysis of heat transfer in CZ crystal growth of oxide, J. Chem. Eng. Jpn., 27, 25-31 (1994).

[8] J. Järvinen, R. Nieminen and T. Tiihonen, Time-dependent simulation of Czochralski silicon crystal growth, J. Cryst Growth, 180, 468-476 (1997).

[9] E. Dornberger, E. Tomzig, A. Seidl, S. Schmitt, H.J. Leister, C. Schmitt and G. Muller, Thermal simulation of the Czochralski silicon growth process by three different models and comparison with experimental results, J. Crystal Growth, 180, 461-467 (1997).

[10] G. Muller, A. Muhe, R. Backofen, E. Tomzig and W. Ammon, Study of oxygen transport in Czochralski growth of silicon, Microelectronic Engineering, 45, 135-147 (1999).

[11] C.J. Jing, N. Imaishi, S. Yasuhiro and Y. Miyazawa, Three dimensional numerical simulation of spoke pattern in oxide melt, J. Crystal Growth, 200, 204-212 (1999).

[12] C.J. Jing, N. Imaishi, T. Sato and Y. Miyazawa, Three-dimensional numerical simulation of oxide melt flow in Czochralski configuration, J. Crystal Growth, 216, 372–388 (2000).

[13] A. Lipchin and R.A. Brown, Hybrid finite-volume/finite-element simulation of heat transfer and melt turbulence in Czochralski crystal growth of silicon, J. Cryst Growth, 216, 192-203 (2000).

[14] Z. Galazka and H. Wilke, Influence of the Marangoni convection on the flow pattern in the melt during growth of Y_{3} Al_{5}O_{12} single crystals by the Czochralski method, J. Crystal Growth, 216, 389–398 (2000).

[15] K. Takano, Y. Shiraishi, J. Matsubara, T. Iida, N. Takase, N. Machida, M. Kuramoto and H. Yamagishi, Global simulation of the CZ silicon crystal growth up to 400 mm in diameter, J. Crystal Growth, 229, 26-30 (2001).

[16] S. Enger, O. Grabner, G. Muller, M. Breuer and F. Durst, Comparison of measurements and numerical simulations of melt convection in Czochralski crystal growth of silicon, J. Crystal Growth, 230, 135-142 (2001).

[17] V.V. Kalaev, I.Y. Evstratov and Y.N. Makarov, Gas flow effect on global heat transport and melt convection in Czochralski silicon growth, J. Crystal Growth, 249, 87-99 (2003).

[18] V.V. Kalaev, D.P. Lukanin and V.A. Zabelin, Prediction of bulk defects in CZ Si crystals using 3D unsteady calculations of melt convection, Materials Science in Semiconductor Processing, 5, 369-378 (2003).

[19] N.G. Ivanov, A.B. Korsakov and E.M. Smirnov, Analysis of magnetic field effect on 3D melt flow in CZ Si growth, J. Crystal Growth, 250, 163-174 (2003).

[20] A.J. Nowak, R. Biaecki, A. Fic and G. Wecel, Analysis of fluid flow and energy transport in Czochralski’s process, Computers and Fluids, 32, 85-95 (2003).

[21] T. Wetzel, J. Virbulis, A. Muiznieks, W. von Ammon, E. Tomzig, G. Raming and M. Weber, Prediction of the growth interface shape in industrial 300 mm CZ Si crystal growth, J. Crystal Growth, 266, 34-39 (2004).

[22] L. You Rong, R. Dengfang and P. Lan, Global simulation of silicon crystal Czochralski growth I.Characteristics of heat transfer and fluid flow, Chinese Journal of Materials Research, 18, 212-218 (2004).

[23] L. Liu, T. Kitashima and K. Kakimoto, Global analysis of effects of magnetic field configuration on melt–crystal interface shape and melt flow in CZ-Si crystal growth, J. Crystal Growth, 275, e2135-e2139 (2005).

[24] S. Shufang, Z. Qingli and S. Jing, Research Progress in Numerical Simulation for Crystal Growth by Czochralski Method, Journal of Synthetic Crystals, 34, 687-693 (2005).

[25] J. Banerjee and K. Muralidhar, Simulation of transport processes during Czochralski growth of YAG crystals, J. Crystal Growth, 286, 350-364 (2006).

[26] V.V. Kalaev, Combined effect of DC magnetic field and free surface stresses on the melt flow and crystallization front formation during 400 mm diameter Si Cz crystal growth, J. Crystal Growth, 303, 203-210 (2007).

[27] M.H. Tavakoli and H. Wilke, Numerical study of heat transport and fluid flow of melt and gas during the seeding process of sapphire Czochralski crystal growth, Crystal Growth Des., 7, 644-651 (2007).

[28] M.H. Tavakoli and H. Wilke, Numerical investigation of heat transport and fluid flow during the seeding process of oxide Czochralski crystal growth Part 1: non-rotating seed, Cryst. Res. Technol, 42, 544-557 (2007).

[29] M.H. Tavakoli and H. Wilke, Numerical investigation of heat transport and fluid flow during the seeding process of oxide Czochralski crystal growth Part 2: rotating seed, Cryst. Res. Technol, 42, 688-698 (2007).

[30] M.H. Tavakoli, S. Omid and E. Mohammadi-Manesh, Influence of active afterheater on the fluid dynamics and heat transfer during Czochralski growth of oxide single crystals, CrystEngComm, 13, 5088-5093 (2011).

[31] A. Raufeisen, S. Jana, M. Breuer, T., Botsch and F. Durst, 3D computation of oxygen transport in Czochralski crystal growth of silicon considering evaporation, J. Crystal Growth, 303, 146-149 (2007).

[32] P. Gunjal, M. Kulkarni, P. Ramachandran, Melt flow simulations of Czochralski crystal growth process of silicon for large crystals, ECS Transactions, 3, 41-52 (2006).

[33] M.H. Tavakoli, Numerical study of heat transport and fluid flow during different stages of sapphire Czochralski crystal growth, J. Cryst Growth, 310, 3107-3112 (2008).

[34] O. Smirnova, N. Durnev, K. Shandrakova, E. Mizitov and V. Soklakov, Optimization of furnace design and growth parameters for Si Cz growth, using numerical simulation, J. Crystal Growth, 310, 2185-2191 (2008).

[35] I.C. Avetissov, E.A. Sukhanova, A.P. Sadovskii, V.A. Kostikov and E.V. Zharikov, Experimental and numerical modeling of Czochralski crystal growth under axial vibrational control of the melt, J. Crystal Growth, 312, 1429-1433 (2010).

[36] F. Mokhtari, A. Bouabdallah, A. Merah, M. Zizi, S. Hanchi and A. lemany, Three-dimensional study of the pressure field and advantages of hemispherical crucible in silicon Czochralski crystal growth, Cryst. Res. Technol., 45, 573-582 (2010).

[37] X. Cen, J. Zhan and Y.S. Li, Large eddy simulation of Marangoni convection in Czochralski crystal growth, Cryst. Res. Technol., 46, 14-22 (2011).

[38] W. Miller, C. Frank-Rotsch and P. Rudolph, Numerical studies of flow patterns during Czochralski growth of square-shaped Si crystals, J. Crystal Growth, 318, 244-248 (2011).

[39] O. Asadi Noghabi, M. M'Hamdi and M. Jomaa, Effect of crystal and crucible rotations on the interface shape of Czochralski grown silicon single crystals, J. Crystal Growth, 318, 173-177, (2011).

[40] J. Chaohua, Effect of crystal and crucible rotations on the interface shape of Czochralski grown silicon single crystals, J. Semiconductors, 34,063005-6, (2013).

[41] M. Kirpo, Global simulation of the Czochralski silicon crystal growth in ANSYS FLUENT, J. Crystal Growth, 371, 60-69 (2013).

[42] A.B. Crowley, Czochralski Crystal Pulling, IMA J. Appl. Math., 30, 173-189 (1983).

[43] F. Dupret, P. Nicodème, Y. Ryckmans, P. Wouters and M. J. Crochet, Global modelling of heat transfer in crystal growth furnaces**, **Int. J. Heat Mass Transfer, 33, 1849**-**1871** (**1990**)**.

[44] N. Van den Bogaert and F. Dupret, J. Cryst. Growth, Dynamic global simulation of the Czochralski process, I. Principles of the method, 171, 65–76 (1997).

[45] N. Van den Bogaert and F. Dupret, J. Cryst. Growth, Dynamic global simulation of the Czochralski process, II. Analysis of the growth of a germanium crystal, 171, 67–93 (1997).

[46] S.V. Bykova, V.D. Golyshev, M.A. Gonik, V.B. Tsvetovsky, I.V. Frjazinov and M.P. Marchenko, J. Cryst. Growth, Features of mass transfer for the laminar melt flow along the interface, 237–239, 1886-1891 (2002).

[47] S. Abbasoglu, Three-dimensional numerical simulation of crystal and crucible rotations during Czochralski growth of Ge_{x}Si_{1−x} single crystals, Turk. J. Elec. Eng. & Comp. Sci., 20, 1343-1358 (2012).

[48] M. Honarmandnia, M.H. Tavakolia and H. Sadeghi, CrystEngComm, Global simulation of an RF Czochralski furnace during different stages of germanium single crystal growth, 18, 3942–3948 (2016).

[49] M. Honarmandnia, M.H. Tavakolia and H. Sadeghi, CrystEngComm, Global simulation of an RF Czochralski furnace during different stages of germanium single crystal growth, part II: to investigate the effect of the crucible's relative position against the RF coil on the isotherms, flow fields and thermo-elastic stresses, 19, 576–583 (2017).

[50] www.ansys.com/Products/Fluids/ANSYS-Fluent

[51] A. Raufeisen, M. Breuer, T. Botsch and A. Delgado, J. Cryst. Growth, Transient 3D simulation of Czochralski crystal growth considering diameter variations, 311 (2009) 695–697.

[52] D.T.J. Hurle, Crystal Pulling from the Melt, Springer-Verlag (1993).

[53] D.T.J. Hurle, Handbook of Crystal Growth. Vol 2, Elsevier (1994).

**Receive Date:**14 February 2016**Revise Date:**30 November 2017**Accept Date:**20 December 2017**First Publish Date:**01 May 2018