Effect of Internal Radiation on Process Parameters in the Global Simulation of Growing Large-Size Bulk β-Ga2O3 Single Crystals with the Czochralski Method

: As a crystal grows, the temperature distribution of the crystal and melt will change. It is necessary to study the dynamic process of single-crystal growth. Due to the relatively low crystallization rates used in the industrial Czochralski growth system, a steady state is used to compute the temperature distribution and melt flow. A two-dimensional axisymmetric model of the whole Czochralski furnace was established. The dynamic growth process of large-size bulk β-Ga 2 O 3 single crystal using the Czochralski method has been numerically analyzed with the parameter sweep method. In this paper, two cases of internal radiation and no internal radiation were compared to study the effect of radiation on the process parameters. The temperature distribution of the furnace, the temperature field, and the flow field of the melt was calculated. The temperature, the temperature gradient of the crystal, the temperature at the bottom of the crucible, and the heater power were studied for the crystals grown in the two cases of radiation. The results obtained in this study clearly show that the loss calculated by including the internal radiation is higher compared to that including the surface radiation.


Introduction
Gallium oxide has five polymorphisms, and other phases transform into the β phase at high temperatures. β-Ga2O3 is the most stable crystal structure, with a melting point of 1800 °C under atmospheric pressure [1,2]. Gallium oxide has a high electrical breakdown field and is an ideal material for the manufacture of high-power diodes. Gallium oxide can be used to manufacture various light-emitting diodes (LEDs) and laser diodes (LDs), new fluorescent materials, and gas sensors based on β-Ga2O3, as well as solar-blind UV detectors and high-temperature, high-frequency, and high-power electronics [3]. Gallium oxide has a good performance at high frequencies and high power and is an ideal material for power electronic equipment [4]. Gallium oxide has a wider band-gap than SiC and GaN, and its unique properties mean that it is used in many devices, including Schottky barrier diodes (SBDs), inverters equipment, circuits used in high-temperature, high-humidity and high-radiation environments, and field-effect transistors (FET) [5]. The application of β-Ga2O3 is expanding rapidly, and cost reduction for this method has become essential. One way to reduce cost is to increase the crystal diameter. However, enlarging the crystal diameter is quite difficult. The company Tamura has successfully prepared two-inch high-quality β-Ga2O3 single crystals using the edge-defined film-fed growth method [6]. Hoshikawa et al. used the vertical Bridgman (VB) method to successfully realize one-inch Ga2O3 single crystals. The crystal quality was good, without crucible pollution [7]. Galazka et al. [8] used the Czochralski method to grow β-Ga2O3 single crystals that were 2-3 inches in diameter. The melting temperature of gallium oxide is very high, and some materials are translucent to the near-infrared (NIR) wavelength range. In these materials, radiative heat transfer is more dominant than thermal transfer. As most oxide crystals have very large radiation losses, it is important to study radiation during growth [9]. The absorption coefficients of doped gallium oxide are different (4-100 cm −1 ). Crystals doped with Mg, Ce, and Al can be grown using the Czochralski method with a good quality, while doping with Cr and Si will lead to the spiral growth of crystals [10][11][12][13].
At present, the mechanism of spiral growth is not clear, and there are many explanations. The authors of [14,15] carried out three-dimensional time-dependent simulations and found that asymmetric temperature fields would cause periodic temperature oscillations. Spiral growth may be affected by asymmetric temperature distribution growth. Uecker [16] insisted that a reduction in radiant heat dissipation will flatten the radial temperature gradient in the melt near the interface, and increases the possibility of breaking symmetry. Enhancing the Marangoni flow also alleviates the formation of a spiral. Moreover, an increase in the radial temperature gradient (or the convexity of the growth interface) can inhibit spiral growth [17]. In addition, impurities and atmosphere affect optical absorption in the crystal and cylindrical growth stability [18,19]. The heat radiation from the iridium crucible reduces the radiant heat transfer through the crystal, resulting in spiral growth [20]. Some authors believed that the surface tension of the melt being reduced by impurities is the cause of spiral growth [21,22]. Schwabe disagreed with this view, believing that the two stages of spiral growth are caused by different reasons. The first stage may be caused by some external asymmetric thermal field, with the second stage of spiral growth being a kind of self-organization [23]. Galazka et al insisted β-Ga2O3 crystals have a tendency to a spiral formation due to free carrier absorption in the NIR wavelength range [10].
The Czochralski method was used to grow large-diameter single crystals of Ge and Si and some oxide crystals [24][25][26]. A small seed was dipped into the slightly overheated melt and pulled upwards from the melt. During the growth process, the crystal heating power, pulling rate, and rotation rate were carefully adjusted to control the shape of the crystal, especially the diameter. The automatic diameter control method for this crystal is called the crucible weighing method and was proposed by Kyle and Zydzik [27,28]. The crystal growth furnace for growing oxide crystals using the Czochralski method includes a heater, a puller, a pulling shaft connected to a weighing cell, a pulling/rotating motor, and a PID (proportional, integral, differential) controller computer. Miller [29] built a 2D global model for a crystal growth furnace to analyze the heat transport simulation and proceeded with 3D heat transport and fluid flow analysis on the basis of 2D simulation. Tang [30] established a 2-dimensional global model, analyzed the flow field and temperature field by steady-state simulation, and calculated the growth of the phase interface of different sizes by a transient method. In previous simulations, the steady-state calculations for a certain crystal length were studied and the effects of radiation in the medium were not taken into account. In this paper, dynamic modeling was used to analyze the effects of radiation on the process parameters during the whole growth process.

