Induction Heating and Heat Transfer Modeling - American Chemical

Design of an RF-Heated Bulk AlN Growth Reactor: Induction Heating and Heat Transfer Modeling Wu,*,†

Noveski,‡

Zhang,†

Bei Vladimir Hui Stephen Beaudoin,# and Zlatko Sitar§

Raoul

Schlesser,§

Subhash

Mahajan,‡

CRYSTAL GROWTH & DESIGN 2005 VOL. 5, NO. 4 1491-1498

Department of Mechanical Engineering, State University of New York at Stony Brook, Stony Brook, New York 11794, Department of Chemical and Materials Engineering, Arizona State University, Tempe, Arizona 85287, Department of Materials Science & Engineering, North Carolina State University, Raleigh, North Carolina 27695, and School of Chemical Engineering, Purdue University, West Lafayette, Indiana 47907 Received January 12, 2005;

Revised Manuscript Received April 29, 2005

ABSTRACT: An induction heating and thermal model was developed to design an inductively heated reaction system for bulk AlN crystal growth by powder sublimation. The electromagnetic field and induction heat generation were calculated using the Maxwell equations, and transient temperature distribution was simulated considering conduction and radiation between various components. The complex dependence of the temperature distribution on the radio frequency (RF) coil current and position was demonstrated by numerical studies. A custom-designed hightemperature reactor was constructed and tested. To provide proper temperature control during envisioned sublimation growth, the response of temperature and temperature difference of the bottom and top external surfaces of the growth crucible to power consumption and coil position was measured. The experimental results corresponded well with the simulated temperature responses. 1. Introduction The nitrides of aluminum, gallium, and indium, and their alloys have recently emerged as materials of choice for applications in short-wave light-emitting devices and high-temperature high-power electronics. This has led to numerous theoretical and experimental studies on structural, electronic, and optical properties of both surface and bulk phases of these materials.1 Because of its high electrical resistivity (1010 Ω m), high thermal conductivity (340 W/m K), wide direct band gap (6.28 eV), high acoustic wave velocity, and excellent lattice and thermal expansion match with GaN, AlN is a serious candidate for high-temperature, high-power, and high-frequency applications. AlN will be used to improve III-nitride devices, which have been already produced on SiC and sapphire but which have suboptimal performance and short lifetimes.2 AlN has the same (2H) hexagonal wurzite crystalline structure as GaN with a lattice mismatch along the c-plane of approximately 2.4%. Furthermore, the thermal expansion coefficient of AlN is smaller than that of GaN at temperatures below 300 °C and larger at temperatures above 300 °C. Consequently, the thermal expansion mismatch approaches zero when integrated from room temperature to approximately 1000 °C, which is highly desirable to avoid cracking of epitaxial films during the cooling process.3 These properties make bulk AlN crystal an ideal substrate for III-nitride devices. Different techniques have been applied to AlN bulk growth.4-6 One of the most promising techniques is sublimation growth.5-9 This method requires a reactor * Corresponding author. Tel.: +1-631-632-1490; fax: +1-631-6328544; e-mail: [email protected]. † State University of New York at Stony Brook. ‡ Arizona State University. § North Carolina State University. # Purdue University.

operating at an extremely high temperature, beyond 1800 °C. Both resistively heated8 and inductively heated reactors9 can be employed for this purpose. The inductively heated reactor features fast and contactless heating/cooling, fast sample change, long lifetime, high operating temperature, high-temperature stability, and fast temperature response to power changes. Low growth rate and parasitic nucleation are the major issues in AlN sublimation growth. Both growth rate and crystal quality can be improved if a higher growth temperature is used. Unfortunately, technical and engineering difficulties have limited the highest achievable temperature in practice. Since a preferred reactor should operate in an oxygen-free ambient condition at temperatures up to 2500 °C, the choice of reactor hot zone materials is a major challenge. There is no existing material that has time-independent properties operating at such elevated temperatures for hundreds of hours. In summary, the reaction systems used in the past have encountered the following problems: (1) highest temperature limit: lower than the one that would provide a commercially viable growth rate (∼1 mm/h) for bulk AlN production and 2) nonreproducibility of the experimental results. Numerical simulation is very useful to predict temperature distribution and subsequent growth kinetics in various vapor growth processes because of the difficulty of in-situ measurements. Multidimensional models have been successfully developed to simulate chemical vapor deposition (CVD) processes.9-11 Simulations have also been performed to study transport phenomena in various SiC sublimation growth systems.12-15 In this paper, an induction heating and thermal model is developed to serve as a guide in the design of an inductively heated sublimation reactor that operates at temperatures up to 2500 °C and pressures up to 800 Torr. Numerical simulations will be performed to predict

