R A D I A T I V E H E A T T R A N S F E R IN G A L L I U M A R S E N I D E L E C C R Y S T A L P U L L E R S By Muna Bakeer B. Sc. (Civil Engineering) B. Sc. (Mechanical Engineering) M. Sc. (Mechanical Engineering) Washington University, St. Louis, Missouri A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF MASTER OF APPLIED SCIENCE in THE FACULTY OF GRADUATE STUDIES MECHANICAL ENGINEERING We accept this thesis as conforming to the required standard THE UNIVERSITY OF BRITISH COLUMBIA November 1990 © Muna Bakeer, 1990 /. In presenting this thesis in partial fulfilment of the requirements for an advanced degree at the University of British Columbia, I agree that the Library shall make it freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may be granted by the head of my department or by his or her representatives. It is understood that copying or publication of this thesis for financial gain shall not be allowed without my written permission. Department The University of British Columbia Vancouver, Canada DE-6 (2/88) Abstract A numerical analysis of radiative heat transfer in a liquid encapsulant Czochralski gallium arsenide crystal puller is developed. The heat transfer and equivilent ambient temper-ature of each surface element axe calculated using the Gebhart radiative model. The effective ambient temperature, to which each surface element is radiating, is found to vary indicating that assuming a constant ambient temperature for all surfaces (simplified radiative model) is incorrect. The importance of including the middle and top cylinders of the growth chamber in numerical analysis of radiative heat transfer in the system is evaluated in the study. The upper section could be replaced by one isothermal surface without significant change of the effective ambient temperature distribution. Fluid flow and heat transfer in the GaAs melt, crystal and encapsulant are calcu-lated using a three dimensional axisymmetric finite difference code which includes the detailed radiative model. The mathematical modelling of the fluid and heat flow de-scribes steady state transport phenomena in a three dimensional solution domain with latent heat release at the liquid/solid interface. The predicted flow and temperature fields using the detailed radiative model differ considerably from the predicted fields using the simplified model. The simplified model shows high axial and low radial temperature gradients in the crystal near the encapsulant region; the axial gradient decreases and the radial gradient increases with increasing distance from the encapsulant top. The detailed model shows a high radial temperature gradient in the crystal near the crystal-encapsulant-ambient junction and nearly flat isotherms in the top half of the crystal. ii Nomenclature A area, m2 ; optical absorbance Ci convective flux from boundary i Cp specific heat, J/gmK Di diffussive flux from boundary i F{j configuration factor from surface i to surface j /,' energy fraction emitted at T< below A = 2p.m fi length ratio Gij Gebhart factor for surface i h' height, m hconvection convective heat transfer coefficient, W/cm2K k thermal conductivity, W/cm.K m-i mass flow rate in face i of the control volume p pressure, N/m2 Pe Peclet number r radius, m Rj, residual of <f> S source term T temperature, K u, v, w x—, T—, 8—direction velocity component, respectively, m/s Up crystal pulling velocity, m/s V volume, m 3 a underrelaxation factor; absorption coefficient c m - 1 8 thermal expansion coefficient, K~l iii Sx, Sr, 69 x—, r—, 9— direction distance between two adjacent grid points Ax, Ar, A0 x—, r—, 8— direction width of the control volume e emissivity T general diffusion coefficient, kg/(ms) 7 approximation coefficient A wavelength, fim fi dynamic viscosity, Ns/m2 p' density, kg/m3 p reflectivity <f> general dependent variable a Planck's constant, 5.729 * lQ-8W/m2K4 Subscripts amb ambient c crucible co outer wall of crucible E neighbour in the positive x-direction, i.e., on the east side e control-volume face between P and E; encapsulant enc enclosure m melt mc middle cylinder wall mcla middle cylinder lower annulus meta middle cylinder top annulus N neighbour in the positive r-direction, i.e., on the north side n control-volume face between P and N P central grid point under consideration iv s 8 St tew tct w w Superscripts * neighbour in the negative r-direction, i.e., on the south side control-volume face between P and S; crystal crystal top top cylinder wall top cylinder top neighbour in the negative x-direction, i.e., on the west side control-volume face between P and W previous-iteration value of a variable v Table of Contents Abstract i i Acknowledgement x 1 I N T R O D U C T I O N 1 1.1 General 1 1.2 The Liquid Encapsulant Czochralski Crystal Growth Process 2 1.2.1 Process History 2 1.2.2 Process Description 2 2 L I T E R A T U R E S U R V E Y 6 2.1 Introduction 6 2.2 Experimental Studies 6 2.3 Theoretical Studies 9 2.4 Scope and Objectives of the Present Work 14 3 Radiative Heat Exchange in L E C GaAs Crystal Pullers 16 3.1 Introduction 16 3.2 Crystal Puller Geometry Assumptions 16 3.3 Mathematical Modelling 17 3.3.1 Configuration Factor Calculation 17 3.3.2 Description of the Analysed Cases 18 3.4 Heat Transfer Model 24 vi 3.4.1 Diffuse-Gray Surfaces 25 3.4.2 Encapsulant Layer Semitransparency 25 3.4.3 Materials' Thermophysical Properties 29 3.5 Temperature Distribution of the Surfaces 36 3.6 Analysis of Results 38 3.6.1 Five Surfaces Enclosure 38 3.6.2 Analysis of the results for case-I 38 3.6.3 Analysis of the results for case-II 40 3.6.4 Analysis of the results for case-Ill 44 3.6.5 One Surface Enclosure 48 3.6.6 Analysis of Case-I 48 3.6.7 Analysis of Case-II 50 3.6.8 Analysis of Case-Ill 52 3.7 Results of the simplified Case-Ill with Cutoff Wavelength 2.5 fim 54 3.8 Comparison of the Results of the Full and Simplified Chambers 56 3.9 Conclusions 58 4 Mathematical Modelling of Fluid and Heat Flow 59 4.1 Introduction 59 4.2 Assumptions 59 4.3 Governing Equations 60 4.4 Boundary Conditions 61 4.5 Solution Procedure 63 4.5.1 Discretization of the Differential Equations 64 4.5.2 Discretization of the Boundary Conditions 68 4.5.3 The Solution Algorithm 69 vii 4.5.4 Solution of the Linear Algebraic Equations 70 4.5.5 Numerical Considerations 71 4.6 Simplified Radiative Model Results 74 4.7 Detailed Radiative Model Results 75 4.8 Comparison between the Radiative Models 75 5 Conclusions and Recommendations 86 5.1 Summary 86 5.2 Conclusions 86 5.3 Recommendations 87 Appendices 89 A Configuration Factors of the System 89 A . l Introduction 89 A.2 Configuration Factor Equations 89 A. 3 Configuration Factor Calculation of Case-I 102 B Gebhart Factor Equations Derivation 106 B. l Introduction 106 B.2 Melt Gebhart factor equations 106 B.3 Gebhart factor equations for the crucible inner wall 136 B.4 Gebhart factor equations for the middle cylinder wall 136 B.5 Gebhart factor equations for the crystal and the crystal top 137 B.6 Gebhart factor equations for the crucible outer wall 137 B.7 Gebhart factor equations for the middle cylinder lower annulus 137 B.8 Gebhart factor equations for the middle cylinder top annulus, top cylinder wall, top cylinder top 137 viii B. 9 Gebhart factor equations for the encapsulant 138 C Programme Flow Chart 139 C. l General 139 D Configuration and Gebhart Factors Calculations Results 144 D. l General 144 E Stainless Steel Emittance 151 E. l General 151 Bibliography 156 ix Acknowledgement I wish to thank Dr. M. Salcudean for giving me the opportunity to work on a very interesting project, and for her financial and moral support during my stay at UBC. I also wish to thank Dr. M. Iqbal for his financial support, Dr. I. Gartshore for his continuous support and encouragement and Dr. P. Sabhapathy for his valuable comments. My thanks go to Mr. A. Steeves who was very helpful in offering assistance with using VAX/VMS system at the Department of Mechanical Engineering. I would like to thank all my friends especially C. Backhouse, V. LeMay and P. Pottier whose continuous support made my path easier to travel. x Chapter 1 INTRODUCTION 1.1 General Gallium arsenide is an electronic material that is emerging, not as a substitute for silicon, but as an important complement to it. Gallium arsenide has many properties that are useful in electronic equipment. Some of those properties include the following: • Speed/power capability that is significantly superior to that of silicon. • Higher radiation resistance of its devices when compared with silicon devices. • Operation at higher temperatures than silicon. • Higher efficiency for solar cells in space applications. • Combining the processing of both light and electronic data on a single chip. • The possibility of being alloyed with other III-V compounds, such as AlAs, InAs, InP and GaP for a broad spectrum of useful properties. Gallium arsenide does however have one major disadvantage - COST. The gallium in a GaAs wafer costs more than ten times as much as the polysilicon in silicon wafers. Therefore, it is imperative to produce crystals with low dislocation densities and point defects in order to improve production yield. GaAs is composed of the low melting point metal gallium and the volatile metal arsenic. Consequently, it shows a very high arsenic 1 Chapter 1. INTRODUCTION 2 dissociation pressure (about 0.1 MPa) at the melting point (1238°C). This fact makes it hard to grow stable crystals [1]. 1.2 The Liquid Encapsulant Czochralski Crystal Growth Process 1.2.1 Process History Metz, Miller and Mazelsky [2] first introduced the liquid encapsulant Czochralski (LEC) method in the early sixties. They pulled PbTe and PbSe crystals using molten boron trioxide (B2O3) as the liquid encapsulant. They chose molten B2O3 because it is less dense than the materials pulled and thus it floats on top of the semiconductor melt. In addition, it has a low vapour pressure which means that it will not evaporate under atmospheric or higher pressures. Also, the liquid metals they used were not soluble in and did not react with the molten boron trioxide. One of the other beneficial properties of boron trioxide is that it is transparent in the visible region which permits direct observation of the melt surface. It is rather viscous and tends to cling to the pulled crystal covering a large fraction of its surface, thus preventing arsenic evaporation from the crystal. However, it is soluble in hot water and may easily be removed from the finished crystals. B2O3 is also used as an encapsulant in the LEC crystal growth of GaAs. 1.2.2 Process Description The LEC growth process is a major means of producing GaAs crystals. The LEC growth chamber is shown in figure 1.1. It consists of a cylindrical crucible that holds the melt and encapsulant. The crucible is heated from the side. The encapsulant covers the melt; its thickness varies depending on (a) the crystal length and diameter, and (b) whether high or low pressure is used (1 versus 4 MPa). A typical thickness is 2-cm. The optimum thickness has not been determined yet, even though full encapsulation seems to allow Chapter 1. INTRODUCTION camera top cylinder middle cylinder — > pulling rod Quartz crystal resistance heaters crucible pedestal crys tal GaAs melt D • 15.24 cm (6in) h = 10.16cm (4in) h°37.25cm (14.7in) D=30.5cm (12/n) JD..'..7.6cm (3in) .Encapsulant . PBN crucible h = 13.8cm (5.4in) D - 15.24 cm (6in) Figure 1.1: Schematic showing GaAs liquid encapsulant Czochralski crsytal puller. Chapter 1. INTRODUCTION 4 the growth of least stressed crystals. Full encapsulation implies that the whole crystal is covered with the encapsulant material as it grows. A seed that is attached to the pulling rod is inserted from the top so that it touches the melt surface. The gas surrounding the crystal - nitrogen or argon - is pressurized up to 4-MPa in order to reduce arsenic evaporation from the crystal and melt. The crucible and crystal are rotated either in the same direction (isorotation) or in opposite directions (counterrotation), at the same or at different speeds. The crucible is rotated to minimize the thermal asymmetry of the system. The crystal is rotated to create a uniform boundary layer at the solid/liquid interface through which heat and mass transfer take place. This boundary layer isolates the growth interface from the velocity fluctuations in the bulk of the melt. Crystal rotation also helps in shaping a cylindrical crystal. The crucible is raised by movement of the pedestal so that the sohdification front remains in a specified region of the heater as the melt level drops. The growth process is viewed through a quartz crystal that is positioned as shown in figure 1.1. Maintaining the same heat transfer conditions throughout the growth of the crystal is complicated by imperfections in the heat transfer environment and by instabilities thought to be inherent in the Czochralski growth configuration. The process is also complicated by the batchwise nature caused by the decreasing depth of the melt in the crucible and the changing heat exchange characteristics between the melt, crucible and ambient gas. The use of B%Oz as an encapsulant introduces the following difficulties: 1. Crystal diameter control becomes difficult because B2O3 which has a large heat capacity, induces a time lag in the temperature control and thereby in the diameter control. Chapter 1. INTRODUCTION 5 2. The large temperature gradient of B2O3, which has a small thermal conductiv-ity, causes a high temperature gradient in the crystal which induces high thermal stresses leading to an increase in dislocation density. The dislocations in turn alter the electrical properties of the GaAs wafers thus degrading device performance. (Dislocations are also caused by propagation and multiplication from an imperfect seed.) 3. Reaction between B2O3 and impurities in the gallium arsenide melt affects impurity concentration in the grown crystal, although B2O3 rarely reacts with the melt itself. Thus, it is important to study the fluid flow and heat transfer in the melt, encapsulant and surrounding gas to understand how they affect the growing crystal and to determine the best conditions for growing GaAs crystals. Also, the information provides data needed for the design of the crystal pullers. Chapter 2 LITERATURE SURVEY 2.1 Introduction The areas of research in the crystal growth field that are of interest here include : • Fluid flow analysis of the melt, encapsulant and gas. • Thermal flow (heat transfer) in the melt, encapsulant, crystal and furnace environ-ment. The literature survey includes both experimental and analytical work. It is found that most studies concentrated on fluid or heat transfer analysis in the melt and did not include the encapsulant as a part of the system, i.e., the analysis was done on a Czochralski - not a liquid encapsulant Czochralski - crystal growth system. However, some recent studies did include the encapsulant region. 2.2 Experimental Studies In 1968 Mullin et al [3] used the liquid encapsulation method to pull GaP and InP crystals. They designed a high pressure chamber which could have been obtained as a standard attachment to commercial versions of the crystal puller available then. Their study was aimed at confirming the viability of the technique and assessing the chemical purity of crystals grown in silica and vitrious carbon crucibles. 6 Chapter 2. LITERATURE SURVEY 7 They concluded that silica crucibles did contaminate GaP and InP melts. However, the vitrious carbon crucibles did not contaminate the melt. They found that the starting material was the most likely potential source of carbon. Lamprecht et al [4] did experiments to investigate the flow and temperature distri-butions in a Czochralski configuration. The experiments were designed to study the importance of buoyacy, thermocapillary and forced convection on the crystal growth. They used NaNOa melt in an iridium crucible. They reported qualitative agreement be-tween their experimental results and numerical results obtained by Langlois in [5, 6, 7]. The agreement showed the importance of thermocapillary forces when considering the problem of fluid flow transition leading to interface transitions. They also emphasized the importance of buoyancy and thermocapillary forces for heat dissipation from the melt surface as well as mass transport through the melt surface. Terashima and Fukuda [8] designed a new magnetic field LEC pulling apparatus. They concluded from their experiments that the higher the magnetic field applied, the lower the temperature fluctuation in the melt and encapsulant. Also, higher pressure in the chamber which caused a large temperature gradient in the melt and encapsulant and increased thermal convection, required a higher magnetic field to suppress the melt and encapsulant temperature fluctuation. Osaka et al [9] grew GaAs crystals using the LEC method in a growth apparatus designed for the study which used a two-heater configuration and controled the main heater power by a computer. They monitored the crystal diameter using a load cell that measured crystal weight during growth. In order to reduce the temperature gradient in the encapsulant, they used a thick layer of boron trioxide. They reported a gradient of 100°C/cm for a 2-cm thick layer and 35°C/cm for an 8-cm thick layer. They con-cluded from their experiments that the yield of single < 100 > oriented crystal reached nearly 100% when the solid/liquid interface profile was controlled such that it was convex Chapter 2. LITERATURE SURVEY 8 towards the melt. Kohda et al [10] grew dislocation-free and striation-free GaAs crystals using fully encapsulated Czochralski method. The crystal diameter was 50-mm. They used a 50-mm B203 layer which allowed full encapsulation of the pulled crystal.. This method reduced stress induced dislocations because arsenic evaporation, which caused surface irregularities, was effectively suppressed from the crystal surface. By using two heaters that were controlled individually, a temperature gradient of 30 — 50°C/cm was maintained along the growth direction in the encapsulant layer. Ozawa and Fukuda [11] observed the solid-liquid interface of GaAs grown by the LEC method. They reported that the region of stable crystal growth was indicated by a meniscus angle 8 in the range of 18° to 45°where 8 is defined as shown in figure 2.2. encapsulant crystal Figure 2.2: Meniscus angle 8 which must be between 18° and 45° for stable crystal growth. Kakimoto et al [12] observed directly the convection in boron trioxide during GaAs crystal growth using x-ray radiography. They studied natural convection in B2O3 with and without crystal in order to seperate the effect of forced convection attributed to Chapter 2. LITERATURE SURVEY 9 crystal and crucible rotation. For that case, they estimated the flow velocity to be about 0.3mm/sec which was two orders of magnitude less than that in the semiconductor melt. In the presence of a crystal, the flow velocity was about 0.03-0.05 mm/sec which was one order of magnitude smaller than that without a crystal. This reduction in flow velocity was attributed to shear viscous flow near the crystal wall. The results of the numerical simulation of the flow in B2O3 agreed very well with the experimental results (0.1-0.5 mm/sec). Their exercise showed that the flow in the encapsulant was steady and axisymmetric. Also, since the flow velocity was two orders of magnitude less than that in the semiconductor melt, conduction may be considered to dominate convection in the encapsulant B2O3. 2.3 Theoretical Studies Ramachandran and Duduhovic [13] used a detailed radiative heat transfer model (Geb-hart model) when simulating the temperature distribution in crystals grown by the Czochralski method. They also used a simplified model (Stefan's model) for the pur-pose of comparison. However, they assumed the melt, exposed crucible wall and the enclosing walls to be at uniform temperatures, i.e., they did not divide those surfaces to elements in order to take into account each surface's temperature gradient. By doing that they reduced the computing time required, however they also lost a more accurate representation of the real system. Atherton et al [14] studied the effect of diffuse-gray radiation on the parametric sen-sitivity and stability of the Czochralski process for growing silicon. They used a detailed radiative model which used the Gebhart method for computing radiative fluxes in diffuse-gray enclosures. For calculating the view factors, they used the programme FACET which was developed at Lawrence Livermore Laboratory to calculate view factors in an arbitrary Chapter 2. LITERATURE SURVEY 10 axisymmetric geometry. They calculated view factors for only four surfaces — crystal, crucible, melt and an enclosure that represented the ambient surface which they divided to 20 subsurfaces. They concluded that the importance of including detailed radiation modelling of the entire growth process was emphasized by the simulations which took into account the decreasing melt volume. Work by Srivastava et al [15] gave similar results as those given in [11] regarding the importance of the meniscus angle. Srivastava et al studied the liquid/solid interface based on a conductive model (conduction being the main mode of heat transfer in the melt and crystal). They used the Gebhart method for calculating radiative heat transfer in the puller. They divided the crystal top, crystal surface and melt surface to Ni elements, and assumed one surface for the crucible wall and one for the enclosing wall. They also assumed the crystal top to lose heat only to the enclosure wall by radiation. They studied the effects of several factors on the interface shape and pulling rate. They found that the effect of the meniscus shape on the interface shape was significant; its incorporation in the solution resulted in an increased concavity of the interface (concavity on the crystal side). The concavity was caused by the reduction in the axial heat transfer from the melt to the interface. Motakef and Witt [16] studied the effect of liquid encapsulation on the thermal stresses induced in the crystal. They found that stresses are strongly effected by the thermal transparency of the encapsulant. Derby and Brown [17] used a conduction dominated heat transfer model to form a simulation of the LEC crystal growth. They considered the extreme cases where the encapsulant was treated as either opaque or transparent. However, neither side of that treatment represented the actual case, since B2Oz was neither opaque nor fully transpar-ent to the energy transferred from the melt. They found that in the case of a transparent Chapter 2. LITERATURE SURVEY 11 encapsulant the direction of heat transfer in the encapsulant was horizontal (radially in-ward) from the wall to the crystal, while in the opaque case it was vertical - from the melt to the ambient (vertical isotherms versus horizontal isotherms). The melt/solid interface shape changed from concave into the crystal for a transparent layer to fully convex from the crystal for an opaque layer. The results of their work indicated that the optical prop-erties of boron trioxide are of great importance in determining the amount of energy that was transmitted from its surface. However, these optical properties were not determind in the temperature range of interest until recently and only by one source. It is worth mentioning here that even though B2O3 is transparent to radiation in the visible region, in the temperature range of interest here, only 1 % of the energy emitted lies within the visible region — assuming blackbody radiation at the GaAs melting temperature. Hicks [18] studied the fluid motion in the encapsulant assuming a temperature depen-dent viscosity and constant thermal conductivity. He introduced the Steffan-Boltzman radiation condition at the encapsulant/ambient interface and assumed a temperature distribution on the crystal-encapsulant, crucible-encapsulant and melt-encapsulant in-terface. He concluded from his results that the flow was predominantly in the azimuthal direction — caused by crystal and crucible rotation — while the flow in the meridional direction was negligible. He also concluded that heat flow was upward and conduction predominated convection in the encapsulant and thus the convective effects could be negelected. The latter result contradicted what Jordan [19] suggested that heat transfer between the crystal and encapsulant was strongly influenced by convective heat transfer within the encapsulant. Thomas et al [20] extended their thermal capillary model to include a low-pressure LEC system for growth of GaAs in an axial magnetic field that was strong enough that convective heat transport became unimportant. They assumed the encapsulant layer to be semitransparent using the data obtained from Ostrogorsky et al [21] to determine Chapter 2. LITERATURE SURVEY 12 that condition. (Semitransparent means that the layer is transparent to radiation be-low a certain wavelength and opaque to radiation above that wavelength.) They also assumed the crystal, crucible, melt and surrounding surfaces to be diffuse gray surfaces for radiative heat transfer calculations. They considered the surrounding surfaces to be at a constant temperature — one surface only. They used Gebhart's method for diffuse gray surfaces to calculate the radiative heat exchange between the surfaces. They found that the crystal shape predicted by assuming a semitransparent encapsulant, different emissivities for the melt, crystal and encapsulated crystal, and an enhanced encapsulant conductivity was the closest to that produced experimentally. Schvezov et al [22] calculated the temperature distribution in the crystal using a finite element model of the heat flow in the LEC growth system. They compared their results with an analytical solution of a simplified similar problem [23] and also with experi-mental results obtained from [24]. They found good agreement between their numerical and experimental results with the numerical results underpredicting the temperature distribution. They also compared the crystal axial temperature gradient and found the measured gradients to be significantly lower than those calculated near the interface, and above the calculated values away from the interface. Sabhapathy and Salcudean [25] studied numerically the effects of melt height and crystal radius on the melt natural convection in LEC growth of GaAs. They also analysed the effect of top and bottom crucible heating on the isotherms. They accounted for top heating by increasing the ambient temperature to which the system is exposed. They used a simplified radiation model in which all surfaces radiate to the same ambient temperature. They investigated the effect of melt height and crystal radius on the critical crucible temperature Tc(above which the crystal starts melting and below which the melt starts freezing at the melt-encapsulant interface). They concluded that by a proper combination of top and bottom heating, the crucible wall temperature could be lowered Chapter 2. LITERATURE SURVEY 13 and nearly flat isotherms at the crystal-melt interface could be obtained. Salcudean et al [26] studied free and forced convection in GaAs melt during the LEC crystal growth. They studied free convection due to melt heating at the crucible wall and forced convection due to crystal and crucible rotation. They found the isotherms in the forced convection case to be nearly equally spaced concluding that conduction dom-inates forced convection (maximum melt velocity in the vertical plane was O.OOAm/sec). However for the free convection case, the isotherms were concentrated near the wall and near the crystal where convective effects were the strongest. For the case of combined forced and free convection, they found the flow to be significantly different from that of the free or forced alone. They concluded that the flow was multicellular with some cells rotating faster than the crucible, and that the flow and temperature fields in the melt beneath the crystal were oscillatory. Sabhapathy and Salcudean [27] studied numerically fluid flow and heat transfer in the melt and encapsulant of GaAs LEC crystal puller. They assumed the encapsulant to be semitransparent to radiation. They also assumed the melt, crystal and encapsulant to be radiating to the same ambient temperature. Thus they used a simplified radiation model. They concluded for the crystal and crucible dimensions they used that the flow in the melt (when natural and forced convection were included) was often multicellular and effects the temperature field significantly. The melt and the encapsulant were much cooler than in the pure conduction case for the same heat input. They also reported that the isotherms in the melt close to the crystal-melt interface were nearly flat when the crystal rotated faster than the crucible. They stated that the flow in the encapsulant was similar to a couette flow between two rotating concentric cylinders, and that the crystal axial temperature gradient increased strongly with a decrease in the ambient temperature. Chapter 2. LITERATURE SURVEY 14 2.4 Scope and Objectives of the Present Work It is clear from the literature survey that the radiative exchange during GaAs crystal pulling is an important factor that must be included when modelling the process. How-ever, the assumptions and simplifications used by the authors influenced their results. The treatment of the encapsulant layer is one area that introduces errors. Another is the assumption that all surfaces radiate to an ambient at a constant temperature. The main objectives of the present work can be listed as follows: • To calculate the radiative exchange including all surfaces. • To obtain the best representative values of the material thermophysical properties. • To evaluate the importance of including radiative exchange of all surfaces versus replacing some of the surfaces by one isothermal surface. • To include the radiative model in the programme that calculates convection and conduction in the various elements of the crystal puller. • To evaluate the difference between the simplified and detailed radiative models in the crystal pulling modelling. In order to achieve those objectives, the radiative view factors which are needed for calculating the Gebhart factors first and heat transfer second, must be calculated. The Gebhart method for calculating heat transfer is a detailed one which includes the effect of multiple reflections in the calculation of heat transfer. A computer code is developed by the author to calculate all the required variables. The programme that models fluid flow and heat transfer in the melt and encapsulant and heat transfer in the crystal has been developed by Salcudean and Sabhapathy [25, 26, 27]. The programme solves the mass, momentum and energy equations using a Chapter 2. LITERATURE SURVEY 15 finite-difference control-volume method. The method i6 described fully in the text. The author's code is included with Salcudean and Sabhapathy's programme and a comparison between the velocity and temperature fields obtained using the simplified radiative model and the detailed Gebhart model is presented. Chapter 3 Radiative Heat Exchange in L E C GaAs Crystal Pullers 3.1 Introduction The position and the shape of the crystal/melt interface during the crystal pulling process and the temperature distribution in the crystal are important factors determining the crystal quality. Also, adverse thermal gradients in the melt may result in solid/liquid interface instabilities or in excessive thermal stresses leading to a poor quality crystal. Therefore, the melt temperature profile, which is coupled with the flow field due to the interaction of natural and forced convection and thermocapillary flow, is an important element in the crystal growth. The convection in the melt is also affected by radiation from its surface. The main objective of this part of the thesis is to determine qualitatively and quan-titatively the radiative heat transfer in the crystal puller. Thus, system geometry and radiative properties are of great importance. 3.2 Crystal Puller Geometry Assumptions The real system will be simplified for the purpose of the present study. However, the simplifications made will be such that the results can still be representative of the original system. The following assumptions are made: 1. The seed rod and the quartz crystal are not included as part of the system since their surface areas are small compared to the middle cylinder area to which they 16 Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 17 will mostly radiate. The seed rod - and quartz crystal - middle cylinder radii ratios are 1:12 and 1:6, respectively. 2. All surfaces are assumed to be diffuse gray ( emissivity is equal to the absorptivity, e = a, and the reflectivity p = 1.0 — a = 1.0 — e). A diffuse gray surface absorbs a fixed fraction of incident radiation from any direction and at any wave length. It emits radiation that is a fixed fraction of blackbody radiation for all directions and wavelengths. 3. The melt, encapsulant and crystal surfaces are assumed to be flat. This assumption is acceptable when calculating radiative view factors and radiative heat transfer. 4. The crucible and crystal are assumed to be cylindrical in shape of radii Rc and R,, respectively. 3.3 Mathematical Modelling 3.3.1 Configuration Factor Calculation The configuration factor is defined as the fraction of the radiation leaving one surface that reaches another surface. The calculations of radiative heat transfer in any system require evaluating the configuration factors of all the system surfaces. In the case of a GaAs crystal puller, the system is composed of the following surfaces: 1. Melt (m). 2. Encapsulant (e). 3. Crucible wall (c). 4. Crystal wall (s). Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 18 5. Middle cylinder wall (mew). 6. Middle cylinder lower annulus (mcla). 7. Crucible outer wall (co). 8. Crystal top (st). 9. Middle cylinder top annulus (meta). 10. Top cylinder wall (tew). 11. Top cylinder top (tct). Surfaces 1-8 will have to be divided into elements each at a uniform temperature for the purpose of configuration factor and later heat transfer calculations. The reason for that is the requirement that configuration factors be calculated between surfaces of uniform temperatures. Since it is desirable to simplify the system as much as possible without losing the heat transfer characteristics of the actual crystal puller, a comparison will be made between the heat transfer calculations obtained before and after simplifying the enclosure of the growth chamber (surfaces 5, 6, 9, 10 and 11) down to one surface. The equations used for calculating the configuration factors between all surfaces are given in appendix A. The reciprocity relationship is UBed for surfaces that are not related by direct equations. 3.3.2 Description of the Analysed Cases The following section describes the cases that are considered in the full chamber and simplified chamber configurations : • System consisting of the melt, encapsulant and upper chamber. Configuration factors and heat transfer from the melt surface first and from the encapsulant Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 19 surface second are calculated . The reason for anlysing this case is to determine the radiative heat transfer that is taking place in the system during the seeding process which is a very crucial step of the crystal growth process. If too much heat is applied then the seed will melt, while not enough heat will cause the melt to freeze on the seed at a fast rate without allowing for good diameter control. Figure 3.3 shows a schematic of the full chamber including the elemental division of the surfaces, while figure 3.4 (a) shows the simplified chamber system. • System consisting of the crystal growing at the rate of lcm/hr past the encapsulant top taking into account the dropping melt level but without allowing the crucible to protrude into the upper chamber. All other surfaces included are as in the previous case. For this case, the crystal is growing at the indicated rate and thus the configuration factors are calculated for crystal height of 3cm, 4cm, 5cm... until the crystal height reaches the crucible edge. Since this growth process is a batch process, the melt height decreases thus exposing more of the crucible wall. The number of surfaces for this case is increased from the previous case by an amount equal to the number of elements on the crystal wall and crystal top. Figure 3.5 shows a schematic of the system and figure 3.4 (b) shows the simplified system. • System consisting of melt, encapsulant, crystal and upper chamber with crystal extending above the crucible edge and the crucible protruding through the upper chamber. The number of surfaces here is increased from the previous case by the number of the elements on the outer surface of the crucible wall. The crucible protrudes through the upper chamber so that the centre of the heating goes through the centre of the melt. Figure 3.6 shows a schematic of the surfaces and figure 3.4 (c) shows the simplified chamber. Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 20 TCT TCW an MCTA -BCL-A MCW B.4. .aJ3. &2. .aJ. hm .hi. GaAs melt MCLA Encapsulant Figure 3.3: Schematic showing the system when there is no crystal (case-I). Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 21 Enclosure ..a or. .aJ. .hi. GaAs melt Enclosure Encapsulant hm Crystal .a2.... jil... .bl... GaAs melt (a) (b) ; Enclosure Crystal ...an. a/izl.. . BJ. ... hm .Al... (c) Figure 3.4: Schematic showing the simplified chamber system: (a)case-I, (b)case-II and (c)case-III. Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 22 TCT TCW an MCTA jan-1. .anr2 MCW hm Crystal .8.4 ...&3.... a2.. ...BL GaAs melt MCLA Encapsulant Figure 3.5: Schematic showing the surfaces under consideration for crystal height any-where between the encapsulant top and the crucible edge (case-II). Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 23 TCT TCW Bn MCTA .aiirl an-2 MCW GaAs melt prystah Ytxl. ...3,4 -.a.3. do... .a2\...dl ...al MCLA Encapsulant Figure 3.6: Schematic showing the surfaces for the crystal growing past the crucible edge (case-Ill). Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 24 The fictitious surfaces 'a', 'b' and'd' shown in figures 3.3 to 3.6 are useful in evaluating the configuration factors from some surfaces to others. A description of configuration factor calculations for some surfaces is given in appendix A. The summation of the configuration factors from each surface to all other surfaces must be equal to one. Any deviation from one represents an error. Appendix D gives the results of configuration factor calculations for the simplified system described above. 3.4 Heat Transfer Model The Gebhart method of determining radiative exchange between surfaces of an enlosure is used here. The method can be used whether heat fluxes or temperature distributions are specified. The general form of the heat equation from [43] obtained by a heat balance on surface k is : Qk = AfcCfccrT/ - X) A&aTfGjk (3.1) i=i where Gjk is the Gebhart factor which is denned as the fraction of the emission from surface Aj that reaches A/, and is absorbed. This includes all the paths for reaching Ak, that is, the direct path, paths by means of one reflection, and paths by means of multiple reflections. In equation form, Gjk is given by: Gjk = Fj-k^k + Fj-ipiG\-k + Fj-2p2G2-k + + Fj-kpkGk-k + •••• + Fj-NpNGN-k (3.2) Appendix B gives the derivation of the Gebhart factor equations for the system surfaces. It also gives the system of equations in matrix form and the method of solution of the Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 25 matrices forming the problem. The flow chart of the programme developed to solve the matrices and calculate the radiative heat transfer is given in appendix C. The net radiation method which was devised by Hottel [40] and later developed by Poljak [41, 42] and the Gebhart method are all equivalent and give the same results. Gebhart's method was chosen in this thesis. The configuration factor calculations were verified by ascertaining that the summation of those factors for each elemental surface is equal to 1. 3.4.1 Diffuse-Gray Surfaces This study will assume all surfaces to be diffuse-gray. According to [43], diffuse signifies that the emissivity and absorptivity do not depend on direction, and gray signifies that the spectral emissivity and spectral absorptivity do not depend on wave length. Thus, a diffuse-gray surface is one which absorbs a fixed fraction of the incident radiation from any direction and at any wavelength, and it emits a fixed fraction of blackbody radiation for all directions and wavelengths. This assumption is in-line with the data available for the emissivity of the surfaces, and until more data is available it will have to be used. As will be found in the section on material thermophysical properties, the data given is of the total normal emissivity of the surfaces as a function of surface temperature. 3.4.2 Encapsulant Layer Semitransparency The absorption of the B2O3 encapsulant layer is very important for the determination of radiative heat transfer in the system. Only one source of experimental data is available (Ostrogorsky et al [21].) Tables 3.1 and 3.2 give the absorbance and absorption coefficient of B2O3 for three layer thicknesses : 1.0cm, 0.4cm and 0.2cm. The absorbance in table 3.1 is obtained from figure 10 of [21] — shown here in figure 3.7. While the absorption coefficient given in table 3.2 is calculated using equation 3.3. Ostrogorsky et al state Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 26 that the indicated layer thicknesses are approximate because they were computed from the removal of weighed quantities of B2O3. They also state that the non-planarity of the melt surfaces represents another uncertainty. Also, the loss of radiation intensity due to scattering of the radiation at the concave melt surface is increased by the long optical path between the sample and detector. Therefore, the uncertainty in this data should be taken into account when choosing a cutoff wavelength for the encapsulant (above which the encapsulant is opaque and below which it is transparent). The absorption coefficient a is calculated by substituting the values of L and A in equation-4 of the same reference (repeated here for convenience), ABIO, = 0.434ai (3.3) where A is the absorbance, a is the absorption coefficient and L is the layer thickness. Also, the optical transmittance r is given by : T = m = e z p ( _ a a j ) (3-4) where a is as defined earlier and x is a position along the path that the radiation is travelling. An average thickness of the encapsulant layer during the growth of GaAs crystals is about 2-cm. However, the data that is given is for a layer that is 1cm thick or less. Therefore, extrapolation from the available data is required (another factor of uncertainty). The transmittance at A = 2.0 and 2.5fim is calculated using equation 3.4. The results are given in table 3.3. It is clear that as the layer thickness increases, r decreases. Also, the values of r are lower for A = 2.5fim than they are for A = 2.0fim. Thus indicating that for layers thicker than 1cm, and for A greater than 2.5 fim, the encapsulant is essentially opaque. Ostrogorsky et al concluded from their data the following: Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 27 WAVELENGTH [vim] 2 2 . 5 3 4 5 ,i i ' » , , 2 1 / ' 1 1 ; 3 / 1 II!' ' ii! il ii il j l • J j i : J hJ • 6000 4000 2000 KAVENUMBER [cm' 1 ] Figure 3.7: Absorbance spectra of B203 at 1250°C for L=0.2cm (1), L=0.4cm (2), and L=1.0cm (3). From Ostrogorsky et al. Table 3.1: B2O3 absorbance for three layer thicknesses obtained from figure-3.7. Wave Length A(/zm) Absorbance 'A' L = 1.0cm L = 0.4cm L = 0.2cm 2.0 0.081 0.032 0.00 2.1 0.10 0.036 0.00 2.2 0.14 0.064 0.00 2.3 0.23 0.097 0.00 2.4 0.48 0.19 0.032 2.5 2.00 0.71 0.29 Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 28 Table 3.2: The absorption coefficient of B2O3 for three different layer thicknesses. Wave Length X(fim) Absorption Coefficient ' a ' L = 1.0cm L = 0.4cm L = 0.2cm 2.0 0.19 0.19 0.00 2.1 0.24 0.20 0.00 2.2 0.34 0.37 0.00 2.3 0.52 0.56 0.00 2.4 1.12 1.12 0.37 2.5 4.61 4.09 3.34 • B2O3 is transparent for A < 2pm. • B2O3 is either transparent or opaque depending on the layer thickness for 1.9pm < A < 2.8/xm. • B2O3 is opaque for A > 2.8um. Since only one source of data is available, and the source describes many areas of uncer-tainty, it seems impractical and unnecessary to divide the wave-spectrum to more than two regions, thus the decision to have one cutoff wavelength. For this work, the cutoff wavelength will be chosen as 2.0 pm. However, calculations using a cutoff wave length of 2.5 pm will be done for comparison. The wavelength spectrum is then divided into two regions in order to consider the absorption property of the encapsulant. The first region is below 2pm where radiative exchange with the melt surface is considered. The second region is above 2pm where radiative exchange with the encapsulant surface is considered. Under that assumption, the radiative heat flux for any surface segment k will be given by: Qk = Akekf'k<rT} - £ A ^ c T J G ^ ^ for A < 2 (3.5) i=i Qk = Akek(l - fk)oT} - £ A i C j ( l - f^T}G^ue for A > 2 (3.6) Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 29 where At, is the surface area, e* is the emissivity of the surface, Gjk is the Gebhart factor and f'k is the fraction of energy emitted at 7* below 2p.m. The heat transfer from each surface is summed up and an effective ambient temper-ature can be calculated as follows : ft- = ( T T E Ajejf'jT}G%an^r'nt + f; Ajej(l - f'j)T}G%>ny/A (3-7) This effective temperature can be compared to a constant ambient temperature that the system would be radiating to if a simplified radiative model is used (Stefan's model). It is also used to produce the correct radiative heat flux from the surface when applying the radiative boundary condition during a complete analysis of heat transfer in the system (including conduction in the encapsulant and crystal, and convection in the melt and surrounding gas). 3.4.3 Materials' Thermophysical Properties The thermophysical properties of the semiconductor GaAs, the encapsulant B2O3, the surrounding gas Ar, and the other chamber surfaces that form the crystal puller are important factors when calculating the fluid flow and heat tranfer in the puller during the crystal growth process. Jordan [19, 23, 32] presents detailed studies of the properties of GaAs and B2O3. He also gives the convective heat transfer coefficient between the encapsulant and crystal, and between the high pressure gas, Ar, and the encapsulant and crystal. The radiative heat transfer coefficient from the GaAs surface is also included. All equations that describe properties that are available from Jordan are included here for the sake of completeness. The equations that describe the heat transfer coefficients are also included. The temperature dependent equations of GaAs density (p'), specific heat (Cp), con-ductivity (k) and thermal expansion coefficient (B) are listed in table 3.4. The radiative Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 30 heat transfer coefficient from the melt and crystal surfaces is also given in the table. While the temperature dependent density and viscosity of the encapsulant B2O3 which are available from [30] and temperature dependent conductivity which is obtained from [31] are tabulated in 3.5. The convective heat transfer coefficient between the encapsulant and the crystal is also given in the table. Table 3.6 gives the properties of high pressure Ar. The data is obtained from [33]. M is the molecular weight of argon, and R is the universal gas constant, p represents the gas pressure in atmospheres; it is introduced in the convective heat transfer equation as a correction factor to take into account the high gas pressure during the Czochralski crystal growth. Properties Used by Different Modellers Table 3.8 presents the physical properties of GaAs (melt and solid), the encapsulant and the chamber surfaces as used by various authors. Some of the references chose to use constant properties while others chose temperature dependent properties. Many references - all not listed in the table, but mentioned in the literature survey - used the values given by Jordan in [19, 23]. It is clear from the table that there is a wide range of values used for the emissivity of GaAs liquid and solid (0.3, 0.5, 0.55, 0.7 and 0.75). The values used for the emissivity of the encapsulant layer also varied from one source to another (0, 0.65, 0.75 and 1). Surface Emissivities Other sources [33, 34, 35] were researched for material properties especially for the emissivity of the surfaces. The emissivity is very important when calculating radiative heat transfer. Therefore, it is important to use the best representative values. The Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 31 Table 3.3: Transmittance for wave lengths of 2 and 2.5 pm for each layer thickness L = 0.2, 0.4 and 1.0cm. X(pm) Transmittance V L=0.2cm L=0.4cm L=1.0cm 2.0 1.00 0.93 0.83 2.5 0.51 0.20 0.01 Table 3.4: Temperature dependent properties of GaAs. PROPERTY EQUATION Density p'(gm/cms) 5.32 - 9.91 * 10" 5 r Specific Heat Cp(J/gm.K) 0.302 + 8.1 *10" BT Conductivity k{W/m.K) 20800T-1-09 Therm. Exp. Coef 0{K-*) 4.68 * 10"6 + 3.82 * 10" 9r h-radiation (cm-1) 2 . 2 7 * 1 0 - 1 1 , rr3 k etIGaA» Table 3.5: Temperature dependent properties of the encapsulant B2O3. PROPERTY EQUATION Density p'(gm/cm3) 2.705 - 2.814 * 10"3T + 2.677 * lO'^T2 -1.204 * 10"9T3 + 2.102 * i o - 1 3 r 4 Them. Exp. Coef. 1 8£ p> BT Dyn. Viscosity p(poise) logioP = -4.343 + 17040T-1 - 1.696 * 10 7T- 2 +7.081 * io 9r~ 3 Conductivity k(W/cm.K) 2.37 *10"3 + 1.1 * i o - 5 r Specific Heat Cp(J/gm.K) 1.83 hconvection (cm-1) ( £ Z ^ ) I / 4 0 . 5 4 8 ^ ( 5 ( - P ' | § ; ^ ) B 3 O 3 ) 1 / 4 Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 32 Table 3.6: Thermophysical properties of high pressure argon as obtained from Touloukian. PROPERTY EQUATION Density p'(gm/cm3) PATKRTA*) Therm. Exp. coef. l/TAr Dyn. Viscosity p,(poise) see table 3.7 Conductivity k(W/cm.K) see table 3.7 Specific Heat Cp{J/gm.K) 0.5204 hconvection (cm-1) ( 1 ^ ) 1 / 4 » o . 5 4 8 ^ ) 1 / a ( ^ ( & ) * ) 1 / 4 Table 3.7: Viscosity and conductivity of .Ar ignoring pressure dependence (from Touloukian). Temp (K) (N.sec/m2) k (W/cm.K) Temp (K) (N.sec/m2) k (W/cm.K) 600 38.2 0.301 * 10"3 760 44.9 0.356 * 10"3 620 39.1 0.308 780 45.7 0.362 640 40.0 0.315 800 46.4 0.369 660 40.9 0.322 820 47.2 0.375 680 41.7 0.329 840 47.9 0.381 700 42.5 0.336 860 48.6 0.387 720 43.3 0.343 . 880 49.4 0.393 740 44.1 0.349 900 50.1 0.398 Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 33 Table 3.8: Materials' properties as reported in the references cited. AUTHOR SURFACE PROPERTY P' (gm/cm3) fi k (W/m.K) Cr (J/g.k) e Kakimoto B203 1.51 3 * 10"B 2.0 1.836 N/A Hicks ° B203 1.52 N/A 1.8 1.830 0.75 Motakef & Witt c B203 N/A N/A see 6 N/A N/A GaAs N/A see d see e N/A 0.5 Derby & Brown GaASmeit 5.71 N/A 14 0.42 0.55 GaAs,oiid 5.17 — 7 0.42 0.55 B203 1.51 N/A 2 N/A 0 & 1 Graphite 1.6 — 0.6 2.1 0.8 Si02liner 2.2 — 0.06 1.3 0.35 Crochet et al GaAsmeit 5.63 1.89 *10"4 15.15 0.507 N/A GaAstoud 5.71 — 15.15 0.507 N/A Meduoye et al / GaAs N/A 1 * 10"B 7.74 N/A N/A Thomas et al GaAsmeit 5.71 N/A 14 0.42 0.3 GaAsnoiij 5.17 — 7 0.42 0.55 & 0.7 B203 1.51 N/A 2 & 8 N/A N/A, A < 2 0.65, A > 2 Dupret et al Steel N/A — 33 N/A 0.6 Graphite N/A — 41.9 N/A 0.81 Gr. Felt N/A — 2.5 N/A N/A Sabhapathy k Salcudean GaAsmeit 5.71 1.87 * 10"4 17.8 0.434 0.75 GaAs goud 5.2 — 7 0.42 0.75 B203 1.5 5 * i<r6 1.85 1.83 0.75 aTcrucible = 1530Ji:,Tomi = 1200JT '2.37 * lO" 1 + 1.1* 10- 3 T eTAr = 700.R-, T„ o J , = 400*, (Bod#). = 0.5, (Rad#)B,0, = 1.0 • k ^ S * 10"6 + 3.82 * 10~ 9T J208 + T - 1 - 0 9 * used temperature dependent hconvection that is available from Jordan. Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 34 emissivity of the melt (liquid GaAs ), the crystal (solid GaAs ), the crucible wall (PBN pyrolytic boron nitride), and the chamber walls (stainless steel) are all required. Only one source [32] gives data for the emittance of n-type GaAs as a function of temperature and the product C<jt, where d is the doping level or impurity concentration in c m - 3 and tis the crystal diameter in cm. (see figure 3.8). The data given in figure 3.8 was obtained by numerical integration of the spectral emittance between the limits 0.5 and 25 pm in 0.1 pm intervals. It is clear from the figure that the total emittance increases as the product dt increases. For the temperature range of interest during crystal growth (1000-1500) K, and for a crystal diameter of 8cm and doping level of l*1016cm~3, the total emittance et ranges from (0.61-0.66). As the impurity concentration increases for the same temperature range, the total emittance varies from (0.65-0.69). With these results for the crystal diameter and doping level chosen, it is not necessary to use temperature dependent emissivity for GaAsloud. Thus, the calculations presented in this work use et = 0.65 for GaAstoud. Lemons and Bosch [38] state in their work that the emissivity of molten silicon is substantially lower when compared to that of solid silicon. The emissivity of the GaAs melt has not been evaluated experimentally but is expected to be different than that of the GaAs solid - in line with silicon. Thomas et al [20] used a value of 0.3 for the emissivity of molten GaAs. With lack of better information, the same value will be used here. For p-type GaAs, Braunstein and Kane [39] determined the absorption coefficient between 2 and 25 pm. The result of their work states that the absorption coefficient increases with wavelength and temperature. Because of the high absorption coefficient, the total emittance is related only to the reflectivity. Therefore, the curve for G^t > 5 * 1017 serves as a reasonable estimate of ej for p-type GaAs. Reference [33] gives data for the normal spectral reflectance and transmittance of Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 35 0.11—i—I 1 i 900 100 TOO tOO 1100 1500 1300 TEMPERATURE (feg-K) Figure 3.8: The emissivity of GaAs as function of temperature and doping density (from Jordan). Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 36 GaAs. However, the data is for room temperature only, which means that it is of no use for calculations of the crystal growth process. All other books on the thermophysical properties of materials surveyed [34, 35] do not provide any information on the emissivity of GaAs. The emissivity of the encapsulant B20z is not available in any of the references researched. Since B20z is basically glass, the emissivity of glass at high temperatures is obtained from the Handbook of Chemistry and Physics and from [33]. Both references give a value of 0.7 for the emissivity of high temperature glass. Therefore, this study will use 0.7 as a representative value of the emissivity of the encapsulant. The emittance of stainless steel is available from many sources. Touloukian et al [33] present an excellent collection of data for oxidized, polished and cleaned stainless steels for a wide range of temperature (320 - 1300)K. The graphs and data are available in appendix E. For the crystal puller chamber, cleaned stainless steel is assumed in the temperature range (450 - 800)K for which the normal total emittance range is (0.2 -0.42). 3.5 Temperature Distribution of the Surfaces The temperature distributions that are used to calcualte the heat transfer from each surface are based on results obtained from numerical and experimental results [12, 22, 26, 27, 25]. The values obtained are analysed to produce a temperature gradient across the surface which is used to evaluate the temperature of each element that the surface is divided to. The temperature equations used are given in table 3.9. The table gives also the temperature gradient across each surface. The gradient is the slope of the temperature equation. 0 Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 37 Table 3.9: Temperature distribution of the surfaces obtained from numerical solutions for cases I, II and III. SURFACE TEMPERATURE DISTRIBUTION (K) (radii and heights in m) TEMPERATURE GRADIENT (K/cm) Melt (I) r m ( t ) = 1511 + 656 * (rm(i) + r m ( t - l))/2 7 Melt (II & III) Tm(i) = 1511 + ' « W - r - W *(rm{i) + rm(i-l))/2' _ Encapsulant (I) Te(i) = 1320 + 2625 * (r.(») + r e ( t - l))/2 26 Encapsulant (II and III) r.(*) = Tc(e) + 2ilfi=2*l * ( r e ( » ) + re(i - l))/2 ' Crucible (I k II) re(i) = 1565 - 1312 * (h'e(i) + h'c(i - l))/2 13 Crucible (III) Tc[i) = 1530 - 2500 * (h'c(i) + h'c(i - l))/2 25 Mid Cylinder Wall (MCW) Tmc(i) = 665 - 671 * (h'mc(i) + h'mc(i - l))/2 7 MCLA Tia{i) = 1465 - 5249 * (r,a(») + rla{i - l))/2 52 Crystal T,{i) = 1511 - 3445 * {h',{i) + h't(i - l))/2 34 Curcible Outer Wall Teo(i) = 1565 - 2500* ((K0(i)+Ko(i-W2 + (H'c-H'co)) 25 Crystal top Ttt = T.(ns) - 787 * (r,t(i) + r,t{i - l))/2 8 MCTA Tmcta = 415 0 TCW Ttcw = 315 0 TCT Ttct = 315 0 Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 38 3.6 Analysis of Results 3.6.1 Five Surfaces Enclosure The results that are given here are obtained using crystal puller dimensions and melt and surface properties that are given in table 3.10. These results are dependent on the temperature distributions which are the best possible estimates. However, some general conclusions can be made regarding the assumptions that are used when analysing global heat transfer in the system. The results of the heat transfer calculations are given in tables 3.11 through 3.15. The calculated effective ambient temperatures for each surface element are plotted versus element number in figures 3.9 to 3.11. The analysis of these results is given next for each case. 3.6.2 Analysis of the results for case-I Table 3.11 gives the heat transfer from each surface to other surfaces in the enclosure. It also gives the total heat absorbed and emitted and the net heat transferred from each surface. The net exchanged in the system is 4.2kW. The melt, crucible, encapsulant and the middle cylinder lower annulus lose heat, while the middle cylinder wall and the other three surfaces absorb the heat lost from the other surfaces. The following notes describe the exchange of each surface: • The melt loses heat mainly to the crucible wall and the middle cylinder wall. It absorbs heat from the crucible. • The crucible exchanges heat mainly with the encapsulant. It loses six times as much heat to the middle cylinder wall as it does to the melt. • The middle cylinder wall absorbs heat from the melt, encapsulant, lower annulus and its surface at the same order of magnitude. However, it absorbs heat from the Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 39 Table 3.10: Dimensions and surface properties of the crystal puller. DIMENSION or PROPERTY SURFACE Melt Encap Crucible Crystal MCW MCLA,-TA Top Cyl Diameter (mm) (inch) — — 152.4 6.0 76.2 3.0 304.8 12.0 (DiyD0) (Dc,Dmc) 152.4 6.0 Height (mm) (inch) 48 1.9 20 0.79 90-130 3.5-5.1 0-203 0-8.0 372.5 .14.67 — 101.6 4.0 e 0.3 0.70 0.5 0.64 0.3 0.3 0.3 Mass (Kg) 5.0 — — — — — — Density (Kg/m*) 5720 — — 5200 — — — Table 3.11: Heat transferred (in Watts) from surface i to surface j for case-I. SURFACE j SURFACE i Melt Crucible MCW MCLA Encap MCTA TCW TCT Melt 13 262 0 2 — 0 0 0 Crucible 294 2232 17 33 1043 1 0 0 MCW 144 1652 335 928 800 19 3 1 MCLA 19 176 70 118 111 2 0 0 Encapsulant — 1064 26 41 146 1 0 0 MCTA 25 199 48 115 125 3 0 0 TCW 12 99 26 74 61 1 3 1 TCT 7 43 11 33 35 0 1 0 Emitted (2) 514 5727 533 1344 2321 27 7 2 Absorbed (1) 277 3620 3882 496 1278 515 277 130 Net (2)-(l) 237 2107 -3349 848 1043 -488 -270 -127 Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 40 crucible that is one order of magnitude higher than the other surfaces. • The lower annulus exchanges heat mainly with the middle cylinder wall. It also absorbs, through multiple reflections, heat generating from its surface. • The encapsulant radiates mainly to the crucible and the middle cylinder wall. It absorbs heat from the crucible and, through multiple reflections, from its surface. • The heat exchanged with the middle cylinder top annulus, top cylinder wall and top is one to three orders of magnitude less than that exchanged between the other surfaces. Figure 3.9 shows the effective ambient temperature of the surface elements versus the element number. The plot shows clearly that assuming one constant ambient temperature for all surfaces - or all elements of a surface - is a poor assumption. The jump in the crucible curve at node-4 is due to the exposure at that node to the encapsulant surface where a higher radiative exchange between the surfaces takes place. The jump in the middle cylinder wall curve is due to the higher exposure of node-2 to the melt and encapsulant than node-1. Eventhough the melt surface temperature is higher than the encapsulant surface temperature, its effective ambient temperature is lower since only the percentage of the heat radiated below 2pm reaches the melt from the crystal, crucible and other surfaces. The effective ambient temperature range for this case is (750-1350) K. 3.6.3 Analysis of the results for case—II The results obtained for this case are calculated for a crystal height of 61-mm (2.4-in). The numerical values of the heat transfer between the surfaces are given in tables 3.12 and 3.13. The net exchange for this case is 2.7kW. The melt, encapsulant, crucible, Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 1400 _ • — 1300 T E M E 1200 R A T U * 1100 1000 900 / .0 « -4 a — a M e l t B---B E n c a p s u l a n t 0 - - - 0 c r u c i b l e H a l 1 4 6 ELEMENT NUMBER 10 1050 T E M P E R A T U R E 1000 J 950 J 900 J 850 800 750 700 a — ^ M i d d l e C y l i n d e r W a l l Q---Q MCLA «---<> MCTA * - - - * TCW TCT 12 14 0 2 4 6 8 10 ELEMENT NUMBER Figure 3.9: Effective ambient temperature versus element number for case-I. Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 42 crystal and the middle cylinder lower annulus lose heat, while the middle cylinder wall and the other three surfaces absorb heat. The following notes describe the exchange of each surface. (The numbers between brackets represent the percentage of heat lost or gained by the surface.) • The melt exchanges heat mainly with the crucible and crystal walls (about 72% of total heat lost and 92% of total heat gained.) • The crucible wall exchanges heat with (loses to and absorbs from) the encapsulant, crystal, middle cylinder wall and itself. It loses heat to the middle cylinder wall. All quantities lost or absorbed are of the same order of magnitude, however, that exchanged with the melt and crystal top are 2 to 10 times lower. • The middle cylinder wall absorbs heat from the lower annulus (35%), crucible (22%), encapsulant (12%), crystal (6%), and from its surface(13%). • The middle cylinder lower annulus loses heat to the middle cylinder wall and top annulus. It absorbs heat from the middle cylinder wall and itself (through multiple reflections). The quantity that is lost to the middle cylinder wall is nine times higher than that lost or absorbed from the other surfaces. • The encapsulant exchanges heat with the crucible, crystal and it loses heat to the middle cylinder wall . The quantity of heat that is absorbed from the crucible is 4-10 times higher than that absorbed or lost to the other surfaces. • The crystal loses heat to the crucible wall, the encapsulant, the middle cylinder wall. It absorbs heat from the crucible and the encapsulant. The amount of heat that is exchanged with the crucible wall is 4-6 times higher than that exchanged with the other surfaces. Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 43 Table 3.12: Heat transferred (in Watts) from surface i to surface j for case-II / i represents major surfaces. SURFACEj SURFACE i Melt Crucible MCW MCLA Encapsulant Crystal Melt 17 157 0 0 — 58 Crucible 159 1204 12 23 527 710 MCW 58 583 334 926 328 153 MCLA 9 86 70 118 51 23 Encapsulant — 960 13 18 127 266 Crystal 82 952 4 8 252 154 Crystal Top 1 40 8 14 8 6 MCTA 4 39 47 115 23 10 TCW 2 22 26 74 13 6 TCT 1 10 11 33 6 3 Emitted (2) 333 4054 525 1330 1333 1389 Absorbed (1) 233 2662 2624 390 1391 1457 Net (2)-(l) 100 1392 -2099 940 -59 -68 Table 3.13: Heat transferred (in Watts) from surface i to surface j for case-II / i represents minor surfaces. SURFACE j SURFACE i Crystal Top MCTA TCW TCT Melt 1 0 0 0 Crucible 26 0 0 0 MCW 219 19 3 1 MCLA 32 2 1 0 Encapsulant 7 0 0 0 Crystal 5 0 0 0 Crystal top 6 0 0 0 MCTA 32 3 0 0 TCW 16 1 3 1 TCT 9 0 1 0 Emitted (2) 353 25 8 2 Absorbed (1) 83 273 164 74 Net (2)-(l) 270 -248 -156 -72 Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 44 • The crystal top loses heat mainly to the middle cylinder wall (62% of total heat lost). • The quantities exchanged between the top cylinder wall and top and all other surfaces are 1-5 orders of magnitude less than those exchanged between the other surfaces. Figure 3.10 shows the effective ambient temperature versus the element number for the system surfaces. The jump in the crucible wall and crystal curves show the effect of higher heat gain due to location above the encapsulant layer. Element 4 on the crystal and crucible surfaces is the first one above the encapsulant top. The general shape of the melt and encapsulant curves is the same - as expected. However, the encapsulant curve is higher than the melt curve - also as expected. Figure 3.10 also shows a jump at element 2 on the middle cylinder wall curve which is due to exposure of that element to the crucible, crystal and encapsulant more than element 1. The range of the equivilent ambient temperature in this case is (650-1400) K. 3.6.4 Analysis of the results for case-Ill Tables 3.14 and 3.15 give the heat transferred between all the surfaces. They also give the total heat absorbed and emitted and the net transferred from each surface. The net heat exchanged in this case is 3.9kW. The crucible, melt, crystal wall and top and middle cylinder lower annulus lose heat. While the encapsulant, middle cylinder wall and top annulus, the top cylinder wall and top absorb the heat lost by the other surfaces. The crystal height in this case is 20.3-cm (8.0-in). The following notes describe the exchange of each surface with the other surfaces. (The percentages in brackets represent percent of total heat that is lost or gained by the surface.) Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers T E M P E R A T U R E 1500 1400 1300 1200 1100 J 1000 900 . . Q -* — A M e l t • - - - Q E n c a p s u l a n t *---<> C r u c i b l e Wal 1 * - - - * C r y s t a l Wal 1 C r y s t a l Top 10 ELEMENT NUMBER 1050 T E M P E R A T U R E 1000 J 950 900 J 850 800 750 700 J 650 -o a — ^ M i d d l e C y l i n d e r W a l l MCLA o--o MCTA * - - - * TCW o—-o T C T 6 8 10 ELEMENT NUMBER 12 14 Figure 3.10: Effective ambient temperature versus element number for case-II. Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 46 Table 3.14: Heat transferred (in Watts) from surface i to surface j for case-Ill / i repre-sents major surfaces. SURFACEj SURFACE i Melt Crucible MCW MCLA Encap Crystal Crucible OW Melt 15 136 2 0 — 65 0 Crucible 184 1970 3 4 728 1827 10 MCW 7 262 290 593 50 685 1460 MCLA 1 42 58 82 6 87 486 Encapsulant — 1078 1 1 133 395 2 Crystal 94 2441 44 65 446 665 140 Crucible OW 1 39 63 276 5 78 346 Crystal Top 0 3 7 7 1 10 13 MCTA 0 18 32 37 5 56 60 TCW 0 19 18 50 3 28 51 TCT 0 7 9 24 1 12 25 Emitted (2) 304 6016 525 1140 1377 3908 2591 Absorbed (1) 217 4718 3465 771 1610 3902 814 Net (2)-(l) 87 1298 -2940 369 -233 • 6 1777 Table 3.15: Heat transferred (in Watts) from surface i to surface j for case-Ill / i repre-sents minor surfaces. SURFACEj SURFACE i Crystal Top MCTA TCW TCT Melt 0 0 0 0 Crucible 1 0 0 0 MCW 99 16 2 1 MCLA 6 2 0 0 Encapsulant 0 0 0 0 Crystal 6 1 0 0 Crucible OW 5 1 0 0 Crystal top 2 1 0 0 MCTA 20 2 0 0 TCW 10 1 3 1 TCT 6 0 1 0 Emitted (2) 155 24 6 2 Absorbed (1) 44 230 184 87 Net (2)-(l) 111 -206 -178 -85 er 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 47 The melt loses heat mainly to the crucible wall (61%) and the crystal (31%). It absorbs heat from the crucible wall (63%). The crucible loses heat to the encapsulant, the crystal and itself by an amount that is one order of magnitude higher than that to the melt or the middle cylinder wall. It absorbs heat from the crystal, the encapsulant and melt. The heat absorbed from the crystal is one order of magnitude higher than that absorbed from the melt. The middle cylinder wall radiates heat mainly to itself through multiple reflections (55%). It absorbes heat from the crucible (50%), the lower annulus and the crystal. The middle cylinder lower annulus radiates mainly to the middle cylinder wall (52%) and the crucible outer wall (24%). It absorbs heat from the crucible outer wall (63%). The encapsulant exchanges heat mainly with the crucible and the crystal (about 85% of heat lost and 82% of heat gained). The crystal radiates to the crucible (46%), encapsulant, middle cylinder wall and itself - through multiple reflections (44% of total heat lost). It absorbs heat from the crucible mainly (63%) and also from the encapsulant. The outer wall of the crucible radiates heat mainly to the middle cylinder wall and the lower annulus (75% of total heat lost). It absorbs heat from the lower annulus (34%) and itself (43%) through multiple reflections. The heat exchanged with the crystal top, middle cylinder top annulus, top cylinder wall and top is one to three orders of magnitude less than that exchanged between the other surfaces. Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 48 Figure 3.11 shows the effective ambient temperature versus the element number for case-Ill. The jumps in the crystal curve at nodes 4 and 18 are due to the increased heat exchange with the crucible wall at node-4 and with the middle cylinder wall at node-18. The general shape of the crucible, melt and encapsulant curves is the same as in the previous cases. The temperature range for these surfaces is (650-1400)K in this case. 3.6.5 One Surface Enclosure The middle cylinder wall, lower and top annuli, and the top cylinder wall and top are replaced by one surface that is referred to as the enclosure. Radiative heat transfer calculations are then performed for the three cases discussed earlier in order to compare the results of the full chamber and the simplified chamber. The numerical values of the heat transferred are given in tables 3.16 to 3.18. The plots of the effective ambient temperature versus the element number are given in figures 3.12 to 3.14. The following is a detailed analysis of each case. 3.6.6 Analysis of Case-I Table 3.16 gives the heat transfer results for case-I. The net heat exchanged in the system is 2.6kW. The crucible, melt and encapsulant lose heat to the enclosure which is at 600K. The following notes describe the exchange of each surface: • The melt loses heat to the crucible wall and the enclosure at the same order of magnitude. It absorbs heat from the crucible wall while not absorbing anything from the enclosure. • The encapsulant exchanges heat with the crucible. It loses heat to the enclosure and absorbs heat originating from its surface through multiple reflections. • The crucible exchanges comprable quantities with all surfaces (including itself). Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers . - B -a — * M e l t • • - - • E n c a p s u l a n t ^ - " • C r u c i b l e Wal 1 C r y s t a l Wal 1 C r y s t a l Top T~ 18 10 12 14 16 ELEMENT NUMBER 20 1^ 22 24 26 28 ELEMENT NUMBER Figure 3.11: Effective ambient temperature versus element number for case-Ill. Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 50 • The enclosure being at a lower temperature absorbs heat from all other surfaces. Figure 3.12 shows the effective ambient temperature versus element number for this case. It is clear that the general shape of the graph is the same as that for the full chamber case. However, the numerical values are (3-8) % higher for the simplified case. 3.6.7 Analysis of Case-II Table 3.17 gives the results of the radiative heat transfer calculations for a crystal height that is less than the crucible wall height (2.6-in). The net heat exchanged in this case is 2.1kW. The crucible, melt, crystal and encapsulant lose heat to the enclosure which is at 600K. The following notes can also be obtained by further examination of the table : • The melt exchanges heat mainly with the crucible and crystal (about 67% of total heat lost and 90% of total heat gained). • The crucible wall loses heat to the encapsulant, enclosure and crystal wall, and absorbs heat from the encapsulant, crystal wall and itself. All those exchanges are of the same order of magnitude. • The encapsulant loses heat to the crucible (36%), enclosure (30%) and crystal wall (16%). It absorbs heat from the crucible wall (65%), the crystal (16%) and itself (14% - through multiple reflections). The exchange between the encapsulant and the other surfaces is one to two orders of magnitude lower. • The crystal wall exchanges heat with the crucible (48% of total heat lost and 64% of total heat gained), encapsulant (18% lost and 17% gained) and enclosure (18% lost and 11% gained). The exchange with the melt is 4-10 times lower than that with the crucible and encapsulant. Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 51 Table 3.16: Heat transferred (in Watts) from surface i to surface j for case-I with, one surface enclosure. SURFACE j SURFACE i Melt Crucible Encap Enclosure Melt 25 317 — 0 Crucible 361 2870 1305 8 Encapsulant — 1323 331 6 Enclosure 145 1675 803 10 Emitted (2) 531 6185 2439 24 Absorbed (1) 342 4544 1660 2633 Net (2)-(l) 189 1641 779 -2609 1400 . 1300 1200 1100 1000 * — & M e l t B---O E n c a p s u l a n t ° - - - o C r u c i b l e Wal 1 E n c l o s u r e 10 E L E M E N T NUMBER Figure 3.12: Effective ambient temperature versus element number for case-I / one surface enclosure. Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 52 • The crystal top loses heat to the enclosure (60%) and to the crucible wall, it absorbs heat from the crucible wall (47%) and the encapsulant (20%). The exchange with the other surfaces is 1-2 orders of magnitude lower. • The enclosure absorbs heat mainly from the crucible, encapsulant and crystal wall. It also absorbs heat from the melt and itself, however that exchange is 5-7 times lower. The effective ambient temperature versus element number is given in figure 3.13. The plot has the same general shape as that for the case of the full chamber. The temper-atures are about (2-7) % higher than those of the full chamber case. The differences are expected since the geometry of the two systems is different. However, since the per-centage difference range is low it can be accepted as a tradeoff to lower computing time. (The CPU time for this case is reduced by 63% from 628sec to 233sec). The temperature range here is (950-1450) K. 3.6.8 Analysis of Case-Ill Table 3.18 gives the numerical results of the heat transfer between the surfaces for case— III. The net heat exchanged in the system is 1.6kW. The crucible, melt, encapsulant and crystal lose heat to the enclosure which is at 600K. The following is a brief description of the exchange of each surface seperately: • The melt loses heat to the crucible (46%) and crystal wall (35%). It absorbs heat from the crucible wall (60%) and crystal (31%). • The encapsulant loses heat to the crucible (39%), crystal wall (33%) and enclosure, while it absorbs heat mainly from the crucible (59%) and crystal (28%) walls . It also absorbs, through multiple reflections, heat generated from its surface. Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 53 Table 3.17: Heat transferred (in Watts) from surface i to surface j for case-II with one surface enclosure. SURFACE j SURFACE i Melt Crucible Encap Enclosure Crystal Wall Crystal Top Melt 21 169 — 0 58 5 Crucible 174 1159 604 23 731 89 Encapsulant — 1072 229 22 271 54 Enclosure 93 978 495 56 277 327 Crystal Wall 84 973 263 8 161 28 Crystal Top 11 142 61 13 35 40 Emitted (2) 383 4493 1652 122 1533 543 Absorbed (1) 253 2780 1648 2226 1517 302 Net (2)-(l) 130 1713 4 -2104 16 241 1500 1400 . 1300 1200 . 1100 1000 900 * — * M e l t Q " " Q E n c a p s u l a n t 0 o C r u c i b l e W a l 1 E n c l o s u r e ° — "° C r y s t a l W a l 1 C r y s t a l T o p 10 E L E M E N T NUMBER Figure 3.13: Effective ambient temperature versus element number for case-II / one surface enclosure. Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 54 • The crucible exchanges heat mainly with the encapsulant and crystal and itself. It also loses heat to the melt and enclosure that is 5-10 times lower than to the other surfaces. • The enclosure absorbs heat from the crucible, encapsulant and crystal wall — all of the same order of magnitude. • The crystal wall exchanges heat with the crucible (37% of total heat lost and 49% of total heat gained) and encapsulant (11% of lost and 15% of gained). It also loses heat to the enclosure (22%) and gains heat originating from its surface (28% - through multiple reflections). • The crystal top loses heat mainly to the enclosure (46%) and it absorbs heat from the crystal wall (49%) through multiple reflections. Figure 3.14 shows the effective ambient temperature versus the element number of each surface. The general shape of the curves is similar to the full chamber case. The temper-ature range is (900 - 1400) K. The effective ambient temperatures for this case are within one percent of those for case-Ill with five surfaces enclosure (except for the crystal top 10 nodes where the differences are higher, 8%). That is explained by the fact that the melt, encapsulant, crystal and crucible walls exchange a higher percentage of heat with each other than with the other surfaces. Again, the differences for this case are accepted as a trade off to lower computing time ( CPU time is reduced by 63% from 704sec to 263sec). 3.7 Results of the simplified Case-Ill with Cutoff Wavelength 2.5 pm The results obtained using a transparent/opaque cutoff wave length A = 2.5pm instead of A = 2.0pm are given in table 3.19. The effective ambient temperature is plotted Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 55 Table 3.18: Heat transferred from surface i to surface j for case-Ill with one surface enclosure in Watts. SURFACE j SURFACE i Melt Crucible Encap Enclosure Crystal Wall Crystal Top Melt 18 122 — 0 64 0 Crucible 178 1093 647 15 1463 15 Encapsulant — 910 176 11 433 10 Enclosure 49 437 252 27 848 75 Crytal Wall 133 1754 546 58 1010 55 Crystal Top 6 54 31 10 106 9 Emitted (2) 384 4370 1652 121 3924 164 Absorbed (1) 204 3411 1540 1688 3556 216 Net (2)-(l) 180 959 112 -1567 368 -52 1500 1400 1300 1200 1100 1000 900 - -Q-<=>—a M e l t • E n c a p s u l a n t o - - - ° C r u c i b l e Wal 1 E n c l o s u r e ° — - ° C r y s t a 1 Wal 1 ^ - - ^ C r y s t a l Top T T 10 12 14 16 18 ELEMENT NUMBER 20 22 24 26 28 Figure 3.14: Effective ambient temperature versus element number for case-Ill / one surface enclosure. Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 56 versus element number in figure 3.15. It is clear that the effect of increasing the cutoff wavelength from 2pm to 2.5pm is the increase of the heat exchange with the melt surface and a decrease of exchange with the encapsulant surface. However, the exchange between the other surfaces did not change significantly (5-10 % difference only). Comparing figures 3.14 and 3.15 that show the effective ambient temperature, indicates that all the curves have the same general shape in both cases. However, the melt and encapsulant curves are closer in figure 3.15 than they are in figure 3.14 due to the increase (decrease) in the melt(encapsulant) heat transfer rate. The effective ambient temperature effects the fluid flow in the melt which also effects the liquid/solid interface. It is therefore important to obtain the correct cutoff wave length for better modelling of the fluid flow and heat transfer in the crystal puller. 3.8 Comparison of the Results of the Full and Simplified Chambers A comparison between the values of the heat transfer and equivilent ambient temper-atures that are obtained for both chambers is done in order to decide the importance of including the full chamber. The reported results show that the melt, encapsulant, crucible and crystal walls exchange heat between each other by higher percentages than they do with the other surfaces. The other surfaces being the middle cylinder wall, top and.lower annuli, and the top cylinder wall and top for the full chamber case. While for the simplified chamber, the other surfaces are simply the enclosure. The comparison between the effective ambient temperatures Tk,a of both chambers shows that those of the simplified chamber are generally 5% higher than those of the full chamber. The slight difference due to chamber geometry difference is expected and accepted. Also, the curves have the same trend as those of the full chamber case. Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 57 Table 3.19: Heat transferred from surface i to surface j for case-Ill with one surface enclosure in Watts for A = 2.5/zm . SURFACE j SURFACE i Melt Crucible Encap Enclosure Crystal Wall Crystal Top Melt 28 189 — 0 102 0 Crucible 274 1251 527 16 1551 10 Encapsulant — 747 138 10 345 6 Enclosure 76 440 190 27 708 49 Crytal Wall 209 1883 439 58 952 36 Crystal Top 10 55 24 10 88 6 Emitted (2) 597 4565 1318 121 3746 107 Absorbed (1) 319 3629 1248 1490 3577 193 Net (2)-(l) 278 936 72 -1369 169 -86 1400 T E M P E R A T U R E 1350 . 1300 1250 1200 1150 . 1100 1050 . 1000 * — * M e l t B---B E n c a p s u l a n t 0 o C r u c i b l e Wal1 *- - -* Enc 1 o s u r e C r y s t a l Wal 1 C r y s t a l Top ELEMENT NUMBER Figure 3.15: Effective ambient temperature versus element number for case-Ill / one surface enclosure and A = 2.5^m. Chapter 3. Radiative Heat Exchange in LEC GaAs Crystal Pullers 58 3.9 Conclusions Eventhough the results in this chapter are obtained using temperature distributions from previous numerical models, they still show the importance of a detailed radiative heat transfer model in GaAs crystal pullers. The effective ambient temperatures obtained for the surface elements show that the use of one ambient temperature for all surfaces introduces errors in the calculation of radiative heat transfer. The comparison between the results of the full and simplified chambers show that it is acceptable to use the simplified chamber as representative of the real crystal puller chamber when calculating heat transfer and fluid flow in the system. The comparison between the results of two cutoff wavelengths A = 2 and 2.5 pm shows that it is important to determine which one is a better representative of the encapsulant absorptive property. Chapter 4 Mathematical Modelling of Fluid and Heat Flow 4.1 Introduction The fluid flow and heat transfer during the crystal growth process are governed by the conservation of mass, momentum and energy equations. Those equations are solved numerically using a control volume method that is described in this chapter. The modes of heat transfer that are present during the crystal growth are (a)forced convection due to crystal and crucible rotation, (b) natural convection due to heating at the crucible wall, (c) conduction, and (d) radiation . This chapter also compares the results obtained using the programme developed by Sabhapathy and Salcudean in [25, 26, 27] when a simplified radiative model is used (radiation to a uniform temperature ambient) versus the detailed radiative model (radiation to effective ambient temperatures) developed in the previous chapter. Only the effect of natural convection is considered here since the main objective of the exercise is the determination of the effect of the detailed radiative model. 4.2 Assumptions The numerical results obtained are dependent on the following assumptions: • Axisymmetrical and laminar flow in the melt and encapsulant. • Cylindrical crucible and crystal whose axes coincide. 59 Chapter 4. Mathematical Modelling of Fluid and Heat Flow 60 • No slip and non-penetration conditions at the crystal and crucible surfaces and also at the melt-encapsulant interface. • Free slip at the encapsulant-gas interface. • Zero radial velocity and zero radial gradient of the temperature at the axis of symmetry. • Negligible convective heat loss from the encapsulant to the ambient gas and from the crystal to the gas. • Planar melt surface and planar encapsulant surface. • Boussinesq's approximation is used. 4.3 Governing Equations The partial differential equations that govern the fluid flow and heat transfer in GaAs, B2O3 and Ar during the crystal growth of GaAs are the mass, momentum and energy conservation equations. The equations can be written in general form for an axisymmetric flow in cylindrical coordinates as follows : where, d> is the dependent variable, is the exchange coefficient for the variable <f>, S(p is the source term. Table 4.3 defines the variables <f>, T^> and 5^ for each equation. The equation that governs the heat transfer in the crystal is : dT. d2T' 1 d 8T Ox ox2 r Or or Chapter 4. Mathematical Modelling of Fluid and Heat Flow 61 4.4 Boundary Conditions The following boundary conditions are imposed on the system (see figure 4.16 for further clarification of the terms) : • Bottom of the crucible: x = 0, r c > r > 0, u = 0, v = 0, w = 0 , § £ = 0 • Crucible wall upto melt surface: hm > x > 0, r - rc , % = 0,v = 0,w = 0, ^ = 0 • Crucible wall between melt surface and encapsulant surface: hm + he > x > hm, r = rc, u = 0,v = 0,w = 0, Tm = Tc • Melt surface between the crystal and crucible : x = hm, rc>r > r„ u = 0, -kms-t = e m < C - ra4mb) • Melt surface at the liquid/solid interface : x = hm, r„ > r > 0, u = 0, v = 0, w = 0, T = Tmeiting • Axis of symmetry upto the melt top : hm>x>0, r = 0, ^ = 0,^ = 0 ,^ = 0, ^ = 0 Chapter 4. Mathematical Modelling o f Fluid and Heat Flow 62 Table 4.20: Definition of the variables in the conservation equations. Equation 4> s* m a s s 1 0 0 x - m o m u r - m o m V 8 — mom w V-Energy T k 0 hs he hm ,CL Crystal B203 GaAs • Melt he rc Figure 4.16: Schematic of the calculation domain showing all variables of the system. Chapter 4. Mathematical Modelling of Fluid and Heat Flow 63 • Axis of symmetry between melt top and crystal top : h m + h, > x > h m , r = 0, Br U • Encapsulant-crystal interface : h m + h e > x > h m , r = ra, u = 0 ) V = Q)W = o, -k& + K*£ = *.<T? - r<m6) • Encapsulant—ambient interface : x = hm + he, re>r> r„ u = 0, £ = 0, £ = 0, -kee-£ = ^(Tt - TQU) • Encapsulant-crucible interface : h m + he > x > h m , r = rc, u = 0, v = 0, w = 0, Te = Tc • Crystal surface that is exposed to the ambient : h m + h, > x > he + hm, r — r„, u = up,v = 0,w = 0, -k,!£- = etcT(T?-TZmb) • Crystal top: x = h, + hm, r, > r > 0, u = up, v = 0, w = 0, -k.^f = e.a(T^-T^mb) • Release of latent heat at the crystal-melt interface : 2ro-{-km^ + k.B£)rdr = upr]LPt 4.5 Solution Procedure The partial differential equations introduced earlier are complex and cannot be solved analytically. The numerical solution procedure used in this work is based on the control Chapter 4. Mathematical Modelling of Fluid and Heat Flow 64 volume method. The control volume method requires that the solution domain be divided into a finite number of six-sided control volumes or cells. A grid-point is placed at the geometric centre of each control volume. This arrangement has the following advantages: • The value of the general variable that is directly available at the centre of the cell is a good representation of the average value over the control volume and can be used without interpolation to calculate the source terms and physical properties. • Discontinuities at the boundaries can be conveniently handled by locating boundary cells where discontinuities occur. A staggered grid arrangement is used since it solves some of the problems introduced by non-staggered grids (like producing non-realistic solutions) [45]. The scalar quantities (pressure and temperature) and the velocity component w are calculated at the geometric centre of the control volume, while the velocity components (u and v) are calculated at the scalar cell surfaces, (see firgure 4.17). In this way, the velocities are directly available for the evaluation of convection through the boundaries of the control volume. Figure 4.18 shows clearly the control volumes used in the analysis. 4.5.1 Discretization of the Differential Equations The finite volume equations are derived by the approximate integration of the general transport equation 4.8. The integral form of the differential equation for steady state conditions over a typical control volume for a general dependent variable <j> can be written as: (4.10) where A and V are the control volume area and volume Chapter 4. Mathematical Modelling of Fluid and Heat Flow 65 n is the outward normal unit vector, is the diffusion coefficient of the variable <j>. The left hand side of equation 4.10 expresses the convective and diffusive inflow/outflow through the area A of the scalar cell. It can be rewritten for a cylindrical coordinate system as follows (see figure 4.19 for definition of the coordinate system and idealized geometry): l(p'u4>-T4)Vd>)ndA= fXn(p'v<f>-TArcdx JA Jx, Or - / " ( A * - TM)..ir (4.11) Jrw OX The first two terms on the right hand side describe the total transport in the radial direction. They can be discretized as follows : {pv<t>e-T,t,—— )-(pv<}>w-T4> — ) (4.12) where 6r is the grid spacing in the radial direction (see figure 4.17) , and e &; w represent the east and west faces of the control volume. The variables that are not available at the scalar cells surfaces must be determind by linear interpolation between neighbouring cells : Pe = feP'p + (1 " f.)p'E (4-13) where fe is a special weighting factor defined by : / . = 0 . 5 | i (4.14) The total transport through the the east cell boundary consists of convective and diffusive fluxes; Ce = (p'v<f>)e and De = (r^|^)e. Due to the staggard grid arrangement, fluxes Chapter 4. Mathematical Modelling of Fluid and Heat Flow 66 are calculated at the scalar cell faces. However, the value of the variable <j> is directly available only at the nodal points. There are a few possible approximation schemes to obtain the fluxes at the scalar cell boundaries. The central difference scheme is one that assumes a piecewise-linear profile of the variable d> between nodal points. The convective flux can thus be expressed by : ft = ( P ' v ) . ^ (4.15) This formulation is second order accurate and gives satisfactory results for flows chara-terized by low Peclet number, Pee — ((p'v)e/-rf ), that is flows dominated by diffusion. However, in cases when convection is much larger than diffusion, the central difference method leads to numerical instabilities and yields unrealistic results. These difficulties can be overcome by using the upwind scheme approximation for convection dominated flows. The value of <f> at the scalar cell face is assumed to be the same as the value of <f> at the upwind grid point. The convective flux at the east scalar cell boundary is evaluated as: Ce = (p'v)e<t>P if ve > 0 (4.16) Ce = (p'v)ad>E if ve < 0 (4.17) The hybrid differencing scheme [45] combines the advantages of central differencing for small Peclet numbers \Pe\ < 2 and upwind differencing for large Peclet numbers \Pe\ > 2. As a result, the numerical stability and accuracy of the solution are enhanced. The discretized term representing convective and diffusive transport across the east boundry Ze can be approximated by : Ze = me(<t> - | ^ ) e = m e(l - 7 e & + 7e<M (4-18) Fe or where me = p'u and the coefficient 7 e is defined for hybrid central/upwind differencing Chapter 4. Mathematical Modelling of Fluid and Heat Flow 67 scheme as : 7e = ((1 - /) + max[-(l - f)Pe, fPe, 1]/Pe)e (4.19) Assuming a linear dependence between the source term and the dependent variable <p, the right hand side of equation 4.10 can be rewritten as follows: / S+dV = /"" r S+rdrdx = S* - SUP (4-20) Jv Jx, JTW where, Sf is the constant part of the linearized source term, Sp is the coefficient of <pp. The final form of of the discretized transport equation can thus be written as follows: me[(l - Ke)<f>E + 1e<t>p] + mw[(l - fw)4>W + lw<t>p\ + mn[(l - l n ) < t > N + LN4>P] + m.[(l - i.)<ps + l.M = Sf - SP<I>P (4.21) For numerical convenience a compact form of the above equation will be introduced : Apcpp = AE<i>E + AW(j>w + AN<f>N + AS(f>s + S* (4.22) which can also be written as : AP4>p = YJAz<pz-rSi (4.23) z where the subscript z denotes the neighbouring grid points, that is E, W, N and S, and the coefficients Az are defined as follows : AE = max[[Ce/2], De] - CJ2 (4.24) Aw = max[[Cw/2], Dw] + Cw/2 (4.25) AN = max[[Cn/2],Dn] - Cn/2 (4.26) As = max[[C./2],D.] + C,/2 (4.27) Ap = J2Az + Sp (4.28) Chapter 4. Mathematical Modelling of Fluid and Heat Flow 68 4.5.2 Discretization of the Boundary Conditions Imposed Temperatures and Velocities The values of the variables that are known or assumed can be directly assigned at the boundary points. However, in the case of an internal grid point, the given value of any variable <f> can be imposed via the source terms Sc and Sp in the following way : Sc = <f>Biven * lO 3 0 (4.29) SP = - l O 3 0 (4.30) Introducing these terms into the discretized governing equation 4.10 makes all other terms negligible and the equation simplifies to: <f>pSP « -Sc (4.31) and this forces the following value to be returned by the numerical procedure: <PP = = (frgiven (4-32) Free Surfacs The free slip condition (zero velocity gradiant) is imposed by cancelling the respective flux terms; e.g. for the east boundary : A* = 0 (4.33) Axis of Symmetry The zero radial gradient of the temperature at the axis of symmetry is imposed in the same manner as the zero velocity gradient at the free surface. Chapter 4. Mathematical Modelling of Fluid and Heat Flow 69 4.5.3 The Solution Algorithm The SIMPLE (Semi Implicit for Pressure Linked Equations) algorithm of Patankar and Spalding [45] is used to solve the system of nonlinear coupled equations. This algorithm uses an iterative solution which involves the following steps : 1. The fields of dependent variables u, v, w, p, T, are guessed for the first iteration loop. 2. The momentum equations in the r—, 6—, x— directions are solved and the first estimates of the velocity field u*, v*, w*, are obtained from the guessed pressure field. 3. A pressure correction equation is solved. 4. Corrections are made to the velocity components and pressure: u = u* + u' (4.34) v = v* + v' (4.35) w = w" + w' (4.36) p=p'+p' (4.37) 5. The equation for the scalar quatity T is solved. 6. The results are checked for convergence. If convergence is not attained then the new values of the dependent variables are used as the starting values, and the procedure from step 2 on is repeated. Chapter 4. Mathematical Modelling of Fluid and Heat Flow 70 4.5.4 Solution of the Linear Algebraic Equations The general finite volume equation 4.10 can be rewritten using the indices as follows +Aiij+1<Piij+1 + Atj.^j-i + S* (4.38) Figures 4.20 to 4.22 show control volume faces for the variables u, v and also the scalar variables defined at the grid node 'P' (1,3). The figures also show other variables as defined in the code. It is possible in principle to solve the above system of equations by matrix inversion. However, this is a computationally costly operation. An alternative solution procedure is to solve the system using an iterative fine by fine method. For the fine of constant i, equation 4.38 can be expressed as : <t>j = O j ^ i + i + bjd>j-i + Cj (4.39) where , *i = Aij+i/Aij (4-40) ^ = Aij.JAij (4.41) cj = (Ai+uh+u + Ai-uitu-ij + StyAij (4.42) The nonzero coefficients of the set of equations defined by equation 4.38 form a tri-diagonal matrix. Such a system of equations can be solved by a Gaussian-elimination method known as tri-diagonal matrix algorithm (TDMA) or Thomas-algorithm. When the values of d> on the fine are found , the same procedure is carried out for all fines in the x-direction and can be repeated for the r- direction as well (see figure 4.23). The convergence of the above method is fast because the boundary conditions information is transmitted at once to the nodal points lying inside the solution domain. Chapter 4. Mathematical Modelling of Fluid and Heat Flow 71 4.5.5 N u m e r i c a l C o n s i d e r a t i o n s U n d e r r e l a x a t i o n In order to handle nonhnearities in the iterative solution of the algebraic equations, it is necessary to underrelax or slow down the iteration process. The dependent variable, i.e. the temperature T, can be underrelaxed by the introduction of an underrelaxation factor a into the general discretized equation : ApTp = JL AZTZ + S^. Thus the discretized equation can be written as follows : ^TP = £ AZTZ + STC + (1 - a)^TP (4.43) a t a where Tp is the value of Tp from the previous iteration. At convergence, Tp = Tp and the equation with the underrelaxation factor is identical to the discretized equation. In the present work, the following underrelaxation factors are used : a u = 0.25 , av = 0.25 , ctu, = 0.25 , a p = 0.50 , CLT = 0.25. The source term can be underrelaxed as follows : Sc = ctSc + (1 - ct)S*c (4.44) where S* is the source term from the previous iteration. Another way of controlling convergence is by the use of the transient approach. The practice of solving a steady state problem via an unsteady formulation can be regarded as a special kind of underrelaxation. C o n v e r g e n c e C r i t e r i a An iteration procedure is normally converged when further iteration does not produce any significant changes in the values of the dependent variables. In addition, it is necessary to have a convergence criterion which measures the degree to which a computed solution Chapter 4. Mathematical Modelling of Fluid and Heat Flow 72 satisfies the original set of transport equations. In the present work, this is done by using the 'residual source' method. The residual Rj, over one cell is denned as : Rt = 52Aa<pM + S*-AP<pp (4.45) X and the sum of the residuals for the whole field is : field At convergence R% becomes 0. For practical purposes, it is sufficient for the sum of the residuals to drop below a specified small number , i.e. : (4.47) where R™ is the reference value for the variable <f>, and A is the level of convergence for all dependent variables (e.g. A = 0.0001 ). The value of R% can be monitored during the iteration process and any unexpected behaviour can indicate the source of error in the programme. For example, a divergent behaviour could mean errors in the implementation of some of the boundary conditions. False Source In order to ensure convergence and stabilize the numerical behaviour of the solution, one may introduce into the finite volume equations the so called 'false source' term Sf defined by: Sf = Cp(<f>x 1 — 4>x) (4.48) where Cp is the mass imbalance for a control volume defined by: (4.49) z Chapter 4. Mathematical Modelling of Fluid and Heat Flow 73 The false source term can be incorporated into the discretized equation as follows : (AP + CP)<f>P = £ Az<t>, + S* + Cpfc1 (4.50) z At convergence, Cp and Sj become 0 and have no effect on the final converged solution. The addition of a false source term allows avoiding the singularity problem which may arise during the iteration process when Ap takes on the value 0 ( corresponding to zero d> transport into the cell). False Diffusion For high Peclet numbers (\Pe\ > 2), the hybrid scheme reduces to the upwind scheme. The upwind scheme, though numerically stable, has the disadvantage of being only first order accurate and introduces a discretization error referred to as false diffusion. The effect of false diffusion is an artificial increase in the (physical) diffusion coefficients which results in smearing of the gradients in the flow field. The errors introduced in the solution can be quite important, particularly when the grid fines are inclined with respect to the local velocity vectors [46]. False diffusion can be limited by using (a) a fine mesh, (b) grid fines that are aligned with the flow direction, (c) discretization schemes accounting for the multidimensional nature of the flow and involving more neighbouring cells, and (d) a curvilinear coordinate system. Numerical Errors and Instabilities The most likely sources of numerical errors are : (a) unrealistic initial values and (b) grids that are too coarse. Numerical instabilities are often due to (a) inappropriate underrelaxation factors (usually too large), (b) too few sweeps in the iterative line-by-line method. The correct evaluation of numerical parameters can be based on previous experience or can be achieved through trial and error. Chapter 4. Mathematical Modelling of Fluid and Heat Flow 74 The errors in the results are usually caused by mistakes in the implementation of boundary conditions or wrong input data (e.g. physical properties or geometry). 4.6 Simplified Radiative Model Results The system is assumed to be radiating to a uniform ambient at 1211.K'. The stream functions and isotherms in the melt, encapsulant and crystal are shown in figures 4.24 and 4.25. Figure 4.24 shows the stream lines in the melt and encapsulant. The melt flows upwards near the crucible wall due to heating at that surface. It is then cooled by the encapsulant and the crystal, falls towards the crucible base and is pushed towards the crucible wall. The flow near the crucible wall, base and near the encapsulant and crystal is one of the boundary layer type where high velocity and temperature gradients occur as can be seen clearly in figures 4.24 and 4.25. The flow velocity at the centre of the melt is lower than at the edges. The flow in the encapsulant moves upwards near the crucible wall and is pushed towards the crystal where it is cooled and sinks downwards near the crystal. The flow velocity in the encapsulant is two orders of magnitude lower than that in the melt. Thus, conduction dominates convection in the encapsulant. The crystal isotherms show a high gradient near the encapsulant region, however, the gradient decreases with increasing distance from the encapsulant top. The isotherms are curved near the crystal surface where radiative effects are greatest but they flatten out towards the centre. The encapsulant isotherms are nearly flat except near the crucible wall where it is heated. At the encapsulant free surface, there is a high radial gradient of the temperature near the crucible wall due to heat gain from the wall and heat loss by radiation to the Chapter 4. Mathematical Modelling of Fluid and Heat Flow 75 ambient. The melt isotherms are highly concentrated near the melt-crystal-encapsulant junction and near the crucible wall where the highest rate of heat transfer takes place. 4.7 Detailed Radiative Model Results The results obtained when a detailed radiative model is used are shown in figures 4.26 and 4.27. The melt flows in a manner similar to that described earlier: upwards near the crucible wall and downwards under the crystal. The central region has a lower velocity since it is away from the direct effects of heating and cooling. The encapsulant streamlines show a circular flow pattern. The crystal isotherms show a high radial gradient near the upper half of the encap-sulant region; they also have a steep slope in that region near the crystal surface. The temperature of the crystal surface decreases near the encapsulant-crystal interface away from the melt surface, however it starts increasing again above the encapsulant top due to exposure to a higher ambient temperature (reflected from higher heat absorption from the crucible wall and the encapsulant). The isotherms in the top half of the crystal are nearly flat. The isotherms in the melt show high concentration near the melt-crystal-encapsulant junction and near the crucible wall where convective effects are greatest. The encapsulant isotherms show the domination of conduction in that layer; they also show a uniform gradient at the encapsulant-ambient interface. 4.8 Comparison between the Radiative Models It is clear from the previous two sections that the temperature fields in the melt, encap-sulant and crystal are affected by the introduction of a radiative model that takes into Chapter 4. Mathematical Modelling of Fluid and Heat Flow 76 account exchange between all surfaces. The field in the crystal is highly affected espe-cially at the crystal surface where heat gain or loss by radiation occurs. It is also clear that the crystal temperature distribution resulting from the detailed model shows higher values in the top half (1336-1211)K than the simplified model (1261-1211)K. That is due to reduced heat loss to the surroundings shown by radiation to a higher ambient temperature. The melt surface loses more heat by radiation when the detailed model is used since it is radiating to a lower ambient temperature distribution than the constant ambient tem-perature. The effect of that on the temperature distribution is clear where the isotherms are not as concentrated near the crystal. (Melt cooling drives the flow faster downwards). The temperature distribution of the encapsulant free surface changes significantly from one model to the other. In the simplified model, the gradient is high near the crucible wall (w 150K/cm) while it decreases drastically to (w 20K/cm) away from the crucible wall. When the detailed model is used, the gradient obtained is more uniform across the surface (w 45K/cm). That is due to including detailed calculations of heat exchange between the encapsulant and the crystal, crucible and ambient. Eventhough, the analysis presented here included only the effects of natural convec-tion, the exercise does evaluate the effect of including a detailed radiative model on the temperature distribution in the crystal, melt and encapsulant. Therefore, it can be con-cluded that for better representation of heat transfer during the crystal growth process, a detailed radiative model be used. Chapter 4. Mathematical Modelling of Fluid and Heat Flow 77 LOCATION VARIABLES STORED t Y///A CONTROL VOLUME FOR J * [SS^J CONTROL VOLUME FOR « | 1 CONTROL VOLUME FOR u Figure 4.17: A finite volume grid in the (x-r)-plane. Chapter 4. Mathematical Modelling of Fluid and Heat Flow Figure 4.18: Control volume and notation for (a, (c) transport of scalar quantity in ^ -direction. ,b) u- and v-momentum analysis, Chapter 4. Mathematical Modelling of Fluid and Heat Flow P, w,T UL Cryatal Fnrapaulant do Alt melt Figure 4.19: The coordinate system and idealized geometry. Figure 4.20: Control volume for all scalar varibles denned at the grid-node (I,J) Chapter 4. Mathematical Modelling of Fluid and Heat Flow 80 © XU( I ) © nvu) YV(J) K - D X E P ( I ) «H • » S E W ( I ) * T » — SEW( I * 1 ) «H Figure 4.21: Control volume for u-velocity denned at the grid node (I,J) 0 © X ( l | *»- SEW(I) -H Figure 4.22: Control volume for v-velocity defined at the grid node (I,J) Chapter 4. Mathematical Modelling of Fluid and Heat Flow 81 \ CHOSEN LINE Figure 4.23: Representation of the line-by-line method. Chapter 4. Mathematical Modelling of Fluid and Heat Flow 82 0.075 0.000 ir=i> X 1.0»107, contour spoclng is 5.0 in the melt X 1.0*10*, contour spacing is 5.0 in the encapsulant Figure 4.24: Stream functions of the melt and encapsulant - simplified radiative model. Ambient temperature = 1211K, crucible wall temperature = 1531K. Chapter 4. Mathematical Modelling of Fluid and Heat Flow 83 0 0 7 5 0.000 contour spacing is 2.0K in the melt contour spacing is 25.0K in encapsulant and crystal Figure 4.25: Melt, crystal and encapsulant isotherms - simplified radiative model. Am-bient temperature = 1211K, crucible wall temperature = 1531K. Chapter 4. Mathematical Modelling of Fluid and Heat Flow 84 i n b i n o o o 0.075 0.000 i>=i> X 1.0*107, contour spacing is 5.0 in the melt if=i> X 1.0*10*, contour spacing is 5.0 in the encapsulant Figure 4.26: Melt and encapsulant stream functions - detailed radiative model. Crucible wall temperature = 1529K. Chapter 4. Mathematical Modelling of Fluid and Heat Flow 85 in ci 6 in o o ^ ^ _. ^ q 0 0 7 5 0.000 contour spacing is 2.0K in the melt contour spacing is 25.0K in encapsulant and crystal Figure 4.27: Melt, encapsulant and crystal isotherms - detailed radiative model. Crucible wall temperature = 1529K. Chapter 5 Conclusions and Recommendations 5.1 Summary This thesis presents a description of the LEC crystal growth process of GaAs. It defines the modes of heat transfer present in the system but concentrates on evaluating the radiative mode through the use of a simplified and a detailed model. It also presents an analysis of the effect each model has on the fluid flow and heat transfer in the whole system. The conclusions drawn from the research and numerical work done are presented in the next section. In addition, some recommendations are suggested for better evaluation and understanding of the heat transfer that takes place during the crystal growing in LEC chambers. 5.2 Conclusions The following points may be concluded from the work given in the previous chapters: • Radiation is an important mode of heat transfer that must be included for proper modelling of the crystal growth process. • The upper part of the growth chamber may be replaced by one isothermal surface for the numerical analysis purposes without significant loss of accurate representation of the real system. 86 Chapter 5. Conclusions and Recommendations 87 • The use of a simplified radiative model where all surfaces are radiating to one am-bient temperature does not evaluate the radiative heat transfer during the growth process correctly. • A detailed radiative model changes the melt, encapsulant and crystal isotherms. The changes in the crystal isotherms are quite significant especially near the en-capsulant top. • The emissivities of GaAs (solid and liquid), the encapsulant (B2O3), the crucible material (PBN) and stainless steel which forms the upper chamber surfaces are important factors in the determination of radiative heat transfer. • Varying the cutoff wavelength of the encapsulant absorptive property from 2/xm to 2.5[im effects the effective ambient temperature distribution of the melt, encapsu-lant and crystal. 5.3 Recommendations The following recomendations are suggested for future numerical analysis work in the field of GaAs crystal pulling : • Evaluation of the emissivity of solid and liquid GaAs so that the sources available now may be verified. • Evaluation of the emissivity of the encapsulant B2O3 which is still not available in the literature. • Evaluation of the absorptive property of the encapsulant B2O3 in order to verify the results available from [21]. Chapter 5. Conclusions and Recommendations 88 • Making temperature measurements in the melt and encapsulant layers as possible, and near the crystal surface and crucible wall surface so that better boundary conditions can be imposed on the system when performing numerical simulations of the growth process. A p p e n d i x A C o n f i g u r a t i o n F a c t o r s of the S y s t e m A . l I n t r o d u c t i o n This appendix gives the equations used to calculate configuration factors for the cases analysed. Configuration factor will be referred to as CF. The relationships are obtained from [44]. A.2 C o n f i g u r a t i o n F a c t o r E q u a t i o n s The following is a list of the equations used in the programme written to calculate the radiative heat transfer in the GaAs crystal puller. 1. CF from the ring element on tube to ring element on coaxial annular element on circular fin. Figure A.28 shows rx, r2, dr and I. With the following definitions : L = l/ru R = r2/ri, Z = 1 + R2 + L2 a i = Z2-AR2, a2 = R2 - L2 - 1 the equation is ,2LR,,„2 , N i / 2 2*a2 lt(Z + 2R)(R- 1 ) , 1 / 2 W , W A " "J) + -gr , <"-'((2-MQ(fi + i)> This equation is used to calculate CF from a crystal element to a melt or encapsulant element, or from the crucible outer wall to a middle cylinder lower annulus element. 89 Appendix A. Configuration Factors of the System 90 2. CF from the outer surface of cylinder to annular disk at end of cylinder. Figure A.29 shows the cylinder and disk and the variables R, r and I. With the following defi-nitions : A = l2-R2 + r2, B = l2 + R2- r2, C = l2 + R2+r2 £ = ( £ " 4 ( ? f ) 1 / 2 , E = C 0 5 _ 1 ( s i ) the equation defining CF is : Fx.2 = ±- * {cos-\^) - {±){D * E + ( £ ) « n - l ( £ ) - ( i r /2 ) (£) ) ) (A.52) Surfaces 1 and 2 may be the crystal element #1 and the melt element #1, re-spectively. They may also be the crucible outer wall element #1 and the middle cylinder lower annulus element #1, respectively. 3. CF from an annular ring on a cylinder to an annular disk at the end of the cylinder. The variables r, R and x are as shown in figure A.30. Therefore, with the following definition of a, b, c, d, and e : o = x2 - R2 + r2 b = x2 + R2-r2 c = x2 + R2 + R2 J c a — (e»-4(ilr)a)1/» e = ra/Rb the CF equations is given by : Fdl-2 = ^-(cos-^a/b) - -(d * cos-\e) - cos^ir/R))) (A.53) Appendix A. Configuration Factors of the System 91 Figure A.30: Schematic showing ring on cylinder and disk connected at end of cylinder. Appendix A. Configuration Factors of the System 92 This equation may be used to calculate CF from a crystal element #(2-ns), to the melt element #1. It may also be used from a crucible outer wall element #(2-nco) to the middle cylinder lower annulus element #1. 4. CF from the inside surface of right circular cylinder to itself. Section and variables are shown in figure A.31. With the following definition : H = h/2r , the CF equation is : F1.1 = 1 + H - (1 + H2f2 (A.54) This equation will be useful in evaluating CF from the crucible wall to itself, from the middle cylinder wall to itself and from the top cylinder wall to itself. 5. CF from the base of a right circular cylinder to inside surface of cylinder. With the variable H = h/2r , where h and r are as shown in figure A.32, the equation is : = 2tf((1 + H2f'2 - H) (A.55) The CF from the crucible wall or the the middle cylinder wall to the fictitious surface (see section 3.3.2) may be calculated using the equation. 6. CF from a disk in cylinder base or top to inside of right circular cylinder. With the variables R = r 2/ri and H = h/r1} where Tx, r2 and h are as shown in figure A.33, the equation for the CF is : Fi_ 2 = (^1 -R2-H2 + ((1 + R2 + H2)2 - IR2)1'2) (A.56) This equation is used to calculate CF from the middle cylinder wall to the crystal top. 7. CF from an annular ring on cylinder base or top to inside of right circular cylinder. The section and variables are shown in figure A.34. With the same definition of H Appendix A. Configuration Factors of the System Figure A.31: Schematic showing surface-1. Figure A.32: Schematic showing surfaces 1 and 2 on a right circular cylinder. Figure A.33: Schematic showing surface-1, disk an base of cylinder, and surface-2 Appendix A. Configuration Factors of the System 94 and R as in the previous point, and the following definition of a, b and c : a = R2 - 1 b = (AR2 + H2)1'2 c = ((1 + R2 + H2)2 - 4R2)1'2 the equation for the CF is : Fx-2 = i ( l + kH*b-c)) l a (A.57) Surface 1 may be a melt, encapsulant or middle cylinder lower annulus element and surface 2 may be the crucible wall or the middle cylinder wall. 8. CF from the interior of right circular cylinder to finite annular ring in base.The variables rj, r 2 , 1 and a are as shown in figure A.35. Therefore, with the following definition of R, L and X : This equation is used to calculate CF from the crucible wall to a melt or encapsulant element, or from the middle cylinder wall to a lower annulus element. 9. CF from the interior of a finite length right circular coaxial cylinder to itself. With the same definition of R and H as in the previous point (section shown in fig-ure A.36), and the following definition of a, b, c, d and e : R = r/o, L = l/a, X = (L4 + 2L2(1 + R2) + (1 - R2)2f'2 The governing equation is : Fx_2 = ±r(X, - X, + (R2)2 - (i2a)2) (A.58) 2 ( f l 2 - l ) » / 3 H (4fl 2+H 3) 1/ 2 H Appendix A. Configuration Factors of the System 95 Figure A.34: Schematic showing surface-1, annular ring on base of cylinder, and surface-2. Figure A.35: Schematic showing a right cylinder and a finite annular ring on its base. Figure A.36: Schematic showing surface-1 for which CF is calculated. Appendix A. Configuration Factors of the System 96 _ 4(R3-D+{H/R)3(R7-2) H the equation governing the configuration factor is : F^ = 1 - i + * ron_1(o) - ^-(b * sin~l(c) - sin~l(d) + ^ * e)) (A.59) Surface 1 may be either the crucible wall or the middle cylinder wall. 10. CF from the interior of the outer right circular cylinder of finite length to exterior of inner right circular coaxial cylinder. The section is shown in figure A.37. With the same definition of R and H and the following definition of a, b, c, d and e : a = H2 - R2 + 1 b = H2 + R2 - 1 c = H2 + R2 + 1 d = a/b e = (c2 - (2J2)2)1/2 / = e * cos-l(d/R) g = a* 5in-1(l/i2) the CF governing equation is : This equation may be used to calculate CF from the crucible or middle cylinder walls to the crystal wall. 11. CF from the interior of an outer right circular cylinder of finite length to the annular end enclosing space between the coaxial cylinders. The section is as shown in figure A.38. With the following definitions of H, R, a through i and kx through (A.60) Appendix A. Configuration Factors of the System 97 h : H — h/r2, R = ri/r2 a = 1 - R2 - H2 b = 1 - R2 + H2 c = l + R2 + H2 d = 2R2 - 1 e = (c2 - AR2)1'2 f = (4 + H2fl2 g = (1 - fl2)1'2 i = R* a/b h = tan-\glH) - tan^^g/H) k2 = sin~1(d) — sin-1(R) ks = | + am" 1^) k4 = | + sin - 1 ( t ) fcs = f + «n"1(fc) the governing equation will be : FX-2 = -(R *kl + ^-*k2 + ^-*k3--^-*k4 + ^-*k6) (A.61) 7r 4 4ii 4ii 4 v ' This equation is used to calculate CF from the crucible or middle cylinder walls to the fictitious surface a,- or the fictitious surface b{. 12. CF from a differential element on annulus between coaxial cylinders to interior of outer cylinder as shown in figure A.39. The equation for this CF will have to be integrated over the whole annulus for it to be useful here. With the following definition of H, R, a, b, c, d and w : H = h/r2 Appendix A. Configuration Factors of the System 98 Figure A.37: Schematic showing surfaces 1 and 2 on two coaxial cylinders. Figure A.38: Schematic showing annular ring between two coaxial cylinders to which CF is calculated. Figure A.39: Schematic showing differential element on annulus between two coaxial cylinders from which CF is calculated. Appendix A. Configuration Factors of the System 99 R = r / r 2 a = 1 + R2 + H2 b = l-R (v = cos 1(ri/r) + cos 1(ri/ra) the equation for the CF can be written as follows : l — i —i & —1 —i w Fdi-2 = —tan (b*tan (—)) + c*tan (d*tan (—)) 7T (A.62) After integrating this equation, it can be used to calculate CF from a melt or encapsulant element to the crucible wall. It can also be used for the configuration factor from a lower annulus element to the middle cylinder wall. 13. CF from the wall of the smaller radius cylinder to the wall of the other coaxial cylinder as shown in figure A.40. The variables R, H\,H2, L\ and i 2 are denned as follows: This equation is used for calculating CF from the top cylinder wall to the middle cylinder wall, and from the crucible wall to the middle cylinder wall. R = r i / r 2 , H2 = h2/ri, L\ = R2 -f H2, X 2 = R? + H\. The configuration factor is denned by : (l-L2 + ((l + L2)2-AR2y/2) (A.63) Appendix A. Configuration Factors of the System 100 14. CF from the interior of a circular cylinder of radius ri to a disk of radius r 2 as shown in figure A.41. With the following definition of R,Hi,H2,Xi and X2 : R = ri/r 2, Hi = hi/r2, H2 = h2/r2, Xx = HI + R2 + 1, X2 = H2 + R2 + 1 the CF equation is X,-X2- (X2 - AR2)1'2 + (X2 - 4R2)1/2 Fl~2 AR(H2-Hl) ( A ' 6 4 ) This equation is used to calculate CF from the middle cylinder wall to the fictitious surface Oj at the crucible edge or at the top cylinder edge. 15. CF from the inside surface of a right cylinder to a coaxial disk of the same diameter seperated from the base of the cylinder as shown in figure A.42 . H\, H2, oi, a2, 0 3 and 0 4 are defined as follows : #i = h/r, H2 = h2/r, a ^ l + ffi/ffj, 02 = ^ + H\ a 3 = (4 + (^1 + ^ 2) 2) 1/ 2 , a 4 = (4+ # 2) 1/ 2 and CF is defined as follows : * a 3 - a 2 - a 4 F i - 2 = (A.65) This equation calculates CF from the middle cylinder wall to a fictitious surface aj or from the crucible wall to a fictitious surface a; or b{. 101 Figure A.40: Schematic showing two coaxial cylinders, one on top of the other for which CF is calculated. Figure A.41: Schematic showing surface-1, circular cylinder, and surface-2, disk of smaller radius. Figure A.42: Schematic showing surfaces 1 and 2 for which CF is calcualted. Appendix A. Configuration Factors of the System 102 16. CF between two parallel disks with centers along the same normal as is shown in figure A.43. The variables r i , r 2 and h are as defined in the figure. With the following definition of Ri, R2 and x : This equation may be used to calculate CF from a fictitious surface O j to another fictitious surface O j . A . 3 Configuration Factor Calculation of Case-I The following is a description of the configuration factor calculation for case-I of the simplified chamber. This case is chosen to illustrate the solution approach because it is the simplest case. All other cases are solved in a similar manner. 1. Equation A.54 is used to calculate CF from the crucible wall element to itself, 2. Equation A.55 is used to calculate CF from the fictitious surfaces 'a' and 'b' to the crucible elements. The reciprocity relationship is then used to calculate CF from the crucible wall elements to the fictitious surfaces. The CF from crucible element C j to element C j + i can be calculated by the following equation : R2 = r2/h, * = ! + (! + R\)IR\ the equation of the CF will be defined as : Fl.2 = ^(x-(xi-A(r2/r1)2)) (A.66) Fa-ci-Fci-ci+i = Fa (A.67) Appendix A. Configuration Factors of the System 103 3. Equations A.57 and A.58 are used to calculate CF from the melt or encapsulant elements to the crucible wall elements. The reciprocity relation is then used to calculate CF from the crucible wall element to the melt or encapsulant elements. 4. CF from from each element to all other elements are summed up. The sum is then subtracted from one to obtain the CF from the each element to the enclosure. Table A.21 gives the summation of CF from surface i to surface j that are obtained using the method described above. The total from each surface is equal to 1. All dimensions of the puller are as given in the main text (table 3.10). Appendix A. Configuration Factors of the System 104 Figure A.43: Schematic showing two disks of different radii for which CF is calcualted. Appendix A. Configuration Factors of the System 105 Table A.21t CF summation from surface i to surface j for case-I of the simplified chamber. ELEMENT # & SURFACE i SURFACEj ENCAPSULANT CRUCIBLE ENCLOSURE ENCAPSULANT 1 — 0.3348 0.6652 2 — 0.3596 0.6404 3 — 0.4123 0.5877 4 — 0.4964 0.5036 5 — 0.6073 0.3927 Crucible 4 0.4756 0.2733 0.2511 5 0.4294 0.2905 0.2802 6 0.3868 0.3009 0.3123 7 0.3478 0.3043 0.3478 8 0.3868 0.3009 0.3868 9 0.2802 0.2905 0.4294 10 0.2511 0.2733 0.4756 ENCLOSURE 1 0.1678 0.1655 0.6667 MELT CRUCIBLE ENCLOSURE MELT 1 — 0.5050 0.4950 2 — 0.5250 0.4750 3 — 0.5644 0.4356 4 — 0.6210 0.3790 5 — 0.6891 0.3109 CRUCIBLE 1 0.4756 0.3438 0.1805 2 0.4294 0.3691 0.2015 3 0.3868 0.3883 0.2249 4 0.3478 0.4011 0.2511 5 0.3123 0.4075 0.2802 6 0.2802 0.4075 0.3123 7 0.2511 0.4011 0.3478 8 0.2249 0.3883 0.3868 9 0.2015 0.3691 0.4294 10 0.1805 0.3438 0.4756 ENCLOSURE 1 0.1273 0.2060 0.6667 Appendix B Gebhart Factor Equations Derivation B.l Introduction The derivation of the Gebhart factor equations of the system surfaces is given in this appendix. A detailed derivation is shown for only one surface being the melt. The equations of the other surfaces will be written in a similar manner as described in each section. The genaral form of the Gebhart factor equation is given in the main text and is given again here for convenience [43]. Gjk = Fj-kCk + Fj-ipiGi-k + Fj-2p2G2-k + + Gj-kpkGk-k + + Fj-Np^Gs-k (B.68) B.2 Melt Gebhart factor equations The following set of equations are the Gebhart factor equations for the melt surface which is divided to nm elements. The variable i goes from 1 to nm in equation B.69 ; from 1 to nc for equation B.70 ; from 1 to nmc for equation B.71 ; from 1 to nla for equation B.72 ; from 1 to ns for equation B.77, from 1 to nst for equation B.76 and from 1 to nco for equation B.78 for each k. k goes from 1 to nm . These equations include the crystal surface elements, the crystal top and the outer surface of the crucible. However, for the case which does not include the crystal, the terms involving the crystal and outer crucible 106 Appendix B. Gebhart Factor Equations Derivation 107 will be eliminated. All other terms remain the same. Gm-m(l, /c) — fm—m{}} ty^m (1|*0 + fm-m(ii2)pmGm-m{2,k) + + fm-m(i,nm)pmGm-m(nm, k) +/m_ e(i, l > c G c _ m ( l , k) + fm_e(i, 2 > e G c _ m ( 2 , k) + + fm-c(i,nc)pcGc-m(nc, k) + fm-co{i, l)PcoGco-m(l, k) + /m-co(i, 2)pcoGco-m(2, k) + + fm-co{i,nco)pcoGco-m(nco, k) +/«-.(*, l)p,G,-m(l, k) + fm-.(h 2)p,G,-m(2, k) + + fm-.(i, ns)ptG,-m(ns, k) 4~/m-mc(^> 1 )pmcGmc_m(1, fc) +fm-me(i, 2)pmcG (2 ,*) + + /m-mc(i) ™ c ) / ' m c G r a c - r a ( i r o C , A) —mc/a —meJa (*) 2)/JmcZaG!mcIa-m(2, fc) + + fm-mclaiiynl^PmclaGmcla-minla, k) + fm (i)PmctaGmcta-m(k) +fm-tcw(i)ptcwGtcw-m(k) + fTn.-tet(i)ptctGtct-Tn(k) +fm-.t(i, l ) p . t G . t - m ( l , *) + / m - r f ( » , 2 ) p . « G r t _ r a ( 2 , fc) + + fm-.t(i, nst)p,tG.t-m{nst, k) (B.69) Appendix B. Gebhart Factor Equations Derivation 108 Gc-m{i,k) = fc-m(i, k)em + / c-m ( i , l)PmG (1, fe) + / e_ m(», 2 > m G m _ m ( 2 , k) + + fc-m{i,nm)pmGm-m(nm, k) +fc-e(i, l)pc<7c-m(l, *) + / c - c ( » , 2)p c G c _ r o ( 2 , fe) 0 + + fe-e{i, nc)pcGc-m(nc, k) +fe-eo{i, l)PcoG co—m (1, fe) + /»_«,(», 2)/> c oG c o _ m (2 , fe) + + fc-co{i, nco)pcoGco.m(nc, fe) + / « - . ( * , l ) p . G . _ m ( l , fe) + / c - . ( i , 2 ) p # G . _ T O ( 2 , fe) + + fc-*(i, ns)p.G,-m(nc, fe) (1 ,*) +fc-mc(i, typmcG (2,fe) + + / c-mc(i)7i7nc)p m cG m c_ m(nmc, fe) + / c —mc/o (*> IJPmc/aGmcJo-mCl) fe) (i> 2)p,„cIaC7mc/a-m(2, fe) + + fc-mcla(i,nla)pmclaGmcla-Tn(nla, fe) + fc —mcta{i') Pmcta Gmcta—m (fe) +fc-tcw(i)PtcwGtcw-m{k) +fc-tct(i)PtctGtct-m(k) + /e-.t(», l > , t G . t _ m ( l , fe) + / e _ . t ( * , 2>. tG . t _ m ( 2 , fe) + + fc-.t{i, nst)pttGet.m(nc, k) (B.70) Gmc-m(i> k) = fmc-m(i) k)em 4"/mc-m(i> typmGm-m (1, fe) + fmc-m(i, 2 > m G m _ m ( 2 , fe) Appendix B. Gebhart Factor Equations Derivation 109 + + fmc-m{i,nm)pmGm-m(nm, fc) +/me-e(», l)PcGc-m(l, k) + fmc-c{h 2)PcGc.m(2, k) + + fmc-c(i,nc)peGc-m(nc, k) + /mc-co(i , l)PcoGeo-m(l, k) + fmc-co(i, 2)pcoG (2,*) + + frnc-co(i,nco)pcoGco-m(nco, k) +fmc-.(i, l)p.G..m(l, k) + fme-.(i, 2)p.G.-m{2, k) + + fmc-,(i,ns)P,Gt-m(ns,k) fmc—mc(^) tyPmcGmc—m (1,*) (2,k) (nmc, k) ~\~fmc—mcla(i> tyPmdaGmcla—m(l> &) ~T fmc-mcla{} •> 2)pmclaGmcla—m(2, fc) + + frnc-rncla{hnla)PrnclaGmcla-m{nla, fc) (i)pmctaGmcta-m{k) +fmc-tcw(i)ptcwGtcw-m(k) -rfmc-tct{i)ptctGtct-m{k) + / m e - . t ( » , l ) p . t G . t - ™ ( l , *) + / m c - . « ( » , 2)PttGtt-m{2, k) + + fmc-,t(i, nst)petG,t-m(nst, fc) (B.71) C m c l o - m ( ' i^ ) = fmcla—m{}jk)€Tn + fmcla-m(i, ^)pmG (1, k) + fmcla-m(i, 2)pmG (2,fc) + + frncia-rn(i,nrn)pmGrn-rn(nm, fc) + / m d Q - c ( i , l > c G c - m ( l , *0 + / m c i o - c ( i , 2 ) p c G c _ m ( 2 , fc) Appendix B. Gebhart Factor Equations Derivation 110 + + fmcia-c(i,nc)pcGc-m(nc, fc) co—m (1, k) + fmda-cb(i, 2)pcx>G (2,*) + + fmcia-co(i, nco)PcoGeo-m(nco) k) + + fmcia-,(i,ns)p,Gt-m(ns, fc) ~T fmcla—nc('i l)PmcC m c _ m (1,*) ~T fmcla—mc('i 2)pmcGmc—m (2,k) + + / m c J a - m c ( i , « " l c ) ^ m c G m c _ r n ( n m C , fc) ~\~ fmcla—mcla{i"> 1 ^ Pmcla Gmcla—m (1,*) —mcio (2,fc) + + fmcla-mcla(i, rila)pmciaGmcla-m(nla, fc) —mcta (i)Pmcta Gmda-m {k) ~rfmcla-tcw(i)PtcwGtcw-m(k) -rfmcla-tct(i)ptctGtct-m{k) +fmcla-tt{i, 1)p»tG,t-m(\, fc) + fmcla-»t{h 2)p,tGBt-m(2, A) + fmcia-,t(i,nst)pttGtt-m(nst,k) (B.72) mcta-m^) — Jmcta—m\™)€m •rfmcta-m{X)pmG (1, fc) + fmcta-m{2)pmGm-m(2, fc) + + fmcta-m(nm)pmGrn-m(nm, fc) + /mda-c(l)PcGc-m(l,*) + fmcta-c{2)pcGc.m{2,k) + + fmcta-c{nc)PcGc-m(nc, fc) "H fmcta - co (1) Pco G c o _ m (1, fc) + /mcfa-co(2)/3coGco-m(2, fc) Appendix B. Gebhart Factor Equations Derivation 111 + + fmcta-co{nCo)pcoGCo-m(nCO, fc) +fmcta-»{l)p,G,-m(l,k) +fmcta-,(2)p.G,-m(2,k) + + fmcta-.(ns)p,Gt-m(ns, k) "f" /mcta—mc(l) PmcGmc—m (1,*) ~\~ fmcta—mc(2)/?Tnc Gmc—m (2,k) (nmc, fc) ~\~fmcta—mclai,]-)PmclaGmcla—m(l,, fc) ~\~ fmcta—mcla(2)pmciaGmcla—m(2, fc) + + fmcta-mcla(nla)PmclaGmcla-m(nla, fc) "T" fmcta—met a —m (fc) ~\~ fmcta-tcw (1) A'tciu Gjciu-m (fc) -rfmcta-tct(l)ptctGtct-m(k) + fmcta-»t{l)p,tG,t-m{l, k) + /mcto-«t(2)p,tG,t_m(2, fc) + + fmcta-.t{nst)P.tG.t-m{nst, fc) (B.73) Gtcw—m{k) = ftcw-m{k)em -rftcw-m (1 )Pm G m—m {l,k) + ftcw-m(2)pmG m—m (2,k) + + ftcw-m{nrn)pmGm-m{nrn, fc) T / ( » - e ( l ) / » . G e - m ( l ) k) + /tC U,-c(2)^cGc-m(2 ) fc) + + ftcw-c(nc)pcGc-m(nc, fc) -\~ftcw-co (1 )Pco Gco-m (1, + /tcu,-co(2)pCOGc0-m(2) fc) + + ftcw-co(nco)pcoGco-m{nco, fc) Appendix B. Gebhart Factor Equations Derivation 112 +ftcW-s{l)p.G.-m{l, fe) + ftcw-.{2)p.G.-m(2, k) + + ftcw-,{ns)p,G,-m(ns, k) "T'/tctu-mc(l)PmcG'mc_m(l, fe) "T'/tciu-mc(2)pmcGm (._m(2, fe) + + ftcw-mc{nmc)pmcGmC-m(nmc, k) ~tftcw—mcla (l)PmclaGmcla-m(l, k) (2)pmclaGmela-m(2, fe) + + ftcw-mcia(nla)pmciaGmcia-m(nla, fe) ~T ftew—meta (1 ~)pmcta Gmcta—m(k) +ftcw-tcw{l)ptcwGtcw-m(k) +ftcw-tct(l)PtctGtct-m(k) -\-ftcw-it{l)p»tGst-m (1, k) + ftcw-.t{2)p.tG,t-rn% fe) + + ftcw-,t{nst)p,tGst-m{nst, fe) (B.74) Gtct-m(k) = ftct-m{k)em +ftct-m{l)pmG m—m (l,k) + ftct-m(2)PmG m—m + + ftct-m{nm)pmGm-m{nm, k) +ftct-c(l)PcGc-m(l, fe) + ftct-c(2)pcGc-m(2, fe) + + ftct-c{nc)pcGc-m(nc, k) -rftct-co{^)pcoGco-m(l, fe) + ftct-co(2)pcoGco~m(2, fe) + + ftct-c0(nco)pcoGco-m(nco, fe) +ftct-.(l)p.G,-m(l, fe) + ftct-.(2)p,G,-m(2, fe) + + ftct-.(ns)p,G,-m(ns, fe) Appendix B. Gebhart Factor Equations Derivation 113 ~\~ ftct—mc (1 )Pmc Gmc—m (1,*) +ftct-mc('2)pmcGrnc-m(2, k) + + ftct-mc(nmc)pmcGmc-m(nmc, k) +/tct {l)PmclaGmcla-m{l, k) +ftct —mcla {2)pmclaGmcla-m(2, k) + + /let —mcla (nla)pmciaGmcia-m('n'la, k) —meta {^)PmctaGmcta —m (k) -tftct-tcw (l)Ptctu Gtcw-m (k) +ftct-tct (1 )ptct Gjct _ m ( fc ) +/*<*-,* ( l )p.tG . t _ m ( l,fc) + /tct-«t(2)p»tC7,f-m(2, A:) + + ftct-.t(nst)p,tG.t-m{nst, k) (B.75) G,t-m(i,k) = f,t-m(i,k)er, m—m + + f,t-m{i,nm)pmGm-m(nm,k) +f.t-c(i, l )p eGc-m ( l , fc) + /.t-e(», 2 )p c G c_ r n(2, fe) + + /«t-c(», n c ) p c G e _ T O ( 7 i c , fe) +f.t-co(i, l)pcoGco-m(l, k) + f.t-co(i,2)PcoG co—tn + + /,t-co(i, nco)pcoGco-m{nco, k) ' + / * - . ( * , l)PmC?.-m(l, *) + / r t _ . ( t , 2)PmG.-m{2, k) + + /,t-,(i, nm)pmG,-m{nm, k) -rf,t-mc{h l)pmcGmc-m (1,*) (2,*) Appendix B. Gebhart Factor Equations Derivation 114 + 4- f.t-mc{i,nmc)pmcGmc-m(nmc, fe) +f.t —mcla (i, l)PmclaG mcla—m (1,*) + fit—mcla{i) 2)p m c/ 0G T O (Ja-m(2, fe) + + fit—mcla (i, nla)pmciaGmcia-m{nla, fe) —meta meta Gmcta—m ( k ) +fit-tew(i)PtcwGtcw-m(k) -rfit-tct(i)PtctGtct-m(k) +fit-it{i, l)pmG.t-m(h k) 4 fit-.t(h 2)PmGtt-m(2, fe) + + f.t-it(i, nm)pmGtt-m{nm, fe) (B.76) Gs-m(i,k) = ft-m(i,k)em +fi-m{i, l ) / 3 m G m _ m ( l , fe) + fi-m{i, 2)pmG m—m (2,k) + + f.-m{i, nm)pmGm-m(nm, fe) + / . - e ( * \ l > c G c - m ( l , * ) + / . _ c ( » , 2)pcGc-m(2, fe) + 4- f.-c(i,nc)pcGc-m(nc, fe) + / . - « , ( » , l ) p C O G c 0 _ m ( l , fe) + f.-co{i, 2)pcoGco-m(2, fe) + 4 fi-eo(i, nco)pcoGeo-m{nco, fe) +/._.(», lKG._ m(l, As) + / . - . ( t , 2)p mG ._ m (2,fe) + + f.-e(i,nm)pmGt-m(nm,k) (1,*) ~V fi—mc('i 2 ) p m c G m c — m (2,*) 4 + f.-mc{i,nrnc)prncGmc-m(nrnc, fe) (1,*) 4"/j—mcio(*> 2 ) p m c f a G m c J a _ m ( 2 , fe) Appendix B. Gebhart Factor Equations Derivation 115 + + fs-mcia(i,nla)pmciaGmda-m{nla, fc) —met a mclo Gmcta—m ( fc ) (fc) +/ < -trf(Opt ctG t r f_ m(fc) +/.-.*(», l K G , . m ( l , fc) + / . _ , t ( i , 2>mt7,_ro(2, fc) + + ft-.t(i,nm)pmGt-m(nm,k) (B.77) Gco—m(*j fc) — /co-m(*i fc)em m—m m—m (2,fc) + + /co-m(i, nm)pmGTn-rn(nm, fc) +/eo-c(»\ l ) p e G e _ m ( l , fc) + /eo-c(t, 2)PcGe.m(2, fc) + 4- fco-c(i,nc)pcGc-m(nc,k) +/»-»(», l > c o G c o _ m ( l , fc) + /«,-»(*! 2)/)coGrco_m(2, fc) + + fco-co(i, nco)PcoGco-m(nco, fc) +/--.(*, l ) p . G . _ m ( l , fc) 4- /«,_.(<, 2>,G,_m(2, fc) + + fco-,{i,ns)P.G.-m(ns,k) 4"/eo-me(*> 1 )PmcGmc—m (1, fc) +/co-mc(i, 2)pmcG (2,fc) + 4- feo-mc{i, nTnc)PmcGmc-Tn{mnc, fc) "("/co—mdo(*'i l)/'mcJoG m c/ 0_ m(l, fc) —mcla (*» 2 ) p m c j a G m c i a _ m ( 2 , fc) + 4" /co-mc/a(i,^a)/'mcJaGfmcia-m(^a, fc) 4/co —mcto (OP mcto Gmct a—m ( fc ) Appendix B. Gebhart Factor Equations Derivation 116 •T fco-tcw(i)ptcwGtcw-m(k) -rfco-tct{i)ptctGtct-m(k) +fco-.t(i,l)p.tG.t-m(l,k) + fee—it (i,2)p.tG.t-m(2,k) + + feo-.t(i, nst)p$tGtt-m(nst, k) (B.78) These equations will be solved simultaneously for the Gebhart factors Gi-j where j repre-sents the melt element k and i represents each other surface. The system that is produced by these equations can be written in matrix form as follows: where [A] is a square matrix of size N * N where N = nm + nc + nmc + ns + nla + nco + nst + 3 if the crystal, crystal top and outer crucible wall are included. If the system is solved for the encapsulant instead of the melt, nm will be replaced by ne , and N will be reduced by twice times the number of elements of the crystal that are below the encapsulant surface. Otherwise, N will be either equal to ( nm + nc+nmc + nla + 3 ) or ( nm + nc + nmc + nla + 3 + ns + nst ). (x) is a column matrix of size N * 1 that represents the Gebhart factors, and [b] is a column matrix of size N * 1 which represents the right hand side of the equation.The matrix is solved using routines available on the VAX-VMS system. The matrices [A], (x) and [b] are given next. The matrix [A] is large and will Appendix B. Gebhart Factor Equations Derivation 117 therefore be written as a 8 * 8 matrix where each element represents a matrix. Al A2 A3 A4 A5 A6 A7 A8 Bl B2 B3 B4 B5 B6 B7 B8 CI C2 C3 C4 C5 C6 C7 C8 DI D2 D3 D4 D5 D6 D7 D8 El E2 E3 E4 E5 E6 E7 E8 Fl F2 F3 F4 F5 F6 F7 F8 Gl G2 G3 G4 G5 G6 G7 G8 HI H2 H3 H4 H5 H6 H7 H8 Bl, B2 ... B8, C1,C2 .. . C8, DI, D2 E8, Fl, F2 ... F8, Gl, G2 ... G8 and Hi, H2 ... H8 are as follows. Al = /f»-m(l> l)Pm — 1 /m-m(l,2)pm / m _ m ( l , nm)^ /m-m(2, l)pm /m-m(2, 2)/J m — 1 fm-m{2,nm)prl ^ fm-m(nm, l)Pm fm-m{nm,2)pm fm-m(nm,nm)pm — 1 ^ A2 = / m - e ( l , l > e /m-c(l,2> c / m _ e ( 2 , l ) p e / m _ e ( 2 , 2 ) p e fm-e(l,nc)pc fm-c(2,nc)pc K fm-c(nm,l)Pc fm-c(nm,2)pc fm-c(nm,nc)Pe } Appendix B. Gebhart Factor Equations Derivation 118 A3 = fm—mc(l) l)Pmc fm—me(l> 2)pn /m-mc(2, l)pmc /m-mc(2, 2)p„ nmcjp n fm—mc{2, Tl ^ fm-mc{nm,l)pmc fm-mc(™rn,2)pmc /m-rac(™, Timc)prnc ^ I A4 = fin—mcla (l,l)pmcla fm-mcla{l,2)Pmcla /m-mcJa(l, nla)pmcla /m-mcio(2, l)pmcia /m-mc/a(2, 2)p m c/ a fm—mcla{2, Tll(l)pmcla ^ fm-mcla{nm, I)pmcla fm-mcla(nm, 2)pmc/a /m-mc/a(""l, nla)pmcla J fm—mcta meta /m-fcuj(l)Ptcn) /m-tct(l)Ptct /m-mcta(2)pmcto /m-tcu;(2)ptcu; fm-tct (2)ptct A5 = ^ fm-mcta{nm)pmcta fm-tcw(nm)ptcw fm-tct{nm)ptct j / / m _ . ( l , l ) p . /m-.(l,2)p /m-,(l,n 5)p, ^ / m _ , ( 2,l)p, /m_,(2,2)p, / m _,(2,ns)p. A6 = K fm-.(nm,l)p. fm-.(nm,2)p, fm..(nm,ns)p, J Appendix B. Gebhart Factor Equations Derivation ( A 7 = / m - c o ( l , l)Pco / m - e o ( l , 2)pc /m-co(2, l)pco /m-co(2, 2)p c / m - c o ( l , W C 0 ) / J c fm-co{'2,nco)p c y fm-co(nm, l)pco fm-co{nm, 2)pco fm-co{nm,nco)pco ^ fm-,t{lA)P>t fm-3t{l,2)p.t fm-,t(l,nst)p,t ^ fm-.t(2,l)ptt fm-.t{2,2)p,t fm-,t{2,nst)p,t A 8 = ^ fm-,t(nm, l)p,t fm-it(nm,2)p,t fm-,t(nm,nst)p,t ) / c - m ( l , l )Pm / c - m ( l , 2 ) p m / c - m ( l , nm)pm fc-m(2,l)pm fc-m(2,2)pm /c_m(2,nm)pm Bl = ^ fc-m(nc, \)pm fc-m(nc,2)pm fc-m(nc, nm)pm J ' / c _ c ( l , l > c - l / c_ c(l , 2 > c / c _ c ( l , n c > c / c _ c ( 2,l)p c / c_ c(2,2)^-1 / c _ c (2,nc> c B 2 = ^ fc-c(nc,l)pc fc^c(nc,2)pc fe-e(nc,nc)pc - 1 J endix B. Gebhart Factor Equations Derivation / c - m c ( l > l)Pmc / c - m c ( l , 2 ) p m c / c - m e ( l , nmc)pme fc-mc{2,l)pmc fe-mc(2, 2)Pme /c-mc(2, nmc)pmc B3= ^ /c-mc(nc, l ) p m c fc-mc(nc, 2)pmc fc-mc{nc,nmc)pmc / c - m c i o ( l ) l)pmcZa /c-mtJo(l) 2 ) /J m c J a fc-Tncla{l-inla)Pmcla fc-mcla(2,l)pmcla fc-mcla(2, 2)/9 m c J a /c-mcia(2, nZa)/3 m c j a B4= ^ fc-mcla(nC, l)pmcla fc-mcla(nC,2)Pmcla fc-mcla(nC,nla)pmcla ^ Jc—mcta mcto fc-tcw{2)ptcw fc-tct{2)ptct B5 = ^ fc-meta(nc)pmcta fc-tcwi^^Ptcw fc-tct(nc)ptct ^ ( fe-.{l,l)p. fc-s(l,2)ps fc-.(l,ns)p. ^ fe-.(2,l)p. /c_,(2,2>, fc-.(2,ns)p. B6 = ^ fc-.(nc, l)p, fc.,(nc,2)p, fc-,(nc,ns)p. ) Appendix B. Gebhart Factor Equations Derivation 121 /c-co(l, l)Pco fc-co(l,2)Pco }c-co{^inco)pco fe-co(2,l)pco /c-co(2,2)pco fc-co(2,nco)pco B7= ^ fc-co(nc,l)pco fc-co{nc,2)pco fc-co(nc,nco)pco ( / c _ , t ( l , l)Plt fc-.t{l,2)p,t fc-.t{L,nst)ptt ^ fc-,t{2,l)p.t fc.tt(2,2)pst fc_.t(2,nst)pet B 8 = ^ fc-»t(nc,l)ptt fc-tt(nc,2)ptt fc-$t{nc,nst)p,t j /mc-m(l, l)/?m /mc-rn(l, 2)pm fme-m(l, nm)Pm fmc-m{2,l)Pm /mc-m(2, 2)pm fmc-m{2, Tim)pm CI = ^ fmc—m ( nmc, l)pm fm.c-m(nmc,2)pm / m c - m (nmc, nm)pm J ( /mc-c(l,l)/»c /mc-c(l ,2)^ c fmc-c(^,nc)pe fmc-c(2,l)pc fmc-c(2,2)pc fmc-c(2,Tlc)Pe C2 = ^ fmc-c(nmc, \)pc fmc-c(nmc,2)pc fTnc-c{nmc,nc)pc J Appendix B. Gebhart Factor Equations Derivation 122 / C3 = /mC-mc(l, l)Pmc — 1 /mc-mc(l, 2)^mc /mc-mc(l. nmc)P„ fmc—mc(2, 1)/Jmc /mc-tnc(2, 2)^Jmc 1 fmc—mc{2)TlTnc)pTl ^ fmc-mc(nmc)l)pmc fmc-mc(nmc, 2)pmc fmc-me{nmc,nmc)prnc — 1 fmc—mcla [l>l)Pmcla fmc-mcla(l,2)pmcla /me-mela(l| ttfo)PmeJa | fmc—mcla (2,l)pmcla fmc—mcla (2,2)Pmcla fmc—mcla C4 = ^ fmc-mclainmc, l)pmda fmc-mda{nmc,2)pmcla fmc-mda(nmc,nla)pmcla j I C5 = fmc—mcta fmc-tcwiX)ptcw fmc-tct(l)Ptct fmc-mcta{2)pmcta fmc-tcw(2)ptcw fmc-tct{2) ptct ^ fmc-mcta(nmc)pmcta fmc-tcw{nmc)ptcw fmc-tct{nmc)Ptct j ( C6 = / m c - , ( l , l)p, fmc-.{l,2)p, fmc-.(2,l)P. fmc-.(2,2)P, fmc-,{l,ns)p, fmc-.(2,ns)p. ^ fmc-,(nmc, l)p, fmc-s(nmc,2)p, fmc-t(nmc,ns)pt / Appendix B. Gebhart Factor Equations Derivation i C7 = /mc-co(l) 1)P c /me-eo(2, l ) p c fmc-co(l, 2)p c /me-co(2, 2)p c c c nmc, l)pco fmc-eo(nmc,2)pco fmc^co(nmc)nco)pco J ( C8 = /mc-«t(l, l)p«t /mc-*t(l ,2)p 4f /mc-»t(2, l)p,t /mc-«t(2, 2)p,t /mc-«t(l,n5t)/J,t fmc-.t{2,nst)p at ^ fmc-st(nmc,l)p,t fmc-tt{nmc)2)p,t fmc-,t(nmc,nst)p3t ) ( DI = /mcia-m(l, l)Pm /mda-m(l, 2)^ (2,1> (2,2> /mcia-m(2,nm)/> n ^ fmda-m.(nla, l)pm fmda-m(nla,2)pm fmcla-m(nla,nm)pm j I D2 = /mdo-c(l) l)Pc fmcla-c(2, l)pc /mcia-c(l, 2)/3 c /mc/o-c(l, nc)pc fmcla-c{2,2)pc fmcla-c{2,nc)pc ^ fmcia-c(nla, \)pc fmcia-c(nla,2)pc fmda-c(nla,nc)pc t Appendix B. Gebhart Factor Equations Derivation 124 D3 = n /me/a-mc(2, l)p n fmcla—mc n fmcla—mc n fmcla—mc (2,nmc)p K fmela-me(nla,l)pme fmda-mc{nlat2)pme fmda-mc(nla,nmc)pmc ( D4 = fmcla—mcla fmcla—mcla{p"i ^ )Pmcla fmcla—mcla (2,2)pmcia — 1 ^ fmcla-mcla{nla, l)pmcla fmda-mda{'nla,2)prncla fmda-mda{^,flla)pmcla fmda-mda{2, nla)pmcla fmda-mda{nla,nla)pmcla — \j fmda-mcta{X)pmcta fmda-tcwiX)Ptcw / m c i a - t c f ( l ) P t e t fmda (2)p mcta /mcfa-tcui(2)ptcu> fmda-tct{2) Ptct D5 = ^ fmcla—mcta (nla)pmcta fmda-tcw(nla)ptcw fmcia-tct{nla)ptct J ^ / m c l o - « ( l , l ) p « / m c / o - » ( l ) 2 ) p , fmda-s{l,ns)P. ^ fmcla—M (2,2)p, fmda-.(2,ns)p. D6 = ^ fmda-,(nla,\)p, fmda-,(nla,2)p, fmcia-,{nla,ns)p, ) Appendix B. Gebhart Factor Equations Derivation ( D 7 = fmcla—co c fmcla-co{2, 1)P c /mdo-co(l,2)/) c /mc/a-co(2, 2)pc fmcla—co (l,nco)p c • •• fmcla-co(2,nco)pc ^ fmcla-co(nla,l)pco fmcla-co{nla,2)pco ••• fmcla-co(nla,nco)pco j ( D 8 = fmcla-tt{l, 1)P«1 /mcio-«t(2, l)p,{ /mda-«t(l,2)p,t fmcla->t(l, nst)Pit fmcla-tt(2,2)p,t fmcla-tt(2, nst)p,t ^ fmcla-tt(nla,l)ptt fmcia-,t{nla,2)p,t fmcia-»t(nla,nst)pat ) E l = fmcta-m{X)pm fmcta-m(2)pm fmcta-m(nm) Pm ftcw-m{l)Pm ftcw-m{2)pm ftcw-m{™m)pm ^ ftct-m(l)pm ftct-m{2)pm ftct-m(™n)pm j E2 fmcta-c(l)pc fmcta—c {2)pc fmcta-c{nc)pc ftcw-c(l)Pc ftcw-c(2)pc ftcw-c{nc)pc ^ ftct-c[X)pc ftct-c{2)pc ftct-c{nc)pc j E3 = ^ fmcta—mc{X)Pn ftcw—mc (1 )Pm< ^ ftct-mc{X)pmc fmcta—mc(2)pmc fmcta—mc{nmc)pTnc ftcw—mc(2)pmc ftcw—mci, nmc)p mc I ftct-mc{2)pmc ftct-mc{nmc)pmc j Appendix B. Gebhart Factor Equations Derivation 126 E 4 = fmcta—mclai^) Pmcla fmcta—mcla[2*)pmda ••• fmcta—mcla(jll&)pmcla ftew—mcla(X)P mcla ftew—mcla$)P mcla ftew—mcla (nla)pmcia ^ ftct-mela{l)Pmcla ftct-mcla{2)Pmcla •• ftct-mcla(nla)Pmcla j E 5 = fmcta—mctaPmcta 1 fmcta—tew P tew fmcta— tctPtct ftew—meta Pmcta ftew—tew P tew 1 ftcw—tctPtct ^ ftct-metapmcta ftct-tcwPtcw ftct-tctptct ~ 1 y E 6 = ' fmcta-i(l)p, fmcta-,{1)pt fmcta-t(ns)p, ftcw-t(l)Ps ftcw-,{2)Pe ftcw-»{ns)pt ^ ftct-.(l)p, ftet-.{2)p, ftct-,{ns)p„ j E 7 = ^ fmcta-coiX)peo fmcta-coitypco fmcta-co{nco)Pco /teu>-eo(l)/>co ftcw-co{1)Pco ftew-eo(nco)pco \ ftct-co{l)Pco ftct-co(2)Pco ftet-co(nco)p CO ) E8 = ( fmcta-3t(l)p,t fmcta-,t{2)p,t fmcta-it{nst)ptt ftcw-Mt{l)Pst ftcw-,t(2)p,t ftcw-tt(nst)p,t ^ ftct-,t{1)pit ftct-,t{2)p.t ftct-,t(nst)p,t j Appendix B. Gebhart Factor Equations Derivation ( F l = / . _ m ( l , l > m /._ m(l , 2 )^ / . _ m ( 2,l)p m f.-m(2,2)p„ / ,_m(l,nm)/j f,-m(2,nm)p K f.-m(ns,l)pm /,_m(ns,2)pm f„_m(nsinm)pm J F2 = /.-c(l,l)Pc /.-c(l ,2)p e /.-c(2,l)pe /.-e(2,2)pe ^ f.-c(ns, l)pe f.-e(ns, 2)Pc f.-e(l,nc)pe f.-c(2, nc)Pc f,-c(ns,nc)pc j ( F3 = J i—mc Tl J *—me (2, l)p m c f.-mc{2,2)pT1 ft—mc ( l,nmc)p n (2,7imc)p ^ f.-mc(ns, l)pmc f,-mc(ns,2)pmc f,-mc{ns,nmc)pmc j ( F 4 = ft—mcia(l) l)Pmeia /#—mcia(l) 2)pmc/a -mcla (2, 2)pmcla (US, l)pmcla f,-mcla{nS,2)pmcla ft-mcla (1, nla)pmcla fs-mcla{2,nla)pmcia fi-Tncia{ns,nla)pmcia j Appendix B. Gebhart Factor Equations Derivation ft-mcta{l)Pmcta ft-tcwiX)ptcw fi-tctiX)Ptct J i—meta meta f»-tcw(2)ptcw ft-tct(2)ptct F5 = V ft-mcta(ns)Pmcta ft-tcw(ns)ptcw f,-tct{ns)ptct J / . _ , ( l , l ) p , - l /.-.(l , 2)p. f.-.{\,ns)p, / , - . ( 2,l>. / ,_,(2,2)p,-1 f,.,{2,ns)Pe F6 = . . . ^ fs-t(ns,l)pt ft-t(ns,2)p, f...(ns,ns)p, - 1 J /,_co(l, l ) p c o f,-co{^,2)pco f,-co(l,nco)pco /«-co(2, l)pco f,-co{2,2)pco f,-co{2,nco)pco F 7 = . . ^ fi-co(ns,l)pco ft-co(ns,2)pco f8-co(ns,nco)pco ^ l)p . t f,-.t(l,2)p,t f.-.t(l,nst)ptt ^ f,-,t(2,l)p,t f,.,t(2,2)p,t f._st(2,nst)p,t F8 ^ f.-it{ns,\)p.t fs-tt(ns,2)ptt f.-.t{ns,nst)p,t ) Appendix B. Gebhart Factor Equations Derivation 129 G l = /co-m(l» l)pm /co-m(l) 2)p„ /co-m(2, l ) p m /co-m(2, 2 ) ^ /co-m(l, nm)p n /co-m(2,nm)/3 ^ fco-m(nCO, l)pm fco-m(nCO, 2)pm fco-m(nCO, Tim)pm ( G2 = /co-c(l, l)Pc /co-c(l, 2)/Jc / c o _ c ( 2 , l > c / c o_ c(2,2> c /«,-e(l|nc)pc /co_c(2,nc)/)c ^ fco-c(nco, l)pc fco-c(nco,2)pc fco-c(nco,nc)pc ) ( G3 = /co—mc(l) l)pVnc fco—mc(1 j 2)pn fco—mc(2] l)Pmc fco—mc(2} 2)pr fco—mc(l j nmc)p /co-mc(2,nmc)/9 fl ^ fco-mc(nCO, l)pmc fco-mc(nCO, 2)pmc /co-rac(7lCO ,?imc) /3 m c ^ /co—mc/a(l> l)Pmc/a /co—mc/a(lj 2)/Jmc/o fco—mcla G4 = fco—mcla (2, l ) p m C / o fco—mcla (2, 2 ) p m c j a /co-mda(2, n/o)/?mc/a ^ fco-mcla{nCO, l)Pmcla fco-mcla(nCO,2)pmcla fco-mcla{nCO,nla)pmcla / Appendix B. Gebhart Factor Equations Derivation ( G5 = /co (1)P meta fco-tcuiiX) Ptcw feo-tetiX)ptct fco-mcta(2)pmcta fcc—tcw(2)Ptcw fco-tct(^)Ptet ^ fco-meta(nco)pmcta fco-tew{nco)ptcw fco-tet{nco)ptet j ^ feo-.(l,l)ps fco-,{l,2)p, fco.t(l,ns)p, ^ / c o _ . ( 2 , l > . / c o_,(2,2)p, feo-.{2,ns)p. G6 = v feo-,(nco,l)ps feo-,(nco,2)p, fco-.(nco,ns)pt J ( G 7 = /co—co(lj l)Pco 1 /co-co( l , 2)p CO fco-co{2} l)Pco /co-co(2,2)pco 1 ^ fco-co{nco, l)pco fco-co(nco,2)pco fco-eo(l,nC0)p c fco-co(2,nco)p c fc -co(nco,nco)pco — 1 j / G8 /co-«t(l, l)Pst /co-*t ( l , 2)/),t fco-»t(2, l)ptt fco-,t(2,2)p,t ^ fco-Mt(nco,l)p,t fco-Mt(nco,2)p,t fco-,t{l,nst)p,t fco-tt(2,nst)plt fco-tt{nco,nst)p,t j Appendix B. Gebhart Factor Equations Derivation 131 HI = / , t _ m ( l , l)pm fat-m(l, 2)p„ f.t-m(2,l)pm f,t-m{2,2)pn H2 = ^ f.t-e(nst,l)pc f.t-c(nst,2)pc n Tl ^ fst~m{nst,l)pm f,t-m(nst>2)pm f,t-m(nst,nm)pm ( fat-c(l, l)pc /.t-c(l,2)/>c fH_c(l,nc)pc ^ f.t-c{2,l)pc f,i-c(2,2)pc fat.c{2,nc)pc fst-e(nst,nc)pc ^ ( H3 = /«f-mc(l, lJPr fst-Tnc{2,l)P Ti TI ftt-mc(2,2)p Tl /*f—mc( l,nmc)p Tl /*t-mc(2, nmc)p Tl ^ f>t-mc(nst, \)pmc f,t-Tnc(nst,2)pmc f,t-Tnc{nst,nmc)pmc t I H4 = /«t-mcio(l> l)Pmc/a /*t—mcia(lj 2)p m c Zo —mcia (2, l ) p m C/o /»(-mcia(2, 2)/>mcia ^ ftt-mcla{nst, l)pmcla ftt-mcla(nst,2)pmcla /«t-mcio(l> nla)pmcia fst-mcia(2, nla)pmcla ftt-Tnclainstjnl^Pmcla ) Appendix B. Gebhart Factor Equations Derivation 132 -meta[^)Pmcta ftt-tcwiX)Ptcw f$t-tctiX)Ptct \ ftt-mcta{2)pmcta fst-tcw{2)Ptcw ftt-tct{2)ptct H5 = ^ f$t-mcta{nst)pmcta f,t-tcw(nst)ptcw ftt-tct(nst)ptct j ( /.t_,(l,l>. /.t-.(l,2)p /.t_.(l,nj)p. > /„_,(2,1>. /<t.J(2>2)p. /.t-.(2,Tw)p. H6 = ^ f,t-.(nst, l)p. f,t-,{nst,2)p. f.t-s(nst,ns)p, } /jt-co(l, l)Pco f$t-co(l,2)pco f,t~a>0-,nco)pco ftt-co{2,i)pco ftt-co(2>2)pco ftt-eo(2,nco)pco H7 = ^ ftt-co(nst, l)Pco f,t-co{nst,2)Pco f,t-co{nst,nco)pco j /.i-.t(l,l)/0#t - 1 /.t-»t(l,2)p,( f.t-,t{l,nst)p.t ^ ftt.,t(2,l)p.t f.t-,t(2,2)p.t - 1 f.t-.t(2,nst)p.t . H 8 = ^ f,t-it(ns,l)ptt ftt-it(ns,2)p,t ftt-tt{nst,nst)p,t - 1 ^ Appendix B. Gebhart Factor Equations Derivation 133 The column matrix (x), as stated earlier, represents the Gebhart factors that must be evaluated, (x) will depend on the particular surface under consideration. It will be given here for a general case and for one specific case, that being the melt. Therefore, in general, (x) will have the form: Gk-i,j \ Gk'j 1 For the specific case of the melt (x) will be a long matrix. Therefore, it will be divided to 8 smaller column matrices that are as follows : Appendix B. Gebhart Factor Equations Derivation ( - ) -x l x2 x3 x4 x5 x6 x7 v x 8 / where x l , x2, x8 are as defined below. x l = ^ C7m_m(l, fe) ^ Gm-m(2,k) "m—m (M) y Gm-m(nm, k) J 1 Gc_m(l,fe) ^ Gc_m(2,fc) ,x2 = ^ Gc.m(nc.fe) j (1,*) 1 Gmc-m(2, fe) ,x3 = ^ Gmc-m(nmc, k) J x4 = Gmcla—jn(11 fe) GmcJa-m(2, fe) ^ Gmcia-m(nla, fe) j ,x5 = "meta—m (fe) * Gtcu,-m(fe) ^ Gtct.m{k) J ,x6 = ^ ^ ( l . f e ) N G,.m(2,k) ^ G,-m(ns, k) J Appendix B. Gebhart Factor Equations Derivation x7 = Gco-m(l) fc) G«,-m(2,fe) ^ Gco-minco,*) j ,x8 = / G r t _ m ( l , * ) ^ G,t-m(2,fc) ^ G, t _ m (nst,A:) y The right hand side [b] is given by: b l b2 b3 b4 b5 b6 b7 b8 b l , b2 b8 are defined next for the the melt surface. b l = /m—m(1 > fc)^n /m—m(2, fc)6fi b4 = fm-m(nm, k)em / m c i a - m ( l ) fc)en fmcla—m{2, fc)Cri ,b2 = / c - m ( l ) fc)en / c_m(2,A:)eT 1 ,b3 = fmcla-m(nla, k) ,b5 / c _ m ( n c , fc)en fmcta—m ( fc ) /tctu-m(fc) em /tct-m(fc) c m fmc—m (1 > fc ) fmc—m(2, fc)6n ,b6 = fmc—m ( nrac, A:) / * — m ( l j fc)^m /« -m(2 , fc)Cm /«-m(n«|fc)Cn Appendix B. Gebhart Factor Equations Derivation 136 /.t -m ( l , k)cm f.t-m(2,k)em f3t-m{nst,k)em _ Since the surface considered here is the melt, then k will vary from 1 to nm. Therefore, the system will be solved nm times to determine the Gebhart factors for all melt elements. The matrix [A] is the same for the whole system. However, (x) and [b] are not. The fact that the matrix [A] is the same for all surfaces reduces the computation time required considerably. The method of solution for the other surfaces is similar to that for the melt and is described next. B.3 Gebhart factor equations for the crucible inner wall The equations for the crucible elements will be similar to the melt elements' equations. The'm' subscript will be replaced by 'c' in the ' G ; _ m ' variable and in the product of em and fi-m- Also, k will vary from 1 to nc, which is the number of elements of the crucible. The subscript 'm' in the matrices (x) and [b] will be replaced by 'c\ B.4 Gebhart factor equations for the middle cylinder wall The'm' subscript in this case will be replaced by 'mc' in lGi-m\ and in the product of em and /,_,„. k will vary from 1 to nmc - the number of elements of the middle cylinder. The subsript 'm' will be replaced by 'mc' in the matrices (x) and [b]. b7 = fco—m(l> k)en /co-m(2, k^Cf, ,b8 = fco-m(nco, k) Appendix B. Gebhart Factor Equations Derivation 137 B .5 G e b h a r t f a c t o r e q u a t i o n s for t h e c r y s t a l a n d t h e c r y s t a l t o p The'm' subscript in this case will be replaced by's' for the crystal case and by 'st' for the crystal top case in 'Gi-m\ and in the product of and / , _ m . k will vary from 1 to - the number of elements of the crystal, and from 1 to nst - the number of elements of the crystal top when calculating the crystal top factors, 'm' will be replaced by's' or by 'st' in the matrices (x) and [b]. B.6 G e b h a r t f a c t o r e q u a t i o n s for t h e c r u c i b l e o u t e r w a l l The'm' subscript in this case will be replaced by 'co' in '(?»_„,' and in the product of tm and fi-m- k will vary from 1 to nco - the number of elements of the outer crucible wall. The subscript'm' will be replaced by 'nco' in the matrices (x) and [b]. B.7 G e b h a r t f a c t o r e q u a t i o n s for t h e m i d d l e c y l i n d e r l o w e r a n n u l u s The'm' subscript in this case will be replaced by 'mcla' in '(j;_m' and in the product of tm and fi-m- k will vary from 1 to nla - the number of elements of the lower annulus. In the matrices (x) and [b], the subscript 'm' will be replaced by 'mcla'. B.8 G e b h a r t f a c t o r e q u a t i o n s for t h e m i d d l e c y l i n d e r t o p a n n u l u s , t o p c y l i n -d e r w a l l , t o p c y l i n d e r t o p The 'm' subscript in this case will be replaced by 'mcta', 'tew' and 'tct' in 'C?i_m' and in the product of em and / i _ m for the middle cylinder top equations, the top cylinder wall equations and the top cylinder top equations, respectively. The variable k will be equal to 1 since there is only one element on each surface. For the case of calculating the Gebhart factors for the top chamber that is composed of mcta, tew, tct, the subscript Appendix B. Gebhart Factor Equations Derivation 138 'm' in both matrices (x) and [b] will be replaced by 'meta ' when calculating the middle cylinder top factors and so on for the other two surfaces. B.9 Gebhart factor equations for the encapsulant The encapsulant surface is considered for the case when the radiation emitted from the melt is at a wave length greater than 2/zm. Since at that wavelength, the encapsulant is considered opaque. Therefore, the Gebhart factor equations are written in the same manner as those for the melt except that the'm' subscript is replaced by 'e' in all the equations and matrices. The subscript'm' in the other surfaces' equations and matrices is also replaced by 'e'. The crucible and crystal elements that are below the encapsulant surface will not be considered as part of the system. Appendix C Programme Flow Chart C . l General This appendix gives the flow chart of the programme that is used in calculating the configuration factors and the Gebhart factors. The flow chart is divided tofour sections and given in figures C.44 to C.47. Figure C.44 gives the flow chart of the main section of the programme in which the case is chosen, the surfaces are divided to elements and then the main subroutine is called for encapsulant surface (cases 12, 22 or 32). Depending on which case is calculated the main subroutine is called again for melt surface calculations (cases 1, 2 or 3). The numbers 1, 2 and 3 refer to the cases described in chapter-3. Figures C.45 to C.47 show the flow chart of the main subroutine. 139 Appendix C. Programme Flow Chart START Define lease 12, 22 or 32 Define Ni and Xi I Initialize heat transfer summation to zero Call Main Subroutine Figure C.45: Main programme flow chart. Appendix C. Programme Flow Chart START ~ r ~ " Enter constants' values Calculate surface areas temperature and energy fraction Calculate H(i-j) Calculate Configuration Factor Depending on case Add Configuration Factor for each surface rint Configuration Factor summation for each surface. Figure C.46: Main subroutine flow chart. Appendix C. Programme Flow Chart Begin Gebhart Factor Calculation i Enter Elements of Matrix 'A' Factor Matrix 'A' by Gaussian Elimination Enter Column Matrix 'b' Solve for Gebhart Factor of surface V Calculate the product Xk - eps*Ak*Tk*Gik*fi Calculate the Sums Yk - Yk •• Xk«sigma Zk " Zk • Xk Calculate the sum of the Gebhart Factors for each surface element Print Gebhart Factor summation of each surfaces Figure C.47: Continuation of main subroutine flow chart. Appendix C. P r o g r a m m e Flow Chut .yes Calculate Heat Transfer and the Equivilent Ambient Temperature of Each Surface Element I *rint: A r e a . "TemPk Equiv. Temp. Q (abs), Q (emit) and Q (net) of each element Calculate the Sums of the Heat Transferred from and to each surface Calculate the Net Heat Transfer from the System Figure C.48: Continuation of main subroutine flow chart. Appendix D Configuration and Gebhart Factors Calculations Results D.l General This appendix gives the results of the configuration factor summation for each element calculated for cases II and III for the simplified case of one—surface enclosure. The results for case-I are given in appendix A. The Gebhart factor summation is one for all elements. Tables D.22 and D.23 give the summation for case-II from the melt, crucible, crystal and crystal top elements, and the ambient to the melt, crucible, crystal, crystal top and ambient surfaces. Tables D.24 and D.25 gives the same information but with replacing the melt surface by the encapsulant. Tables D.26, D.27 and D.28 give the same information stated in the previous para-graph for case-Ill (exchange with the melt). Tables D.29, D.30 and D.31 fist the data for exchange with the encapsulant for case-Ill . 144 Appendix D. Configuration and Gebhart Factors Calculations Results 145 Table D.22: Configuration factor summation from surface j to surface i for case-II (ex-change with the melt). SURFACEj SURFACE i # Melt Crucible Crystal Wall Crystal Top Enclosure Melt 1 — 0.4021 0.1818 — 0.4162 2 — 0.4668 0.2422 — 0.2910 3 — 0.5151 0.1915 — 0.2934 4 — 0.5586 0.1516 — 0.2898 5 — 0.6950 0.1207 — 0.1843 Crucible 1 0.4556 0.2050 0.2543 0.0000 0.0851 2 0.3728 0.2027 0.3251 0.0000 0.0994 3 0.3020 0.2269 0.3540 0.0000 0.1170 4 0.2448 0.2477 0.3685 0.0000 0.1390 5 0.1998 0.2651 0.3685 0.0000 0.1666 6 0.1646 0.2797 0.3540 0.0000 0.2018 7 0.1370 0.2911 0.3251 0.0000 0.2468 8 0.1150 0.3266 0.2543 0.0000 0.3040 9 0.0974 0.2682 0.2050 0.0164 0.4130 10 0.0831 0.2740 0.1307 0.0459 0.4663 Table D.23: Configuration factor summation from surface j to surface i for case-II (ex-change with the melt). SURFACEj SURFACE i # Melt Crucible Crystal Wall Crystal Top Enclosure Crystal 1 0.3198 0.5269 — — 0.1534 2 0.2643 0.6754 — — 0.0603 3 0.2096 0.7434 — — 0.0470 4 0.1601 0.7870 — — 0.0528 5 0.1208 0.8081 — — 0.0711 6 0.0913 0.8081 — — 0.1006 7 0.0696 0.7870 — — 0.1434 8 0.0537 0.7437 — — 0.2029 Crystal Top 1 — 0.0390 — — 0.9610 2 — 0.0414 — • — 0.9586 3 — 0.0469 — — 0.9531 4 — 0.0571 — — 0.9429 Ambient 1 0.0706 0.1493 0.0277 0.0792 0.6732 Appendix D. Configuration and Gebhart Factors Calculations Results 146 Table D.24: Configuration factor summation from surface j to surface i for case-II (ex-change with the encapsulant). SURFACE j SURFACE i # Encapsulant Crucible Crystal Crystal Top Enclosure Encapsulant 1 — 0.3153 0.1772 — 0.5075 2 — 0.3780 0.2180 — 0.4040 3 — 0.4349 0.1589 — 0.4062 4 — 0.4930 0.1156 — 0.3913 5 — 0.6133 0.0850 — 0.3017 Crucible 4 0.4556 0.1834 0.2220 — 0.1390 5 0.3728 0.1813 0.2793 — 0.1666 6 0.3020 0.2068 0.2894 — 0.2018 7 0.2448 0.2291 0.2793 — 0.2468 8 0.1998 0.2742 0.2220 — 0.3040 9 0.1646 0.2240 0.1820 0.0164 0.4130 10 0.1370 0.2365 0.1142 0.0495 0.4663 Table D.25: Configuration factor summation from surface j to surface i for case-II (ex-change with the encapsulant). SURFACEj SURFACE i # Melt Crucible Crystal Wall Crystal Top Enclosure Crystal 4 0.3198 0.4940 — — 0.1863 5 0.2643 0.6295 — — 0.1061 6 0.2096 0.6788 — — 0.1116 7 0.1601 0.6954 — — 0.1445 8 0.1208 0.6788 — — 0.2004 Crystal Top 1 — 0.0390 — — 0.9610 2 — 0.0414 — — 0.9586 3 — 0.0469 — — 0.9531 4 — 0.0571 — — 0.9429 Ambient 1 0.0977 0.1292 0.0249 0.0792 0.6690 Appendix D. Configuration and Gebhart Factors Calculations Results 147 Table D.26: Configuration factor summation from surface j to surface i for case-Ill (exchange with the melt). SURFACEj SURFACE i # Melt Crucible Crystal Wall Crystal Top Enclosure Melt 1 — 0.4937 0.1849 — 0.3213 2 — 0.5630 0.2617 — 0.1753 3 — 0.6047 0.2220 — 0.1733 4 — 0.6366 0.1912 — 0.1722 5 — 0.7446 0.1671 — 0.0883 Crystal Top 1 — — — — 1.0000 Ambient 1 0.0434 0.1222 0.3274 0.0833 0.4237 Table D.27: Configuration factor summation from the crucible elements to the other surfaces for case-III (exchange with the melt). Crucible SURFACE i Element # Melt Crucible Crystal Wall Crystal Top Enclosure 1 0.4556 0.2370 0.2763 — 0.0310 2 0.3728 0.2382 0.3543 — 0.0347 3 0.3020 0.2655 0.3935 — 0.0389 4 0.2448 0.2888 0.4225 — 0.0439 5 0.1998 0.3074 0.4432 — 0.0496 6 0.1646 0.3214 0.4576 — 0.0563 7 0.1370 0.3311 0.4677 — 0.0642 8 0.1150 0.3368 0.4748 — 0.0734 9 0.0974 0.3387 0.4798 — 0.0842 10 0.0831 0.3368 0.4833 — 0.0969 11 0.0713 0.3311 0.4858 — 0.1119 12 0.0615 0.3214 0.4874 — 0.1296 13 0.0534 0.3074 0.4886 — 0.1507 14 0.0465 0.2888 0.4892 — 0.1755 15 0.0407 0.2655 0.4894 — 0.2044 16 0.0358 0.2382 0.4890 — 0.2370 17 0.0316 0.2370 0.4806 — 0.2508 Appendix D. Configuration and Gebhart Factors Calculations Results 148 Table D.28: Configuration factor summation from the crystal wall elements to the other surfaces for case-Ill (exchange with the melt). Crystal SURFACE i Element # Melt Crucible Crystal Wall Crystal Top Enclosure 1 0.3198 0.5475 — — 0.1327 2 0.2643 0.7025 — — 0.0332 3 0.2096 0.7795 — — 0.0110 4 0.1601 0.8358 — — 0.0041 5 0.1208 0.8749 — — 0.0077 6 0.0913 0.9009 — — 0.0131 7 0.0696 0.9173 — — 0.0201 8 0.0537 0.9262 — — 0.0290 9 0.0419 0.9291 — — 0.0290 10 0.0331 0.9262 — — 0.0407 11 0.0265 0.9173 — — 0.0562 12 0.0215 0.9009 — — 0.0776 13 0.0176 0.8749 — — 0.1075 14 0.0146 0.8358 — — 0.1496 15 0.0122 0.7795 — — 0.2083 16 0.0103 0.7025 — — 0.2872 17 0.0087 0.5475 — — 0.4437 18 0.0076 0.4439 — — 0.5484 19 0.0070 0.3060 — — 0.6873 20 0.0059 0.2319 — — 0.7622 21 0.0052 0.1740 — — 0.8208 22 0.0047 0.1304 — — 0.8649 23 0.0042 0.0982 — — 0.8976 24 0.0037 0.0747 — — 0.9216 25 0.0034 0.0575 — — 0.9392 26 0.0030 0.0447 — — 0.9522 27 0.0027 0.0352 — — 0.9620 28 0.0025 0.0281 — — 0.9694 Appendix D. Configuration and Gebhart Factors Calculations Results 149 Table D.29: Configuration factor summation from surface j to surface i for case-Ill (exchange with the encapsulant). SURFACE j SURFACE i #• Encapsulant Crucible Crystal Crystal Top Enclosure Encap 1 — 0.4661 0.1849 — 0.3490 2 — 0.5334 0.2611 — 0.2055 3 — 0.5766 0.2210 — 0.2024 4 — 0.6114 0.1898 — 0.1988 5 — 0.6856 0.1653 — 0.1491 Crystal Top 1 — — — - — 1.0000 Ambient 1 0.0525 0.1152 0.3266 0.0834 0.4223 Table D.30: Configuration factor summation from the crucible elements to the other surfaces for case-Ill (exchange with the encapsulant). SURFACE j SURFACE i Element # Encapsulant Crucible Crystal Crystal Top Enclosure 4 0.4556 0.2245 0.2760 0.0000 0.0439 5 0.3728 0.2237 0.3539 0.0000 0.0496 6 0.3020 0.2486 0.3930 0.0000 0.0563 7 0.2448 0.2691 0.4219 0.0000 0.0642 8 0.1998 0.2844 0.4425 0.0000 0.0734 9 0.1646 0.2944 0.4569 0.0000 0.0842 10 0.1370 0.2993 0.4668 0.0000 0.0969 11 0.1150 0.2993 0.4738 0.0000 0.1119 12 0.0974 0.2944 0.4786 0.0000 0.1296 13 0.0831 0.2844 0.4819 0.0000 0.1507 14 0.0713 0.2691 0.4841 0.0000 0.1755 15 0.0615 0.2486 0.4855 0.0000 0.2044 16 0.0534 0.2237 0.4859 0.0000 0.2370 17 0.0465 0.2245 0.4781 0.0000 0.2508 Appendix D. Configuration and Gebhart Factors Calculations Results 150 Table D.31: Configuration factor summation from the crystal elements to the other surfaces for case—III (exchange with the encapsulant). SURFACE j SURFACE i Element # Encapsulant Crucible Crystal Crystal Top Enclosure 4 0.3198 0.5427 — — 0.1375 5 0.2643 0.6964 — — 0.0393 6 0.2096 0.7717 — — 0.0188 7 0.1601 0.8257 — — 0.0142 8 0.1208 0.8616 — — 0.0176 9 0.0913 0.8832 — — 0.0254 10 0.0696 0.8933 — — 0.0371 11 0.0537 0.8933 — — 0.0530 12 0.0419 0.8832 — — 0.0749 13 0.0331 0.8616 — — 0.1053 14 0.0265 0.8257 — — 0.1478 15 0.0215 0.7717 — — 0.2068 16 0.0176 0.6964 — . — 0.2860 17 0.0146 0.5427 — — 0.4427 18 0.0125 0.4399 — — 0.5476 19 0.0107 0.3027 — — 0.6866 20 0.0093 0.2292 — — 0.7615 21 0.0081 0.1717 — — 0.8203 22 0.0071 0.1284 — • — 0.8645 23 0.0062 0.0966 — — 0.8972 24 0.0055 0.0733 — — 0.9212 25 0.0049 0.0562 — — 0.9389 26 0.0044 0.0437 — — 0.9520 27 0.0039 0.0343 — — 0.9618 28 0.0035 0.0273 — — 0.9692 Appendix E Stainless Steel Emittance E.l General This appendix presents the emittance of various stainless steels — oxidized, polished and cleaned as obtained from [33]. Figures E.48 to E.50 show the analysed grouping of the data which give a quick reference range for the emittance of a certain type of steel for a temperature or temperature range. Figure E.51 gives the exact curves obtained from experimental data. The specification and data tables which are not included here for the sake of brevity, describe the specimen condition and exact experimental values of the emittance. 151 !~' Temperature (k) CO ~ Temperature (K) Figure E.52: Normal total emittance of various stainless steels. Bibliography [1] A.K. Willardson, Advances in GaAs Crystal Growth, Advanced Materials and Processes, 1986. [2] E.P.A. Metz, R.C. Miller and R.Mazelsky, A Technique for Pulling Single Crystals of Volatile Matter, J. of Applied Physics, 33-1962-6. [3] J.B. MuUin et al, Liquid Encapsulation Crystal Pulling at High Pressures, J. of Crystal Growth, 3,4(1968)281-285. [4] R. Lamprecht et al, Experiments on Bouyant, Thermocapillary and Forced Convection in Czochralski Configuration, J. of Crystal Growth, 65(1983)143-152. [5] W.E. Langlois, Conservative Differencing Procedures for Rotationally Sym-metric Flow with Swirl, Computer Methods Appl. Mech. Eng., 25(1981)315. [6] W.E. Langlois, Digital Simulation of Czochralski Bulk Flow in Microgravity, J. of Crystal Growth, 48(1980)25. [7] W.E. Langlois, A Parametric Sensitivity Study for Czochralsi Bulk Flow in Silicon, J. of Crystal Growth, 56(1982)15. [8] K.terashima and T.Fukuda, A New Magnetic-Field Applied Pulling Apparatus for LEC GaAs Single Crystal Growth, J. of Crystal Growth, 63(1983)423-425. [9] J. Osaka et al, Growth of Semi-Insulating GaAs Single Crystal by LEC Method, Review of the Electrical Communications Laboratories, Vol 33 #1, 156 Bibliography 157 1985. [10] H. Kohda et al, Crystal Growth of Completely Dislocation- Free and Striation-Free GaAs, J. of Crystal Growth, 71(1985)813-816. [11] Ozawa and Fukuda, In-Situ Observation of LEC GaAs Solid-Liquid Interface with Newly Developed X-Ray Image Processing system, J. of Crystal Growth, 76-1986-323. [12] K. Kakimoto et al, Direct Observations by X-Ray Radiography of Convection of Boric Oxide in the GaAs Liquid Encapsulated Czochralski Growth, J. of Crystal Growth, 94(1989)405-411. [13] P.A. Ramachandran and M.P. Dudkovic, Simulation of Temperature Dis-tribution in Crystals Grown by Czochralski Method, J. of Crystal Growth, 71(1985)399-408. [14] L.J. Atherton et al, Radiative Heat Exchange in Czochralski Crystal Growth, J. of Crystal Growth, 84(1987)57-78. [15] R.K. Srivastava et al, Interface Shape in Czochralski Grown Crystals - Effect of Conduction and Radiation, J. of Crystal Growth, 73-1985-487. [16] S. Motakef and A.W. Witt, Thermoelastic analysis of GaAs in LEC growth cofiguration, J. of Crystal Growth, 80-1987-37. [17] J.J. Derby and R.A. Brown, Thermal Capillary Analysis of Czochralski and Liquid Encapsulant Czochralski crystal Growth, J. of Crystal Growth, 75-1986-227. Bibliography 158 [18] T.W. Hicks, Fluid Motion in the encapsulant region of the LEC Growth Sys-tem, J. of Crystal Growth, 84-1987-598. [19] A.S. Jordan, An Evaluation of the thermal and elastic Constants Affecting GaAs Crystal Growth, J. of Crystal Growth, 49-1980-631. [20] P.D. Thomas et al, Dynamics of Liquid-Encapsulated Czochralski Growth of Gallium Arsenide: Comparing Model with Experiment, J. of Crystal Growth, 96(1989)135-152. [21] A.G. Ostrogorsky et al, Infrared Absorbance of B2O3 at Temperatures to 1250°C, j ; of Crystal Growth, 84-1987-460. [22] C. Schvezov et al, Mathematical Modelling of the Liquid Encapsulated Czochralski Growth of Gallium Arsenide, J. of Crystal Growth, 84(1987)2-218. [23] A.S. Jordan et al, Bell System Tech. J., 59(1980)593. [24] I. Grant, B. Smith and D. Rumsby and R.W. Ware, Cambridge Instrument Co., Rustat Road, Cambridge, UK. private communication of author, 1985. [25] P. Sabhapathy and M. E. Salcudean, Natural Convection in Liquid Encapsu-lated Czochralski Growth of Gallium Arsenide, Canadian Journal of Chemical Engineering, Volume 68, 243-249, April, 1990. [26] M.E. Salcudean et al, Numerical Study of Free and Forced Convection in the LEC Growth of GaAs, J. of Crystal Growth, 94(1989)522-526. [27] P. Sabhapathy and M.E. Salcudean, Numerical Analysis of Heat Transfer in LEC Growth of GaAs, J. of Crystal Growth^ 97(1989)125-135. Bibliography 159 [28] J.J. Derby, R.A. Brown, F.T. Geyling, A.S. Jordan, and G.A. Nikolakopoulou, Finite-Element Analysis of a Thermal Cappilary Model for LEC Czochralski Model, J. Electrochem. Soc. 132(1986)470. [29] J.J. Derby and R.A. Brown, Thermal Capillary Analysis of Czochralski Crys-tal Growth, J. of Crystal Growth, 74(1986)605. [30] Napolitano et al, Viscosity and Density of Boron Trioxide, J. of American Ceramics Society Vol48 no. 12 ,1965, pp.613-616. [31] G.K. Greffield and A.J. Wickens, J. of Chemical Engineering Data, 20-1975-223. [32] A.S. Jordan, Determination of the Total Emittance of n-type GaAs with Ap-plication to Czochralski Growth, Journal of Applied Physics, 51(4)2218-2227, April,1980. [33] Y.S. Touloukian et al, Thermophysical Properties of Matter, Volumes # 2, 3, 6, 7, 8, 9 & 11 (IFI/Plenum, New York, 1970). [34] A. Goldsmith, T.E. Waterman & H.J. Hirschhorn, Handbook of Thermophys-ical Properties of Solid Materials, Vol. 1-5, the Macmillan Company, New York, 1961. [35] Thermophysical Properties Research Literature Guide, Y.S. Touloukian, J.K. Gerritsen & N.Y. Moore, Vol. 1-3, Plenum 1967; Y.S. Touloukian, J.K. Ger-ritsen & W.H. Shafer, Supplement-I Vol. 1-6, IFI/Plenum 1973; J.K. Gerrit-sen, V. Ramdas & T.M. Putnam, Supplement-II Vol 1-6, IFI/Plenum 1977. [36] G.O.Meduoye et al, The Minimization of Thermal Stresses During the Growth of GaAs Crystals, J. of Crystal Growth, 88(1988)397-410. Bibliography 160 [37] M.J. Crochet at al, Finite Element Simulation of Czochralski Bulk Flow, J. of Crystal Growth, 65(1983)153-165. [38] R.A. Lemons and M.A. Bosch, Microscopy of Si Films During Laser Melting, Applied Physics Letters, vol.40, No.8, 1982. [39] R. Braunstein and E.O. Kane, The Valence Band Structure of the III-V Com-pounds, J. Phys. Chem. Solids, vol.23, ppl423-1431, 1962. [40] H.C. Hottel, Radiant-heat Transmission, in William H. McAdams (ed.), Heat Transmission, 3rd ed., chap. 4, McGraw-Hill Book Company, 1954. [41] G. Poljak, Analysis of Heat Interchange by Radiation between Diffuse Sur-faces, Tech. Phys. USSR, voll, nos 5, 6, pp555-590, 1935. [42] M. Jakob, Heat Transfer, vol. II, John Wiley k Sons, Inc., New York, 1957. [43] R. Siegel and J.R. Howel, Thermal Radiation Heat Transfer, second edition, Hemisphere Publishing Corporation, 1981. [44] J.R. Howel, A Catalog of Radiation Configuration Factors, McGraw Hill, Inc, 1982. [45] S.V. Patankar, Numerical Heat Transfer and Fluid Flow, McGraw-Hill (Hemisphere), New York, 1980. [46] G.D. Raithby, A Critical Evaluation of Upstream Differencing Applied to Problem Involving Fluid Flow, Computer Methods in Applied Mechanics and Engineering, 9(1976)75-103.
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Radiative heat transfer in gallium arsenide lec crystal...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Radiative heat transfer in gallium arsenide lec crystal pullers Bakeer, Muna 1990
pdf
Page Metadata
Item Metadata
Title | Radiative heat transfer in gallium arsenide lec crystal pullers |
Creator |
Bakeer, Muna |
Publisher | University of British Columbia |
Date Issued | 1990 |
Description | A numerical analysis of radiative heat transfer in a liquid encapsulant Czochralski gallium arsenide crystal puller is developed. The heat transfer and equivilent ambient temperature of each surface element are calculated using the Gebhart radiative model. The effective ambient temperature, to which each surface element is radiating, is found to vary indicating that assuming a constant ambient temperature for all surfaces (simplified radiative model) is incorrect. The importance of including the middle and top cylinders of the growth chamber in numerical analysis of radiative heat transfer in the system is evaluated in the study. The upper section could be replaced by one isothermal surface without significant change of the effective ambient temperature distribution. Fluid flow and heat transfer in the GaAs melt, crystal and encapsulant are calculated using a three dimensional axisymmetric finite difference code which includes the detailed radiative model. The mathematical modelling of the fluid and heat flow describes steady state transport phenomena in a three dimensional solution domain with latent heat release at the liquid/solid interface. The predicted flow and temperature fields using the detailed radiative model differ considerably from the predicted fields using the simplified model. The simplified model shows high axial and low radial temperature gradients in the crystal near the encapsulant region; the axial gradient decreases and the radial gradient increases with increasing distance from the encapsulant top. The detailed model shows a high radial temperature gradient in the crystal near the crystal-encapsulant-ambient junction and nearly flat isotherms in the top half of the crystal. |
Subject |
Gallium arsenide crystals |
Genre |
Thesis/Dissertation |
Type |
Text |
Language | eng |
Date Available | 2010-11-10 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
DOI | 10.14288/1.0098437 |
URI | http://hdl.handle.net/2429/29916 |
Degree |
Master of Applied Science - MASc |
Program |
Mechanical Engineering |
Affiliation |
Applied Science, Faculty of Mechanical Engineering, Department of |
Degree Grantor | University of British Columbia |
Campus |
UBCV |
Scholarly Level | Graduate |
AggregatedSourceRepository | DSpace |
Download
- Media
- 831-UBC_1991_A7 B34.pdf [ 7.55MB ]
- Metadata
- JSON: 831-1.0098437.json
- JSON-LD: 831-1.0098437-ld.json
- RDF/XML (Pretty): 831-1.0098437-rdf.xml
- RDF/JSON: 831-1.0098437-rdf.json
- Turtle: 831-1.0098437-turtle.txt
- N-Triples: 831-1.0098437-rdf-ntriples.txt
- Original Record: 831-1.0098437-source.json
- Full Text
- 831-1.0098437-fulltext.txt
- Citation
- 831-1.0098437.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
https://iiif.library.ubc.ca/presentation/dsp.831.1-0098437/manifest