Simulation of Effect of Ar Flow Rate on Silicon Ingot Growth in Directional Solidification System

1College of Artificial Intelligence, Yango University, Mawei District, Fujian 350015, China 2Department of Mechanical Engineering, National Kaohsiung University of Applied Science, Kaohsiung 807, Taiwan 3Department of Chemical and Materials Engineering, National University of Kaohsiung, Kaohsiung 811, Taiwan 4Department of Aeronautical Engineering, Chaoyang University of Technology, Taichung 413, Taiwan


Introduction
Recently, many reports have focused on the reduction of threading dislocations for the practical production and growth of high-quality polycrystalline silicon from a silicon-based solution, because high-quality silicon ingots can be used to fabricate solar cells with higher solar photovoltaic conversion efficiency and greater cost-effectiveness. (1,2) Growing silicon ingots in a furnace to find the optimal growth parameters in a directional solidification system (DSS) is an expensive and unfeasible method. During the growth of polycrystalline silicon in a DSS, the impact of the Ar flow rate is less than that of radiation heat transfer, but it has a large effect on the melted silicon/gas interface. (3) Gao et al. focused on the effects of the Ar flow rate and pressure on the distributions of oxygen and carbon in a DSS. (4) They found that increasing the Ar flow rate can significantly change the fluid convection and temperature distribution in a DSS. In this study, we performed simulations to investigate the effects of the Ar flow rate, thermal pads, heaters, and insulated side panels on the growth of a silicon ingot in a DSS.
There are five processes in silicon ingot growth: heating, melting, growth, annealing, and cooling. Heating is performed to raise the temperature to about 1450 K in a short time, causing the solid silicon in the quartz crucible to reach the pre-melted condition. During the early stage of the melting process, the temperature is kept at around 1450 K, giving the silicon source the same temperature as the graphite cooling plate. During this process, moisture and impurities are exhausted, and the impurities can be pumped out of the furnace by Ar. Then, the temperature is increased at a pre-determined rate to around 1800 K, and this temperature is maintained to melt the silicon source. Because the crucible walls are close to the heater, the temperature at the periphery of the melted silicon is higher than that at the center, resulting in the generation of convection as a disturbance in the melted silicon. Therefore, we need to view the process from the peephole of the Ar inlet to determine whether the silicon floating on the melted silicon is completely melted.
Because the effect of the Ar flow rate on the internal temperature field is difficult to determine, in this study, we used a numerical method to analyze the effect of its variation. (5,6) The finite volume method is a discretization technique used to calculate fluid dynamics, in which the analyzed area is divided into several control volumes. Each control volume has a main node, and all the control volumes can be integrated into governing equations with various forms. Then the governing equations are used to combine grids with the corresponding initial value and boundary conditions to derive the discretization equations. If ϕ expresses an arbitrary physical quantity, to discretize the governing equations, we first need to express the relationship between an arbitrary control volume and ϕ as an integral. Solvers used in the finite volume method can be divided into density-based solvers and pressure-based solvers. A density-based solver is a coupled solution, and it can solve the discretization equation set and solve all variables simultaneously. A pressure-based solver can solve the variables of the governing equations in order, but the governing equations are mostly nonlinear and coupled; thus, we need to use an iterative method to find the numerical solutions. (7)(8)(9)(10)(11) In this study, we used the boundary conditions and related parameters in a DSS furnace to simulate the growth of polycrystalline silicon ingots under different Ar flow rates. The novel feature of this study is that we used these methods to simulate and find the optimum Ar flow rate for silicon ingot growth. In contrast to other methods, this simulation method can find the optimum Ar flow rate in advance. This decreases the likelihood of failure in silicon ingot growth, reducing the cost of growing silicon ingots in a DSS furnace.