10.1021/cg050014m CCC: $30.25 © 2005 American Chemical Society Published on Web 06/11/2005

1492

Crystal Growth & Design, Vol. 5, No. 4, 2005

Wu et al. Table 1. Thermophysical Properties Used in the Simulations16,18,19 properties thermal conductivity of AlN (crystal) (1000 K) thermal conductivity of nitrogen (1000 K) thermal conductivity of insulation (300 K) thermal conductivity of insulation (300 K, in a vacuum) thermal conductivity of graphite (300 K) density of insulation density of dense graphite density of AlN (crystal) density of nitrogen (1000 K) porosity of AlN powder heat capacity of graphite (300 K) heat capacity of AlN (crystal) heat capacity of nitrogen (1000 K) electrical resistivity of insulation electrical resistivity of crucible

Figure 1. Schematic diagram of inductively heated reactor used for heat transfer simulations and growth experiments.

the temperature distribution and electromagnetic field in the reactor. The predicted dynamic temperature response will be used to design the reactor operating procedure, with special focus on the protocol for controlling the reactor temperature, and an experimental dynamics study in the constructed reactor will provide data for validation of the model predictions.

values 40 W/(m K) 0.0647 W/(m K) 0.347 W/(m K) 0.252 W/(m K) 115 W/(m K) 170 kg/m3 1760 kg/m3 3250 kg/m3 0.3368 kg/m3 0-40% 720 J/(kg K) 1080 J/(kg K) 1167 J/(kg K) 2.59 × 10-3 (Ω m) 1.17 × 10-5 (Ω m)

The transient temperature distribution is calculated by solving transient energy equation. The nitrogen is approximated to be stagnant. The order of magnitude analysis reveals that the maximum nitrogen flow rate of 100 sccm through the 6-7 in. diameter tube has insignificant contribution to heat transfer comparing with radiation. The energy equation can then be written as

(Fcp)eff

∂T ′′′ + qlatent ′′′ + ∇‚q′′radi ) ∇‚(keff∇T) + qeddy ∂t

2. Mathematical Models 2.1. Inductive Heating and Heat Transfer. The inductively heated sublimation reactor shown in Figure 1 is simulated in this paper. The electromagnetic field in the azimuthal direction of the coil induces eddy currents in the conducting susceptor (crucible) placed inside the coil. The electromagnetic field and Joule heating produced by the induction heating are calculated by solving Maxwell’s equations, as described below. The temperature distribution inside the growth chamber is simulated using an energy conservation equation that accounts for heat generation and conductive and radiative heat transfer. The azimuthal component A0 (A0 ) Ar + iAi) of the magnetic vector potential is calculated for the inductive coils. The Maxwell equations can be expressed in terms of Ar and Ai by the classical electromagnetic theory:16

(

)( ) )( )

∂2 1 1 ∂ ∂2 Ar - 2+ 2 + + mω2Ar + ωσcAi ) -J0 2 r ∂r µ ∂r r ∂z m (1)

(

∂2 1 1 ∂ ∂2 Ai + + + mω2Ai - ωσcAr ) 0 ∂r2 r ∂r r2 ∂z2 µm (2)

where the real part, Ar, is the in-phase component, the imaginary part, Ai, is the out-of-phase component, µm is the magnetic permeability, ω ) 2πf is the angular frequency of the AC current, σc is the electrical conductivity, m is the electrical permittivity, and J0 is the current density in the azimuthal direction in the coil. The computations are performed in an axisymmetric domain and are large enough so that the boundary values of both the in-phase and out-of-phase components can be set to zero, i.e., Ar ) Ai ) 0 at r f ∞, z f -∞ and z f +∞.

(3)