Geometric Model
Czochralski is one of the melt growth methods for large-size and high-quality β-Ga2O3 crystal growth, which has a low cost and a high crystal growth speed [31]. The geometry of Czochralski method is depicted in Figure 1. The growth furnace designed for β-Ga2O3 includes crucible support rods, crucible trays, heaters, shields, thermocouples, zirconia insulation ceramics, seed rods, and an iridium crucible and PID controller for automatic crystal diameter control. Zirconia ceramics are very good conductors at high temperatures, so in this model, resistance heating is used instead of electromagnetic in-duction heating. Gallium oxide melt is highly stable up to about 1200 °C and will decompose at high temperatures. Therefore, the Czochralski growth process requires a mixed atmosphere containing oxygen to inhibit the thermal decomposition of gallium oxide [32]. The growth chamber is the core of the growth furnace and consists of a metal crucible containing melt and insulated material surrounding the crucible. Crystals are pulled out from the crucible and an observation window (not shown) is on the side. After the raw material is melted, the directional seeds attached to the lifting shaft are dipped in the melt, and then slowly lifted while rotating. The seed diameter (usually a few millimeters) is then expanded to a pre-defined cylindrical diameter and grows at a constant rate of 1-5 mm/h. Once the cylinder reaches its final length, it is separated from the melt and slowly cooled to room temperature. The upper part of the furnace shell is supplied with cold air and the lower part with cooling water. From an engineering point, the Czochralski method is a dynamic phase transition problem with internal moving boundaries and free surface, including heat transfer by the melt and ambient, internal radiation in the crystal, surface-to-surface radiation between the boundaries and ambient, flow and heat transfer coupling in the melt and phase transition in the solid-liquid interface.

Mathematical Model
This model provides a solution to the axisymmetric heat transfer problem that accounts for the heat transfer via conduction and heat exchange between solid surfaces by radiation. The melt flow effect can also be accounted for in the model. The convection heat transfers in the melt and gas are included. The surface-to-surface radiation between the boundaries and ambient and heat exchange via convection are calculated. As using the P1 method can significantly reduce the calculation cost, this method is used to solve the 3dimensional internal radiation problem in the crystal medium with emission and absorption. The steady state is applied to compute the furnace temperature, phase interface, and flow field during crystal growth using the with parameter sweep method [33]. The model is established with parameters. The crucible moves upward and the melt quality decreases as the crystal length changes. By scanning the parameter of crystal length, there is no need to modify the geometry manually, and a steady-state calculation can be carried out for each scanning.
The following assumptions are made to simplify the model: (1) The crystal is fully cylindrical and isotropic (2) The melt is considered as an incompressible liquid. The temperature field in the crucible is axially symmetric.
(3) Due to the relatively low crystallization rates (2 mm/h), the growing process can be seen as quasi-static.
(4) The phase interface shape has little effect on the overall temperature distribution and heater power. It is convex over the entire growth process with slight changes slightly in thickness. Therefore, at this stage, we have neglected the changes in the phase interface.

Thermophysical Properties
The thermophysical properties of the Ga2O3 melt are given in Table 1.
The radiative portion of heat transfer is [41]: where k, Cp, and Q represent the thermal conductivity, the specific heat capacity of the material, and the external heating source, respectively. ε is the emissivity, G is the irradiation, σ is called the blackbody radiation constant, and J is radiosity. According to the following equation, the liquid portion is described by the fluid velocity u and pressure p to simulate the laminar flow. The governing equations are as follows [42]: where ρ denotes the density and µ is the viscosity. For an axisymmetric flow, the equations reduce to: − + + = 1 − + + , Here, is the radial velocity, υ is the rotational velocity, and ω is the axial angular velocity. In the model, Fr and are 0. Different from reference [30], we added the phase material module into the heat transfer calculation. In the module, the phase transition temperature, phase transition interval, and latent heat were set. The solid fraction (the ratio of solid mass to the sum of solid and liquid mass) was used to represent the solidification degree of the melt in the phase transformation zone. For example, the solid fraction of solid is 1, the liquid fraction is 0, and the material in the transition interval is 0-1.
When the crystal is pulled out from the melt, the volume of melt in the crucible decreases. In order to maintain a horizontal position in the phase interface, the crucible moves upward at a certain speed. The crucible rising rate is related to the pulling rate, which is derived by the mass conservation of the crystal and the melt.
where denotes the crystal density, denotes the melt density, and the crystal radius and melt radius are known. is set to 2 mm/h.