Simulation Parameters
The internal state of an ingot growth furnace is changing over time; thus, the time-related conservation equations used to analyze the silicon crystal growth process are set to reflect the transient state. The boundary conditions of the proposed simulation of the silicon ingot growth in a DSS furnace are shown in Fig. 1, which modeled a system manufactured by GT-Solar Ltd.
The Ar inlet and outlet flow rates were both set at 20 liters per minute (LPM) to keep the pressure of the crystal growth furnace at 200 mbar. Because the ingot growth furnace had a cooling system, the inner furnace walls were kept at 1000 K and the total growth time was about 40 h. In the early stage of growth, the temperature recorded by sensor TC1 shown in Fig. 1 was 1758 K, and by the end of the ingot growth, the temperature had dropped to 1689 K. In the numerical analysis of the heater control function, we used the actual temperature at point TC1. After that, we used the curve fitting method to segment the temperature curve into the five-segment functions, allowing us to control the heater output and thus control the temperature-time relationship. Figure 2(a) shows the relationship between the actually measured temperature curve and the fitting results of the five-segment functions, and the results for each segment are shown in Figs. 2(b)-2(f). (1) When a fluid is in a turbulent flow, it is in a fluctuating condition, and it cannot be expressed exactly using a mathematical model. In this study, we chose the most representative k-ε model, which was based on the turbulent flow energy equation and the dissipation rate of the turbulent flow energy, to find the k value, ε value, and diffusion coefficient of the turbulent flow at different positions simultaneously. The turbulent flow energy equation is The equation for the dissipation rate of turbulent flow energy is The k-ε model is more suitable for use with a fluid with a high Reynolds number, and the important assumption of our method is that the fluid undergoes isotropic diffusion. Using this assumption in the calculation, we can control both the complexity and computational efficiency at the same time; therefore, we used the k-ε model for the simulation. The basic assumptions of the boundary conditions were as follows: 1. a transient analysis was performed, 2. the geometric model was flat and two-dimensional, 3. we ignored chemical reactions, 4. the melted silicon satisfied the Boussinesq approximation, 5. the system had no moveable boundary, and 6. the melted silicon system was recognized as grey body radiation.
We mainly discuss the effect of the Ar flow rate in the crystal growth furnace on the heat transfer of various parts and the phase change of melted silicon. There are three different types of heat transfer: heat conduction, heat convection, and heat radiation. Heat conduction is heat transfer from a material with a high temperature to a material with a low temperature and involves interactions between different molecules. The heat transfer equation is where k is the thermal conductivity coefficient (W/mK), which affects the heat transfer capacity. Heat convection can be divided into natural convection and forced convection. Natural convection is caused by the temperature difference or gradient of liquid fluid, which causes the gas density to change. Because part of the fluid is heated, its volume increases and its density decreases, the gas rises, and loop circulation occurs in this way. Forced convection is caused by introducing an external force (a pump or fan) to induce fluid movement. The heat convection equation is where h is the average heat convection coefficient, which is 5-25 W/m 2 K under the natural convection condition. Thermal radiation is the radiation from a material with temperature higher than absolute zero in the form of electromagnetic waves, which does not require a medium. Thus, a material can transfer heat even in vacuum, where the radiation intensity depends on the body temperature and surface characteristics. The equation for thermal radiation is where ε is the radiation coefficient, whose value is between 0 and 1, and σ is the Boltzmann constant (5.674 × 10 −8 W/m 2 K 4 ). The phase change has four stages, condensation, evaporation, boiling, and solidification/melting, and in this paper, we discuss the effect of the Ar flow rate on the solidification process.