where q′′radi is the radiative heat flux on the surfaces of inner enclosure, top enclosure, top and bottom windows ˙ is the latent heat (see Figure 1), q′′′ latent ) ∆Hvs˘ M absorbed/released during the dissociation/deposition. The effective heat capacity (Fcp)eff and thermal conductivity keff for the porous medium are estimated from ref 17. Within a nonporous solid or the gas, (Fcp)eff and keff are the heat capacity and thermal conductivity of the solid or the gas, respectively. In Figure 1, we illustrate schematically that radiative heat transfer is important in the top enclosure and inside the crucible. The heat generated by the eddy current in the graphite susceptor 2 2 2 can be expressed as q′′′ eddy ) 1/2 σcω (Ar + Ai ). The properties of each component in the growth system are listed in Table 1.16,18,19 2.2. Numerical Methods. Nitrogen gas species are considered totally transparent, and all the surfaces are assumed to be gray, opaque, and diffusive. The grid-togrid gray-diffusive method is used in this paper. The radiation surfaces are divided into rings, each with uniform properties. The rings coincide with the grids used for temperature calculation. The view factor between each pair of rings is calculated.14 The total exchange factor method (total exchange factors between grids) is used to relate temperature with the heat flux.14 q′′radi is calculated by solving a radiative equation for gray and diffusive surfaces:20

q′′radi,j j

N

-

∑

k)1

Fj,k

N 1 - k q′′radi,k ) σTj4 Fj,kσTk4 k k)1

∑

(4)

where Fj,k is the view factor from ring j to ring k, is the emissivity of the radiative surface, and σ is the Stefan-Boltzmann constant. The electromagnetic and

Design of an RF-Heated Bulk AlN Growth Reactor

Crystal Growth & Design, Vol. 5, No. 4, 2005 1493

Figure 2. Double-walled quartz tube design: water jacket and cooling baffles.

energy equations including radiation are solved iteratively using an in-house developed code.13 Different grid requirements are needed to calculate the magnetic potential and temperature fields. For example, accurate prediction of the magnetic potential requires fine grids in the vicinity of the coil and near the outer surface of the graphite susceptor (the skin depth in the crucible) as well as a large computational domain to satisfy the zero boundary conditions. The skin depth of the graphite susceptor is about x2/2πfµmσc ) 1.7 cm. The geometry inside the crucible is not important for magnetic potential calculation since the electromagnetic field is unable to penetrate into the crucible. On the other hand, the thermal field is substantially affected by the detailed configuration of the reactor, i.e., each component and detailed design of the reactor. However, a smaller computational domain is needed for temperature simulation. A bilinear interpolation method is used to interpolate data between the two grid systems. The magnetic potential is calculated first using the magnetic grid system, and volumetric heat generation on the thermal grid is evaluated by interpolation of data obtained in the magnetic grid system. The energy equation is then solved. Iterative approximations are needed if electric properties are dependent on temperature. 3. Reactor Design 3.1. Water Cooling. Efficient cooling of the radio frequency (RF)-heated reactor is mandatory because all reactor components must be maintained at temperatures far below their melting points. Figure 2 shows the location of three parallel cooling circuits in a crosssection of the reactor design. Chilled water flows in an upward direction through the double-walled quartz tube to extract the heat transferred from the susceptor in a

Figure 3. Reactor cross-section, thermal insulation (crosshatched), and other internal elements.

predominantly radial direction. In addition, customdesigned cooling baffles are used to prevent overheating of the bottom and the top reactor flanges. The baffles also prevent deposition of crystalline AlN on the internal reactor flanges. The heat flow removed by the cooling water can be calculated as Fcp∆TV/∆t. The density and heat capacity of water near room temperature are F ) 1000 Kg/m3, and cp ) 4180 J/Kg K, respectively. A water temperature change of ∆T ) 15 K is assumed for normal operation. From the thermal model, it is found that the current density of 800 A, at power output of about 8 KW, will provide temperatures at the susceptor in the range of 2100-2200 °C. To safely estimate the required cooling, it is simply assumed that all the input power must be removed by water cooling. It is further assumed in the calculation that the required heat removal is in a range from 6 to 10 KW in all operations. The chilled water flow rate needed in the double-walled quartz tube can then be estimated to be 1.5-2.5 gal/min, i.e., 9.5 × 10-5-1.6 × 10-4 m3/s. 3.2. Thermal Insulation. Custom designed, modular graphite foam insulation elements are used to tailor the temperature gradient inside the crucible, as shown in Figure 3. By changing the insulation thickness below and above the hot zone, control over the heat transferred in the axial direction at steady state can be achieved. Similarly, the radial heat transfer can be controlled by changing the thickness of the insulation around the crucible. For example, to achieve higher temperatures at the bottom of the crucible the insulation thickness below the crucible will be larger, and hence the heat transferred is lower in that particular direction. Model calculations showed that reducing the thickness of the