Mesh Model and Boundaries
The inner diameter of the crucible was 228 mm, the radius of the furnace was 780 mm, and the height was 1350 mm. There were 11,204 grids with 2478 elements in the melt. The seed was rotated at a speed of 2 rpm and the crystal and crucible rotate at a speed of 6 rpm in opposite directions. The sliding boundary condition was used on the seed and crucible and the axisymmetric boundary condition was adopted. The initial temperature in the system was set to be 300 K. Time iteration step was 0.1 s. Use adaptive mesh refinement and the minimum grid size was 2 mm.

Results and Discussion
The global temperature distribution and isotherms in the process system are shown in Figure 2. There are obvious high-temperature zones and cold zones in the furnace. When the temperature of the phase interface reached 1820 °C, the heather power was about 20 kW which corresponds to the actual situation. This means that the furnace structure and insulation accessories design were reasonable. Based on the accurate calculation of the whole temperature field, the coupling calculation of heat transfer and flow in the crucible was carried out, taking into account the natural convection. The temperature and velocity fields are given in Figure 3a. The solid-liquid interface was calculated by adding the phase change material into the heat transfer module, as shown in Figure 3b. The global simulation during the dynamic growth process is given in Figure 4. During the dynamic growth process, the crystal grows and the material mass decreased, while the relative position of the crucible and heater also changed. The relationship between the crystal growth rate and the crucible rising rate is given by Equation (5). In Figure 4, the surface-to-surface radiation is considered, excluding the internal radiation of the medium. The temperature gradient at the middle section is uniform and it is larger at the crystal shoulder, which agrees well with the reality. In this paper, the effects of the surface-to-surface radiation and the participating medium radiation on the temperature field and heater power are shown in Figures 5-8. As shown in Figure 5, the trend of the heater power curve was completely different during crystal growth. The heater power, when considering the participating medium radiation, was higher than that when considering the surface to surface radiation as a whole. In the case of the surface to surface radiation, the heater power gradually decreased as the length of the crystal grew. This shows that crystal length increased the thermal resistance and reduces heat dissipation. If the participating medium radiation is considered, the power increased greatly with the increase in crystal length. The translucent crystal increases the heat dissipation of the furnace [14]. The power of the heater changed with the crystal growth and the temperature at the bottom of the crucible also changed. As shown in Figure 6, the temperatures at the bottom of the crucible almost stayed constant with the crystal growth in the case of surface-to-surface radiation. Meanwhile, when considering participating medium radiation, the temperature kept rising and the difference was about 60 K. The participating medium radiation has a greater influence on the bottom temperature. The participating medium radiation increased the heat flux at the solid-liquid interface and the temperature at the bottom of the crucible. The temperature at the bottom when considering the participating medium radiation is at least 30 K higher than that when considering the surface-to-surface radiation as a whole. The participating medium radiation of the crystal also has a great influence on the temperature distribution and temperature gradient. In the process of growth, the working point and dynamic characteristics of CZ systems will change permanently. Due to this situation and the nonlinear characteristics of the CZ system, the practicability of a linear PID controller is limited. This is especially true for growing opaque oxide crystals with a large diameter. This problem can be solved using empirical parameters. However, changes in the thermal environment or the enlargement of the diameter necessitate extraordinary effort to be made to adjust all experimental parameters. Figures 5-6 are helpful to understand the variation of power and the temperature at the bottom. As shown in Figures 7-8, when considering participating medium radiation, the temperature gradient of the crystal is lower than that in the case of the surface-to-surface radiation. The temperature gradient is uniform in the middle of the crystal and higher at the shoulder of the crystal. The temperature gradient drops from 20 K/cm to 5 K/cm. In the case of the participating medium radiation, the temperature gradient difference between the middle section and shoulder is larger. It is necessary to improve the thermal radiation of the crystal shoulder to reduce the temperature gradient and make the temperature gradient of the crystal more uniform.    Surface to surface radiation Considering the participating medium In this paper, two models of surface radiation and internal radiation are compared to study the effect of radiation on the process parameters. The power and temperature trajectories at the bottom of the crucible are obtained. The temperature distribution of the furnace, the temperature field, and the flow field of the melt are calculated. The temperature, the temperature gradient of the crystal, the temperature at the bottom of the crucible, and the heater power are studied when the crystal is growing in the two cases of radiation. The results are helpful to adjust the parameters for use in industry. The results obtained in this study show that the participating medium radiation has a great influence on the power and temperature trajectories. The internal radiation influences the dynamic characteristics of CZ systems and reduces the temperature gradient of the crystal.

Conflicts of Interest:
The authors declare no conflict of interest.

Cp
the Specific Heat Capacity Fr force in r direction force in the circumferential direction force in z-direction G irradiation J radiosity k thermal conductivity p pressure Q external heating source r radial coordinate of the melt/crystal T temperature u fluid velocity ε emissivity υ rotational velocity ω axial angular velocity ρ density µ viscosity σ blackbody radiation constant