Simulation Processes and Results
The main purpose of introducing Ar into the silicon ingot growth furnace is to load the impurities and exhaust them from the furnace at a high temperature. This can prevent the impurities from diffusing into the melted silicon, increasing the purity of the grown silicon ingot. The carbon impurity can react with Si to form SiC, which has very high hardness, increases the processing risk of a growing silicon ingot, and decreases the silicon yield. The formation of SiC also leads to a drop in yield and reduces the transfer efficiency of the silicon wafer. Under the same operating conditions, we compared the effect of different Ar flow rates (10, 20, 30, 40, and 50 LPM) on the silicon growth process. Figure 3 is a schematic diagram showing all the data extraction points and Fig. 4 shows the temperature-time curves for all the data extraction points as a function of time. With increasing flow rate of Ar, the cooling rates of the melted silicon in the central section and the furnace wall section increase. A change in the Ar flow rate had a larger effect at the top (point C) of the melted silicon, and the cooling rate of the top section was higher than that of the bottom section (point A). The reason for this is that as Ar is introduced into the silicon ingot growth furnace, it first comes in contact with the top of the melted silicon. Table 1 shows the different growth heights of silicon ingots and Table 2 shows the different curing fractions of silicon ingots with different Ar flow rates and different growth times. For the same growth time, the growth height increased with the Ar flow rate. According to Table 2, for the same growth time, the cured fraction of the silicon ingot increased with the Ar flow rate. The results in Tables 1 and 2 also suggest that the Ar flow rate has a large effect on the growth mechanism of silicon ingots.
Concavity is defined as the average slope in the solid-liquid interface of the crucible-wall section, and convexity is calculated from the average slope of the solid-liquid central section calculated from the lowest and highest points. Concavity is defined as a negative value, and a  Convexity dy dx = .
(12) Figure 5 shows solid-liquid phase distribution maps at 30 h under different Ar flow rates. These maps showed that as the Ar flow rate increased, the silicon ingot at the central section exhibited a faster curing rate. The silicon ingot growth height and curing fraction as functions of cooling time are shown in Table 3. The silicon ingot growth height increased up to a time of 25 h, whereas the curing fraction showed no apparent trend with increasing time. When the cooling time was increased, the silicon ingot growth height first increased and reached a maximum value at 25 h; the curing fraction first decreased as the cooling time increased from 5 to 10 h, and then it increased with cooling time, as shown in Table 3. Even when the cooling time increases, the rates of increase of the silicon ingot growth height and curing fraction are inconsistent. However, they tend to increase with the cooling time, and the convexity of the silicon ingot at the central section shows a larger increase with the cooling time. Thus, in the central section, the crystal growth direction of the silicon ingot changed from inward to outward, as shown in Figs. 5(e) and 5(f). When a lower Ar flow rate was used, the convexity of the silicon ingot at the central section was increased, the crystal growth direction of the silicon ingot was more vertical, and the quality of the silicon ingot was therefore higher, as shown in Fig. 5(c). As the Ar flow rate increased, the radial temperature difference at the central section increased. This was caused by the exhaustion of the heat energy at the crucible walls by the Ar gas. Thus, the temperature at point B* dropped quickly, then the growth of crystal silicon occurred at the side walls. As shown in Fig. 6(a), when the Ar flow rate increased, the axial temperature difference in the first 15 h of growth increased slightly, and after 15 h of growth, the axial temperature difference increased markedly. When the silicon at point C solidified, a large amount of heat was exhausted by the Ar gas; therefore, the temperature at point C dropped  quickly, and the axial temperature difference dropped. Figure 6(b) shows the variation of the axial temperature at the central section as a function of the Ar flow rate, which indicates that the axial temperature changes with the Ar flow rate. Figure 7(a) shows the variation of concavity as a function of the Ar flow rate. When the growth time of the silicon ingot was less than 25 h, the concavity did not change with the Ar flow rate. When the growth time was over 25 h and the Ar flow rate increased, the top left and right sections of the melted silicon had a larger heat dissipation rate, the temperature at the crucible walls was lower, and the concavity increased. Figure 7(b) shows the variation of convexity as a function of the Ar flow rate. The convexity increased with the Ar flow rate and was greatest at the Ar flow rate of 50 LPM. The simulation results show that the Ar flow rate had the least effect on the variation of the temperature at the bottom silicon material. The kinetic energy of the distant Ar gas was absorbed by the top melted silicon; therefore, its effect on the shape of the solid-liquid interface was also reduced. The Ar flow rate had the greatest effect at the top of the melted silicon, because the Ar flow rate was highest and the cooling was fastest. This caused the solid-liquid interface of the silicon at the top-central section to undergo crystallization in the outward direction, which was unfavorable for crystallization in the vertical direction. Because the Ar gas pushed down on the silicon ingot, the impurities in the melted silicon were not easily exhausted and they accumulated in the silicon ingot. Because impurity accumulation in the polycrystalline silicon has a negative effect on the transfer efficiency, a larger Ar flow rate should be used at the early stage of growth to increase the speed of the flow of gas on the surface of the melted silicon. This should promote the exhaustion of impurities to outside the furnace and reduce the concentration of impurities in the melted silicon. At the late growth stage, a lower Ar flow rate can be used to prevent the central section of the solid-liquid interface from being excessively convex and to increase the stability and uniformity of the temperature field. The Ar gas allows the silicon ingot to grow in a more vertical direction, and then the air flow can decrease the impurity accumulation at the bottom of the silicon ingot.

Conclusions
We simulated the growth of silicon ingots in a DSS furnace under different Ar flow rates. For the same growth time, the growth height and curing fraction of the silicon ingot increased with the Ar flow rate. This increased the radial temperature difference at the central section; in the first 15 h of the growth, the axial temperature difference increased slightly; after which, it increased markedly. At the early stage of growth, a larger Ar flow rate increased the flow speed of gas on the surface of the melted silicon, thus exhausting the impurities to outside the furnace more effectively. At the late growth stage, a lower Ar flow rate was used to prevent the central section of the solid-liquid interface from being excessively convex and to increase the uniformity and stability of the temperature field. When the growth time was more than 25 h and the Ar flow rate increased, the top left and right sections of the melted silicon had a larger heat dissipation rate, the temperature at the crucible walls was lower, and the concavity increased. When the Ar flow rate increased, the convexity increased, and when the Ar flow rate was 50 LPM, the convexity was greatest. From the simulation results, we demonstrated that the Ar flow rate can affect the temperature distribution in a DSS furnace and thus affect the properties of grown silicon ingots.