1494


Wu et al.

Figure 4. Magnetic field simulation results for a coil current of 980 A and frequency of 10 kHz.

insulation above the reactor hot zone to approximately one-third of the insulation thickness below the hot zone yields a temperature drop along the crucible of 100 °C, which is required for typical sublimation growth. During growth experiments or in the thermal response experiments described below, the insulation is not mobile, and dynamic changes in the temperature gradient can only be achieved by adjusting the position of the coils relative to the crucible. 4. Results and Discussion 4.1. Global Reactor Heat Transfer Simulations. Figure 4 shows the distributions of the in-phase (Ar) and out-of-phase (Ai) components of the magnetic field at coil current of 980 A and frequency of 10 kHz, which is the typical operating condition. The top of the coil is lower than the top insulation disk by 1.3 in. The in-phase component reaches its largest field strength around the coil, with the contour shape being obviously influenced by the geometry of the coil, while the maximum absolute value of the out-of-phase component concentrates in the graphite crucible. Interestingly, their absolute values are about the same order of magnitude (∼10-4 Wb/m) in the crucible. Therefore, both the in-phase and the outof-phase components need to be considered in the heat generation calculations. The axial position of the coils relative to the crucible is chosen so that more heat is generated in the bottom portion of the crucible than in the top portion. An axial temperature gradient is needed in the sublimation growth system for providing a driving force for mass transfer. The temperature distribution of the entire growth system is presented in Figure 5. The highest temperature is located at the lower portion of the

Figure 5. Temperature distribution in the growth system for coil current of 980 A and frequency of 10 kHz.

crucible wall, at an axial position corresponding to the center of the induction coil. The temperature difference between the source surface and the crucible lid is about 40 K. The axial temperature gradient is closely related to the axial position of the coil and the diameter of the top window from which heat is transferred out of the chamber. The radial temperature difference between the center and the periphery of the crucible is very small due to the thick insulation layer and comparatively small diameter of the crucible.



Table 2. Effects of Coil Current and Frequency on Heating Power and Maximum Temperature coil current (A) 800 850 850 850 900 frequency (KHz) 10 8 10 12 10 calculated power (W) 1879 1367 2122 3032 2380 2 2 R ) I f /P 34061 33826 34048 34314 34034 R ) I2f 1.5/P 10771 11959 10767 9905.6 10762 R ) I2f/P 3406 4228 3404 2859 3403 R ) I2/P 340.6 528.5 340.5 238.3 340.3 maximum temp (K) 2156 1970 2263 2674 2387

Table 3. Effects of Frequency on Maximum Temperature intended power (W) intended maximum temp (K) frequency (KHz) coil current (A) I ) (RP/f 2)1/2 R ) 34048 I from calculationa numerical predicted power (W) numerical predicted maximum temp (K)

2122 2263 8

8

1062.5 2138 2278

12

12

708.3 1058.5 2122 2271

2104 2237

711.4 2122 2246

a Current I is calculated iteratively to make the heating power reach 2122 W.

The effects of coil current and frequency on maximum temperature and heating power are presented in Tables 2 and 3. From the heat generation expression, q′′′ heater ) 1/2 σcω2(Ar2 + Ai2) and the skin depth formulation, ds ) x2/2πfµmσc, it can be derived that the heat generation should be proportional to frequency with the power in the range of 1.5-2.0. For the cases we simulated, R ) I2f b/P, is calculated with different b values, where I is the current and R is a constant. The R value is almost a constant when b ) 2 is selected. The heating power generated by induction coils is therefore proportional to the square of the current and almost square of the frequency, i.e., I2f 2. R is 34048 for the power of 2122 W with current of 850 A and frequency of 10 kHz in Table 2. This case is considered as the baseline case to study the effects of the frequency. Table 2 also shows that both power and maximum temperature will increase with current and/or frequency. Table 3 shows that when the frequency is increased from 8 to 12 kHz while the current is decreased from 1058.5 to 711.4 A, the maximum temperature is decreased from 2271 to 2246 K even under the same power of 2122 W. This is because the skin depth decreases with the frequency. 4.2. Temperature Response to Heating Power and Coil Position. Dynamic studies are performed in the designed reaction system to understand the temperature response to the power changes and temperature gradient response to the coil position changes. Dynamic response studies are necessary for process control optimization. The time dependence of the temperature response, as shown in Figure 6, does not have an inflection point, which indicates that the system response can be approximated as being of first order. The temperature change at the center location on the top surface of the crucible, [T(t) - T(0)], is defined as a deviation variable of the temperature, which quantitatively represents the dynamic changes of the temperature from its initial steady state before the power changes. If the system is approximated as a first order, its temperature dynamics can be described by a firstorder differential equation. However, the coefficients in the first-order differential equation appear to be slightly

Figure 6. Dynamical response of temperature to power changes.

power dependent, resulting in nonlinearity to the firstorder approximation. Experiments and simulations both show that the equal power change at a higher power level results in a faster temperature response but a smaller magnitude of the temperature change, as shown in Figure 6. When the power supplies are 7, 8, and 9 KW, the temperatures at steady state, T(0), are 1740, 1842, and 1937 °C, respectively. The experimentally estimated time constant is 15-20 min, and the gain is about 100 to 85 °C. The time constant and the gain are estimated from the system response graph: time constant as one-third of the time period required for the system to reach 95% of its final steady-state value, while the gain as a difference of the final and initial steadystate temperature. Assuming the first-order system, the gain-time constant form of the system transfer function G(s) in the Laplace domain can be written:

T(s) k ) G(s) ) τs + 1 P(s)

(5)

Empirically, the functions of the gain (k) and the time constant (τ) of the power can be determined from measurements at different power changes. Figure 7 shows the temperature distributions of the entire growth system at different time periods (5, 10, 20, and 50 min) when the system is heated from room temperature for coil current of 900 A and frequency of 10 kHz. The initially large temperature increase observed in the insulation layers may be explained by the fact that the graphite foam insulation has low thermal mass and, consequently, weak coupling with the RF field is sufficient to induce the observed temperature increase. Simulation results also show that for a sudden power change the system reaches 95% of its final steadystate value after about 50 min. The quick thermal response of the system is due to the small thermal mass of the hot zone components. Figure 8a,b shows the vertical center-line temperature profiles and the horizontal temperature profiles cross the highest temperature. The highest temperature is determined by the position of the coils (Figure 8a). The temperature in the source material is almost

1496


Wu et al.

Figure 7. Predicted temperature distributions at coil current of 900 A and frequency of 10 kHz for different heating time periods: (a) 5 min; (b) 10 min; (c) 20 min; (d) 50 min.

constant (Figure 8b) due to the large effective thermal conductivity of AlN powder. The effective thermal conductivity of the porous AlN powder is the interplay of a low conductive heat flux and a high radiative flux.17 At 5 min, Figure 8b shows that the second highest temperature peak is located inside the insulation material because the insulation material is also electrically conductive and heat can be generated. However, heat generation in the insulation material is much weaker than that in the graphite crucible. This will result in a large temperature drop between the crucible and the outside furnace wall. The temperature in the furnace changes significantly from 5 to 20 min. After 20 min, the temperature in the furnace will vary slowly. After

50 min, the highest temperature reaches 95% of the final value. Likewise, the response of the temperature difference of the bottom and top of the crucible to coil position change shows a nonlinearity with power, albeit of much smaller magnitude, so that it can be assumed to be power-independent. The temperature difference has a much faster response to the coil position changes than to the power changes. Also, even for a large position change, the system will reach a new steady state in approximately 15 min. Figure 9 demonstrates the response of the temperature difference as the coil movement of 2, 4, or 6 mm at 10 KW. Since a combination of adjusting coil position and heating power allows



Figure 9. Dynamical response of temperature difference changes to coil position movements.

coil current and coil position is demonstrated by numerical simulations. The numerical model is validated by comparison of the predictions with the experimental measurements. Numerical simulations have demonstrated that the heating power is approximately proportional to I2f 2. Temperature dynamics study shows that the system responds rapidly to the heating power and coil position changes. The desired temperature distribution can be established by adjusting the coil current and position adequately using a feed forward control scheme. Time constant and gain can be obtained by either simulations or experiments.

Figure 8. Predicted temperature profiles at different times along (a) the center line; (b) the height of the highest temperature.

one to obtain the desired temperature distribution in the growth system, the furnace temperature and its gradient can therefore be controlled during growth. The initial test of the reactor reveals that the temperature stability of the system is excellent and a feed forward control scheme of the heating power and coil position is sufficient to attain and maintain the desired crucible temperature distributions. In particular, at 2000 °C, the process temperature varied only within (3 °C over a period of 50 h. 5. Conclusion Thermal modeling of an inductively heated sublimation system for bulk AlN growth is performed, and temperature distributions are predicted. Simulation has been used as a guide in the design of the inductively heated reactor. The cooling system consisting of a double-walled quartz tube water jacket and water-cooled baffles on the top and bottom of the reactor provides safe and stable high-temperature operation, up to 2500 °C. The dependence of temperature distribution on the

Acknowledgment. We would like to acknowledge the financial assistance received through the DoD Multidisciplinary University Research Initiative (MURI) program administered by the Office of Naval Research under Grant N00014-01-1-1-0716 monitored by Dr. Colin E. Wood. References (1) Kandalam, A. K.; Pandey, R.; Blanco, M. A.; Costales, A.; Recio, J. M.; Newsam, J. M. J. Phys. Chem. B. 2000, 104, 4361-4367. (2) Liu, L.; Edgar, J. H. Mater. Sci. Eng. R 2002, 37, 61127. (3) Rojo, J. C.; Slack, G. A.; Morgan, K.; Raghothamachar, B.; Dudley, M.; Schowalter, L. J. J. Cryst. Growth 2001, 231, 317-321. (4) Bliss, D. F.; Tassev, V. L.; Weyburne, D.; Bailey, J. S. J. Cryst. Growth 2003, 250, 1-6. (5) Slack, G. A.; McNelly, T. F. J. Cryst. Growth 1977, 42, 560563. (6) Balkas, C.; Sitar, Z.; Zheleva, T.; Bergman, L.; Nemanich, R.; Davis, R. F. J. Cryst. Growth 1997, 179, 363-370. (7) Schlesser, R.; Dalmau, R.; Sitar, Z. J. Cryst. Growth 2002, 241, 416-420. (8) Wu, B.; Ma, R.-H.; Zhang, H.; Dudley, M.; Schlesser, R.; Sitar, Z. J. Cryst. Growth 2003, 253, 326-339. (9) Noveski, V.; Schlesser, R.; Mahajan, S.; Beaudoin, S.; Sitar, Z. J. Cryst. Growth 2004, 264, 369-378. (10) Iwanik, P. O.; Chiu, W. K. S. Numer. Heat Transfer, Part A 2003, 43, 221-237. (11) Yoo, H.; Jaluria, Y. J. Heat Trans. 2002, 124, 938946. (12) Ma, R.; Zhang, H.; Prasad, V.; Dudley, M. Cryst. Growth Des. 2002, 2, 213-220.

1498


(13) Ma, R.-H.; Zhang, H.; Ha, S.; Skowronski, M. J. Cryst. Growth 2003, 252, 523-537. (14) Chen, Q.-S.; Zhang, H.; Prasad, V.; Balkas, C. M.; Yushin, N. K. J. Heat Transfer 2001, 123, 1098-1109. (15) Klein, O.; Philip, P. J. Cryst. Growth 2003, 249, 514-522. (16) Jackson, J. D. In Classical Electrodynamics; John Wiley & Sons: New York, 1998. (17) Kaviany, M. In Principles of Heat Transfer in Porous Media; Springer, New York, 1995.

Wu et al. (18) Incropera, F. P.; De Witt, D. P. In Fundamentals of Heat and Mass Transfer, 3rd ed.; John Wiley & Sons: New York, 1990. (19) Slack, G. A., Tanzilli R. A., Pohl R. O., Vandersande, J. W. J. Phys. Chem. Solids 1987, 48, 641-647. (20) Modest, M. F. In Radiative Heat Transfer; McGraw-Hill: New York, 1993.

CG050014M

Induction Heating and Heat Transfer Modeling - American Chemical

Recommend Documents