Scaleup of Cu(InGa)Se2 Thin Film Coevaporative Physical Vapor

Ind. Eng. Chem. Res. 2009, 48, 5975–5991

5975

Scaleup of Cu(InGa)Se2 Thin Film Coevaporative Physical Vapor Deposition Process, 1. Evaporation Source Model Development Kapil Mukati and Babatunde A. Ogunnaike* Department of Chemical Engineering, UniVersity of Delaware, Newark, Delaware 19716

Erten Eser, Shannon Fields, and Robert W. Birkmire Institute of Energy ConVersion, UniVersity of Delaware, Newark, Delaware 19716

Even though rapid advances have been made in improving Cu(InGa)Se2 thin-film-based solar cell efficiencies, the breakthroughs have been limited to the laboratory scale. Most commercially viable thin-film technologies reside at the premanufacturing development stage and scaleup has proven to be much more difficult than expected. Elemental in-line evaporation on flexible substrates in a roll-to-roll configuration is a commercially attractive process for the manufacture of large-area CuInSe2-based photovoltaics. At the University of Delaware’s Institute of Energy Conversion (IEC), such a process is being investigated, at the pilot scale, for a polyimide web substrate. The process works well for 6-in.-wide substrates and for short deposition runtimes. However, a commercially viable process is required to produce large-area, high-quality films at a muchhigher throughput. Specifically, the desired film thickness (∼2 µm) and composition uniformity must be achieved continuously and reproducibly on large-area (12-in.-wide) substrates at translation speeds (∼1 ft/ min) much higher than those currently used in pilot-scale processes. Although achievement of the desired film thickness and composition setpoints is best addressed by proper control system design, the film thickness uniformity is determined by the source design and individual nozzle effusion rates. The nozzle effusion rates, in turn, are dependent on the melt surface temperature profile and, thus, on the source design itself. Therefore, proper source design is critical in achieving film thickness uniformity. We have identified two modeling requirements for effective thermal evaporation source design: (i) a detailed three-dimensional thermal model of the evaporation source, to predict the melt surface temperature accurately, and (ii) a nozzle effusion model, to predict the effusion rate and the vapor flux distribution for a given nozzle geometry (length and diameter), melt surface temperature, and evaporant. To meet the first requirement, we have developed a first-principles three-dimensional electrothermal model using COMSOL Multiphysics software; the Direct Simulation Monte Carlo (DSMC) method has been used for effusion modeling. In this paper, which is the first of a two-part series, we present the details of the two models and their experimental validation. 1. Introduction Elemental in-line evaporation on flexible substrates in a rollto-roll configuration is a commercially attractive process for the manufacture of large-area CuInSe2-based photovoltaics. Cu(InGa)Se2-based photovoltaic technology has all the advantages of thin film processes: lower cost, monolithic module integration, and high-throughput deposition over large areas, when compared to silicon-based photovoltaics. More importantly, Cu(InGa)Se2 solar cells have demonstrated very high efficiencies at both the cellular level (∼19.5%) and the module level (∼13%);1 furthermore, they have shown excellent long-term stability in operation,2 as well as superior robustness to radiation3,4 (compared to Si wafer-based devices), and they can provide flexibility as well as very high power-to-weight ratios if deposited on lightweight substrates such as metal foils and plastics. Although rapid advances have been made in improving Cu(InGa)Se2 thin-film-based solar cell efficiencies, the breakthroughs have been limited to the laboratory scale. Most potentially commercial thin-film technologies reside at the premanufacturing development stage and scaleup has proven to be much more difficult than expected.5 At the University of Delaware’s Institute of Energy Conversion (IEC), copper-indium-gallium diselenide thin films are * To whom correspondence should be addressed. E-mail address: [email protected].

deposited onto a flexible substrate (or “web”), using a modified three-stage coevaporative physical vapor deposition (PVD) technique, in a roll-to-roll pilot-scale process (shown schematically in Figure 1). The film is deposited by thermal evaporation from a series of elemental sources (with two nozzles each) located sequentially in a vacuum deposition chamber. More specifically, copper, gallium, and indium are delivered to the substrate in amounts determined by their effusion rates, while selenium vapor, which is provided in excess, is incorporated into the film as determined by the stoichiometry. The desired elemental effusion rates are achieved by controlling the individual source temperatures, which, in turn, are controlled by manipulating the power supplied to the resistive heaters. The desired web temperature is achieved with two large resistively

Figure 1. Schematic of the pilot-scale inline coevaporative Cu(InGa)Se2 thin-film physical vapor deposition (PVD) process.

10.1021/ie8015957 CCC: $40.75  2009 American Chemical Society Published on Web 05/11/2009

5976

Ind. Eng. Chem. Res., Vol. 48, No. 13, 2009

heated and individually controlled steel plates placed near the top of the web. 1.1. Scale-up Objectives. The IEC’s PVD process works reasonably well at the pilot scale for substrates that are up to 6 in. wide, providing acceptable film thickness uniformity and composition uniformity across the substrate. The quality of the Cu(InGa)Se2 films has been demonstrated to give an average efficiency of 9.2% for 0.56 cm2 cells for a population of 44 cells distributed over a 60-in.-long and 6-in.-wide polyimide web.6 However, the substrate translation speed is 50%. In addition, the process should be designed such that the source operating temperatures are as low as possible for ease of operation, energy efficiency, and to maintain material integrity. The PVD system at the IEC is designed to accommodate substrates up to 12 in. wide. For this reason, 12in.-wide substrates are considered in the present study. Nevertheless, the design methodology presented here should be scalable to wider substrates. Achieving scale-up objectives that have been mentioned previously require knowledge of the source temperature profile, as well as accurate estimation of nozzle effusion rates and vapor flux distribution; therefore, the development of accurate evaporation source thermal and effusion models is essential. In this paper, we present a finite-element electrothermal model to predict source thermal profile, and a Direct Simulation Monte Carlo (DSMC) model to accurately predict effusion rate and vapor flux distribution. The paper is organized as follows: In Section 2, effusion and deposition models are presented for the inline PVD process. In Section 3, two modeling requirements are identified that are essential for proper design of evaporation sources. Sections 4 and 5 present, in detail, the development and validation of the


electrothermal model and DSMC model, respectively. The paper is summarized and concluded in Section 6. 2. Effusion and Deposition Subsystem Models The physics of the PVD process suggests a natural division into two distinct subsystems: (i) the effusion subsystem and (ii) the deposition subsystem. The effusion subsystem consists of the elemental sources for copper, gallium, and indium (selenium is provided in excess and is assumed to react stoichiometrically). Electric power supplied to the source resistive heaters is the input, and the measured outputs are the source temperatures at a specified location (usually the source bottom). The critical but unmeasurable process variables are the melt surface temperature profile, elemental nozzle effusion rates, and vapor flux distributions. The deposition subsystem consists of a moving substrate onto which elemental fluxes from the effusion sources are deposited. The web speed coupled with the incident elemental fluxes determine the film composition and thickness. The work on the PVD process modeling has been ongoing at the IEC for several years; this study builds on the work of Junker,9 Hanket,10 and Rocheleau,11 which is summarized in the following subsections. 2.1. Effusion Subsystem Model. Effusion subsystem model can be divided into two parts: (i) a source thermal model, which involves determination of the melt surface temperature given the electrical power input supplied to the resistive heaters; and (ii) an effusion model, which determines source effusion rate and vapor flux distribution based on the melt surface temperature obtained from the source thermal model. The source thermal model is presented in detail in Section 4. 2.1.1. Effusion Model. The effusion phenomenon in a thermal evaporation source with a rate-restricting nozzle can be described as a combination of two processes:10 evaporation from the melt surface and vapor effusion through the short nozzle. The net mass rate of evaporation (m ˙ evap) from the melt surface is proportional to the difference in the saturation pressure Psat and the gas pressure p above the evaporant surface. It is given by the Hertz-Knudsen equation:12 Psat - p m ˙ evap )η Asurf √2πRT/M

(5)

where Asurf is the melt surface area, T the melt surface temperature (in Kelvin), M the molecular weight, R the ideal gas constant, and η the evaporation coefficient, which is defined as the ratio of the real evaporation rate to the maximum theoretical value.13 The η value is dependent on the extent of evaporant surface contamination (low values for dirty surfaces and unity for clean surfaces). In the current work, η is assumed to be 1, because very-high-purity materials are used. Vapor effusion through the rate-restricting nozzle is dependent on the vapor flow regime and the nozzle geometry. Three flow regimes can be defined, based on the value of the Knudsen number (Kn): Kn )

λ 2Γ

(6)

where Γ is the nozzle radius and λ is the mean free path obtained from the kinetic theory of gases.14 λ)

kbT

√2πσ2p

(7)

Here, kb is the Boltzmann constant and σ is the collision diameter. When Kn > 1, the flow is free molecular, i.e.,

5977

molecule-wall collisions determine the flow behavior; when Kn < 0.01, intermolecular collisions dominate and the flow is said to be Viscous in nature; and when 0.01 e Kn e 1 the flow is transitional, where both molecule-wall and intermolecular collisions define the flow characteristic. In the IEC’s PVD process, we mainly encounter the transitional flow regime. Effusion flow rate in the transitional flow regime can be described by solving the following two equations:11,15 m ˙ ent ) BAnoz m ˙ noz ) Anoz

( )(

p - pent

(8)

√2πRT/M

)[

pent2 - pvac2 Γ2 2 λ 1+4 -1 16µL RT/M f Γ

(

)]

(9)

where m ˙ ent is the rate at which mass is entering the nozzle and pent is the pressure at the nozzle entrance. The entrance losses are taken into consideration by a correction factor B. The parameter f represents the fraction of molecules diffusively reflected from the nozzle walls, and µ is the viscosity, which is estimated using Chapman-Enskog theory.16 Using the mass continuity relationship, m ˙ evap ) 2m ˙ ent ) 2m ˙ noz

(10)

˙ noz, pent, and eqs 5, 8, and 9, the four unknown variables m ˙ ent, m and p can be determined. The vacuum chamber pressure pvac is neglected. Note that the parameter B is dependent on the source geometry and f on the material characteristics. Therefore, these parameters must be fitted experimentally. In the past evaporation source models, a value of f ) 0.9 was assumed based on the data reported in the literature.17 Junker9 obtained a value of B ) 1.09 for the copper source experimentally; this same value was used for the gallium and indium sources for the lack of experimental data. 2.2. Deposition Subsystem Model. The deposition subsystem model determines the evolution of the film thickness and composition on the moving substrate given the elemental evaporation rates and their spatial flux distribution. The local spatial flux distribution on the substrate is dependent on the vapor flux distribution emanating from the source nozzles. (The vapor flux distribution is also referred to as the molecular beam distribution or vapor flux profile.) 2.2.1. Vapor Flux Distribution. The flux profile of the vapor effusing from a cylindrical nozzle, with diameter 2Γ and length L, of an isothermal enclosure is described by the modified cosine law of emission:18 dm ˙ )m ˙ evap

( n +2 1 ) cos (θ)πrcos(φ) dA n

2

(11)

where dm ˙ d is the rate of deposition on an area element dAe that lies at a distance r from the effusion orifice (see Figure 2) and n is the beam collimation parameter. Higher values of n represent a greater degree of beam collimation. For free molecular flows, n has been determined to be proportional to the nozzle aspect ratio L/(2Γ). A value of n ) 3 was determined experimentally for the IEC’s copper source;10 the same value also was used for the other sources. If the substrate is oriented at an angle φ relative to the area element dAe, then the flux rate on the area of interest dA is given by cos(θ) cos(φ) dA (12) πr2 2.2.2. Deposition Model. The deposition model predicts the film growth by estimating the local elemental flux at discrete surface elements on the moving substrate.19 The model is dm ˙ )m ˙ evap

5978


where an empty web strip enters the deposition zone and is sequentially coated with the elemental mass as the strip moves through the deposition zone. The mass gain at a particular area element ∆A in time ∆t that is due to a single nozzle is given by ∆mS )

( n 2π+ 1 )( cos r (θ) )m˙ n+1 2

noz∆A∆t

(17)

where S ) Cu, Ga, In, or Se. Therefore, the combined mass gain due to all the sources at a particular substrate position (i,j), where i ) (1, 2, ..., Nx) and j ) (1, 2, ..., Ny), is obtained as

∑

∆mi,j )

∆mSi,j

(18)

S)(Cu,Ga,In,Se)

Figure 2. Knudsen cell effusion geometry.

where the selenium mass gain is determined stoichiometrically as NSe i,j )

Ga In NCu i,j + 3(Ni,j + Ni,j) 2

Se Se ∆mi,j ) NSe i,j Mw

(19) (20)

S is the molar gain due Here, Mw is the molecular weight and Ni,j to a particular element S (where S ) Cu, Ga, or In), which is given by

S ∆mi,j

S ) Ni,j

Figure 3. Geometry for the deposition model.

simplified by neglecting the reaction kinetics and material diffusion in the film. The sticking coefficient of the constituent elements is assumed to be unity and the excess selenium is assumed to be incorporated in the film stoichiometrically. A well-defined deposition zone of length of Ld and width W is obtained by placing steel plates between the sources and the web. This is the substrate area onto which the vapor strikes as the web moves over the sources. A global Cartesian coordinate system (x,y,z) is defined for the deposition system as shown in Figure 3. The nozzle-to-substrate distance is H, and the web speed is V in the positive y-direction. The substrate is discretized into Nx and Ny finite elements across the width and along the length, respectively yielding small-area elements of dimension ∆x ) W/Nx and ∆y ) Ld/Ny. Because vapor flux distribution is described for individual nozzles, another Cartesian coordinate system (ξ,ζ,z) is defined locally for every nozzle (xnoz,ynoz) such that ξ ) x - xnoz

(13)

ζ ) y - ynoz

(14)

Noting that the substrate plane is parallel to the nozzle plane (i.e., φ ) θ), for the specific nozzle depicted in Figure 3, the local flux rate at the area element ∆A ) ∆x∆y is therefore given by

)(

)

∆m ˙ n + 1 cos (θ) ) m ˙ evap ∆A 2π r2

(

n+1

(15)

j 2 + ξ2 + ζ2]1/2 and cos(θ) ) H/r. where r ) [H Junker19 simulated the web motion, using the cinematic approach with discrete time steps: ∆t )

∆y V

(16)

(21)

MwS

These calculations were repeated for all (i,j) at the next time instant k, by shifting the web by ∆y and, thus, obtaining the following relationship: Mk ) shift(Mk-1) + ∆mk

(22)

where M is the mass matrix of dimension (Nx,Ny). While computationally expensive, this formulation was necessary to obtain a state-space formulation of the deposition model that could be used for control simulations. However, the work in this thesis involves process design, and it is assumed that the source temperatures and the web speed are constant. Consequently, the mass gain ∆mi,j at a given position (i,j) is the same for all times. Therefore, the algorithm is simplified by calculating ∆mi,j only for one time-step ∆t, and then summing the film thicknesses of the individual finite elements in the y-direction to obtain the film thickness at the deposition zone exit (i.e., j ) Ny): Ny

ti,exit )

∑ ∆m

i,j

j)1

FCIGS∆A

(23)

where FCIGS is the film density. The film thickness nonuniformity is given by tnu )

tmax,exit - tmin,exit 2tmax,exit

(24)

The number of moles of a particular metal at the exit is obtained in a similar manner, as Ny

NSi,exit

)

∑N

S i,j

(25)

j)1

Therefore, the final film composition, as described using the copper ratio and the gallium ratio, is


Rcgi i,exit )

Rggi i,exit )

NCu i,exit NGa i,exit

+

NIn i,exit

NGa i,exit In NGa i,exit + Ni,exit

(26)

(27)

3. Model Requirements The most important variable of interest in the PVD process is the melt surface temperature under the nozzles. The effusion rates and, hence, the final film thickness and composition values are determined by this variable. Direct and reproducible measurement of the melt surface temperature is very difficult, if not impossible. In the IEC’s PVD process, the melt surface temperature is estimated indirectly from the source temperature. Such an arrangement suffices for the process control, especially if an outer-level multivariable control loop is implemented, which provides the temperature setpoints for the inner “effusion loops’’ to achieve desired film thickness and composition. In the absence of process output feedback, trial-and-error experiments provide the desired source operating temperatures (a strategy used for the pilot-scale process). However, for process design purposes, accurate knowledge of the temperature profile of the evaporation source is essential. Therefore, the first modeling requirement is the development of a detailed threedimensional thermal model of the evaporation source for the following reasons: (1) to ascertain the effect of heater assembly on melt temperature profile; (2) to accurately predict melt surface temperature, which is required to estimate the nozzle effusion rate; (3) to determine the relationship between source temperature (the measured variable) and melt surface temperature (the variable of interest) (this will be helpful in determining the best location for the source temperature measurement); and, finally, (4) to evaluate the effect of melt depletion on the melt temperature profile. Junker20 developed a one-dimensional source thermal model that was numerically tractable and thus easily incorporated in a model predictive controller. This simple model is not useful for the purposes listed above, because it is based on the assumption that the melt temperature profile is uniformsthe very issue that is under investigation. In fact, such an assumption may lead to the erroneous conclusion that the scale-up issues are the result of an ineffective process control strategy, rather than a faulty design. A more-detailed source thermal model was developed by Doody,21 using Fluent, Inc.’s FIDAP software. However, the model was preliminary in nature, and with several strong assumptions, such as replacing the insulation assembly with an effective emissivity, neglecting the effect of power leads, and assuming constant material properties. In addition, the temperature dynamics were neither investigated nor validated with the experimental data. In the current paper, a three-dimensional electrothermal model is developed for the entire evaporation source assembly, using the COMSOL Multiphysics v3.2 software package. The details are presented in Section 4. Evaporation sources in a commercial process, where the web speeds are much higher, must operate at higher effusion rates to achieve the desired film thickness. This increase in effusion rate may be achieved by operating the source at higher temperatures and/or by increasing the nozzle diameter. However, the former may lead to operational difficulties such as material breakdown, polymer substrate degradation, undesirable chemical

5979

reactions, thermocouple failure, and high power consumption. Therefore, the required effusion rate should be achieved primarily by increasing the nozzle diameter. Thus, the second model requirement is the prediction of effusion rate and vapor flux distribution for a given nozzle geometry (length and diameter) and melt surface temperature. The theoretical work of Rocheleau11 presented in Section 2 is ideally suited for effusion modeling, because of its mathematical tractability. However, the model requires two parameterssthe entrance loss factor B and fractional diffusivity fswhich must be determined experimentally. Therefore, the model’s validity is restricted to the particular source design and evaporant on which the experiments were performed. Obtaining the optimum nozzle geometry that will produce the desired effusion rate will require numerous experiments and, hence, a considerable investment of time and money. To reduce the development time and cost, a modeling tool is needed that will provide accurate predictions of effusion rate and flux distribution for any given nozzle geometry and melt temperature. In the current work, this objective is accomplished using the Direct Simulation Monte Carlo (DSMC) method, which is discussed in detail in Section 5. 4. Source Thermal Model A detailed three-dimensional finite-element thermal model is developed for the evaporation source assembly using the commercially available COMSOL Multiphysics software package. This type of detailed electrothermal model is fairly complex and computationally expensive. Therefore, several simplifying assumptions were made, so that the model could be solved on a personal computer (Intel Pentium, 2.8 GHz, 1 GB RAM): (1) It is assumed that all the contacting surfaces are in good thermal contact with each other. Surface roughness, and thus the possibility of radiative heat transfer between the surfaces in contact, is neglected. (2) Viscous dissipation in the melt and free-convective heat transfer with the source walls are assumed to be negligible, because the melt volume is small and the melt velocity is very low. (3) The interior of the source boat is assumed to be rectangular, to ease the memory requirement and computational load of mesh generation. In reality, the bottom part of the source boat is curved to facilitate the removal of the metal charge upon solidification. (4) All surfaces are opaque and radiate diffusively. (This is a common assumption in radiative heat transfer of real surfaces.22) (5) The contribution of the vapor to the heat transfer is assumed to be negligible, because of its very low density. (6) The outer stainless-steel cage that provides support to the insulation is replaced with its effective emissivity. This assumption is valid because the steel thickness is very small (∼0.025), relative to the thickness of other structures. (7) The chamber is modeled as a blackbody at the uniform temperature of 300 °C. The following heat and electrodynamic equations constitute the model: Heat Conduction Equation: FCp

∂T + ∇·q ) Q ∂t

(28)

where T is the temperature, F the density, Cp the heat capacity at constant pressure, q the heat flux, and Q the heat source.

5980


Figure 5. Source assembly geometry used in the model.

Figure 4. Example geometry for surface-to-surface radiation.

Because the heat transfer is via conduction only, the heat flux is given by q ) -k∇T The general boundary condition equation is

(29)

-n·q ) q0 + ε(G - σT4)

(30)

where q0 is the constant heat flux and the second term on the right-hand side is the radiative heat flux, assuming ideal gray surfaces; G is the irradiation (the total arriving radiative flux on a surface). For surfaces radiating to an ambient temperature Tamb, the irradiation term G ) σTamb4. For surface-to-surface radiation, G is given by G ) Gm + FambσTamb4

(31)

where Gm is the mutual irradiation arriving from other surfaces in the modeled geometry and Famb is the ambient view factor (the portion of view covered by the ambient conditions). Gm and Famb are obtained by solving the following surface integrals (see Figure 4):22

∫

(-n′·r)(-n·r) J′ dS (32) S′ πr4 (-n′·r)(-n·r) Famb ) 1 - S′ dS (33) πr4 Here, J is the radiosity, the total outgoing radiative flux: Gm )

∫

J ) (1 - ε)G + εσT4 Electrodynamic Equation:

(34)

-∇·(σe∇V - Jc) ) Qj

(35)

Here, V is the electric potential, σe the electrical conductivity, Jc the externally generated current density, and Qj the current source. The boundary condition at one end of the current lead is V ) V0, where V0 is the constant voltage applied to the heater, and the other current lead is grounded (i.e., V ) 0). All the other surfaces are electrically insulated (n · J ) 0). The resistive

heating term (in units of W/m3) is given by Q ) σe|∇V|2, which is used in eq 28 as the heat source.23 COMSOL’s heat transfer module is used to solve the heat equation, and the conductive-media DC module is used to solve the electrodynamic equation. The input to the model is the voltage applied across the molybdenum current leads. The geometry of the entire source assembly used in the model is shown in Figure 5. While the thermal model developed here is applicable to all three effusion sources (copper, gallium, and indium), the simulations and experimental validation are performed only for the copper source. 4.1. Material Properties. The materials that constitute the source assembly are given as follows: boron nitride (BN), for the source boat; graphite, for the heater element and the boat support; molybdenum (Mo), for the power leads; alumina, as the insulation; stainless-steel (SS), for the insulation support; and copper, gallium, or indium, which is used as the evaporant. Because the evaporation source temperature spans a wide rangesfrom room temperature (27 °C) up to 1500 °Csthe temperature dependence of the material properties is taken into account whenever the data are available and is summarized in the Appendix. The electrical conductivity of graphite (σGe ) and the outer insulation (stainless-steel) emissivity (εSS) are estimated experimentally, because of the nonavailability of reliable data in the literature. The following linear σGe vs T relationship is obtained by fitting the experimentally obtained steady-state source temperature to that predicted by the thermal model: σGe ) σmax -

σmax - σmin (T - T1) T1 - T2

(36)

where T1 ) 293 K, T2 ) 1773 K, σmax ) 105(Ω m)-1, and σmin ) 7.1 × 104(Ω m)-1. Figure 6 shows the σGe vs T plot. Values of the current that have been determined experimentally at different source temperatures are also shown for comparison, which indicate that the variation in resistance (and, hence, electrical conductivity) is linear with temperature for a constant voltage input. 4.2. Solver Parameters and Solution. The COMSOL Multiphysics software uses the finite-element method to solve the model equations. The parameters used to generate the threedimensional (3-D) tetrahedral mesh elements are chosen such that an acceptable resolution of the temperature profile is achieved using a minimum number of nodes. The mesh used in the model is shown in Figure 7. The UMFPACK solver is used to solve the finite-element problem. The solution is in the form of a 3-D temperature profile of the evaporation source,


Figure 6. Electrical conductivity of graphite.

Figure 7. Mesh elements in the source geometry.

which includes the copper melt. Figure 8 shows one such solution when a voltage of 11 V is applied across the molybdenum current leads; the resulting current is 65 A and power dissipation is 719 W. 4.3. Experimental Validation. An experimental validation study is conducted for both temperature dynamics and melt temperature profile predictions. The experiments are performed on the copper source, which has been placed in a separate vacuum bell-jar experimental setup to avoid disrupting the normal functioning of the pilot-scale system. The copper source is selected because it offers the worst-case design and control scenario, because of its very high operating temperature, compared to other elements in the system (gallium, indium, and selenium) and questionable melt temperature uniformity. 4.3.1. Temperature Dynamics. Experiments were performed to obtain the copper source temperature response when the following step changes are made in the voltage applied across the power leads: (i) 0 V to 11 V and (ii) 0 V to 12 V. The experimental data are then compared with the predictions from the thermal model. The plots are shown in Figures 9 and 10, respectively, which demonstrate that the thermal model is able to predict the source temperature dynamics with reasonable accuracy. The thermocouple measurements become very noisy near the copper melting point, because of some unresolved hardware issues in the National Instruments’s thermocouple module. (This difficulty is most probably due the interaction between the thermocouple module and the power supply.) In the figures, for clarity, such noisy measurements are not plotted.

5981

4.3.2. Melt Temperature Gradient. The two-dimensional (2-D) melt surface temperature contour predicted by the thermal model for a lead voltage of 12 V is shown in Figure 11a. Figure 11b shows the temperature variation along the centerline of the melt surface. The source bottom-center temperature is ∼1400 °C, and the melt temperature difference between the far-side and lead-side nozzles is determined to be ∆T ≈ 25° C. Because experimental data for the melt surface temperature cannot be obtained, the model predictions are validated in an indirect manner, using the normalized film thickness profile. The ∆T ) 25° C value given by the thermal model is used in the effusion-deposition model to obtain the copper film thickness profile, which is then compared in Figure 12 to the experimental film thickness data. The comparison shows that the model predictions match very well with the experimental data. The thermal model, therefore, predicts the melt surface temperature profile fairly accurately. 4.4. Symmetric Heater Design Validation. Experiments have shown that the effusion rates from the two source nozzles are not equal. This mismatch has been attributed to the use of an asymmetric heater (i.e., a “U”-shaped heater element with both power leads on one side). We propose that a symmetric heater assembly should be used to resolve the issue of unequal nozzle flow rates. This design modification is validated using both the thermal model as well as experiments. To perform the experiments, a symmetric graphite heater was constructed, and the boron nitride (BN) support and end insulation were modified to accommodate a molybdenum power lead on each side of the source; the electrical feedthroughs were rearranged accordingly. The melt surface temperature profile predicted by the thermal model is plotted in Figure 13, which shows that the melt temperatures under the two nozzles are equal. Therefore, the two nozzle effusion rates should also be equal. This is confirmed by the experimental film thickness results given in Figure 14, because the film thickness peaks above the nozzles are of equal height. 4.5. Best Location for Control Thermocouple Placement. For effective effusion rate control, accurate measurement of the melt surface temperature below the nozzle (Tsurf) is required. Because direct measurement of the melt surface temperature is difficult (if not impossible), the temperature of a solid source wall is measured instead, using a thermocouple. Therefore, the thermocouple location should be where the source temperature varies closely with Tsurf. Three possible locations are shown in Figure 15, where 2 represents Tsurf: 1, which is the location at the source lid just above 2, location 3 is at the source bottom directly below 2, and location 4 is located at the source-bottom center. (Location 4 is currently being used in the inline pilot-scale system for source temperature control.) To determine the best placement for the control thermocouple, the pulse temperature response of the source at the above four locations is determined using the thermal model under the following conditions: from a steady-state melt surface temperature of Tsurf ≈ 1400 °C and a corresponding heater-voltage value of 7.5 V, the applied voltage is increased to 8 V for a period of 30 min, and then reduced back to the initial value of 7.5 V. The resulting temperature dynamics are shown in Figure 16. From the figure, we can clearly observe that the temperature of location 3 most closely tracks the melt surface temperature Tsurf, and, therefore, should be the thermocouple location. 4.6. Effect of Melt Depletion on the Melt Surface Temperature. Simulations are conducted on a source with a symmetric heater for the following two melt levels to evaluate the effect of melt depletion on the melt surface temperature: (i)

5982


Figure 8. Copper source thermal modeling results.

Figure 9. Copper source temperature dynamics for 11 V voltage input: the temperatures of the source bottom below the far-side and lead-side nozzles are shown.

an almost-full source, where the melt level is 0.025 in. below the nozzle entrance, and (ii) an almost-empty source, where this distance is 0.6. The maximum possible melt height is 0.625. A voltage of 7.5 V is applied across the heater element to obtain melt temperature near 1400 °C. The resulting melt surface temperature profiles are shown in Figure 17, which indicate that the effect of melt depletion on the melt surface temperature above the nozzles is negligible ((3 °C) for the pilot-scale twonozzle sources. However, the change in melt temperature at the center is ∼17 °C. Such a variation in temperature may

significantly affect the nozzle effusion rate in a commercialscale source with more than two nozzles (e.g., a three-nozzle source with a nozzle placed at the center). Therefore, any source design study should take melt temperature variation with melt depletion into account. 5. DSMC Effusion Model We use the Direct Simulation Monte Carlo (DSMC) method to predict the effusion rate and vapor flux distribution for a given


5983

Figure 10. Copper source temperature dynamics for 12 V voltage input: the temperatures of the source bottom below the far-side and lead-side nozzles are shown.

Figure 13. Melt surface temperature profile for symmetric heater predicted using a thermal model. Figure 11. Thermal modeling: melt surface temperature results.

Figure 12. Experimental and model-predicted film thickness profiles.

nozzle geometry and melt surface temperature. The DSMC method, which was first proposed by Bird,24 is a computational algorithm for the stochastic simulation of rarefied gas flows. The method is particularly efficient for the transition regime (0.1 < Kn < 10), which occurs inside the evaporation source at the typical operating temperatures. 5.1. DSMC Fundamentals. The fundamental idea behind the DSMC method is to statistically represent the large number of actual atoms or molecules (Nreal) by a few simulation particles (Nsim). The ratio FN ) Nreal/Nsim may be as large as 1020. Such a considerable reduction in the number of simulation particles makes this algorithm remarkably efficient, with regard to the

Figure 14. Experimental normalized film thickness profile for a symmetric heater.

computational time and memory requirements. The particles’ motion and interactions (interparticle collision or wall collision) are used to modify their positions and velocities while enforcing the mass, momentum and energy conservation laws.25 In the DSMC method, therefore, the flow develops from a set of prescribed initial conditions and subsequently evolves by

5984


Over this time step, the particle motions are calculated deterministically, while the collisions are treated statistically, based on the gas kinetic theory. The physical flow space is divided into virtual cells with dimensions lc not greater than the mean free path (λm). This cell network facilitates the selection of collision pairs and the sampling of the macroscopic flow properties. One may take advantage of the symmetries in physical space to reduce the dimension of the cell network and also the number of cells needed for the simulation. Irrespective of the dimensionality, collisions are always treated as 3-D phenomena, and only the particles in the same spatial cell are allowed to collide. The DSMC algorithm may be divided into five basic steps: initialization, particle motion, particle indexing, interparticle collision, and flow sampling. 5.1.1. Initialization. In this first step, the cell length lc < λm, and the time step ∆t < τm are chosen. The typical values are lc ) λm/3 and ∆t ) τm/5. FN is chosen such that at least 20 to 30 simulation particles are present per λm3: this is to reduce the statistical scatter when the average macroscopic properties are sampled for each cell. The particles’ position and velocity are prescribed based on the entrance boundary conditions and the specified velocity distribution. 5.1.2. Particle Motion. In the second step, all the particles (i) are moved deterministically through distances (r) according to their velocity components (v) and the discrete time step (∆t):

Figure 15. Possible locations for control thermocouple placement.

ri(t + 1) ) ri(t) + vi(t)∆t Figure 16. Source temperature response at various locations for a pulse power input.

(39)

Appropriate action is taken if the particle intersects the boundaries representing solid surfaces, symmetry, or the flow exit. Collisions with surfaces can be treated as being either fully specular, fully diffuse, or some combination of the two. Specular collisions are modeled as a simple reversal of the particle velocity component normal to the incident surface. The other velocity components remain the same. Diffuse collisions cause a random reorientation of the reflected molecule, where the postcollision velocity is related to the temperature of the surface. The distribution function for a thermal velocity component u in an equilibrium gas is given by:24 fu )

( )

1 u2 exp - 2 Vp Vp√π

(40)

where Vp is the most probable velocity, given as Vp ) Figure 17. Melt surface temperature profiles for an almost-full source and an almost-empty source.

tracking the particle positions and velocities while accounting for the collisions and boundary interactions. The particles constantly enter and leave the computational domain, as determined by the specified flow conditions at the boundaries. The essential approximation of DSMC is that the particle motion and interparticle collisions can be decoupled over small time intervals, provided the simulation time step ∆t is less than the mean collision time τm: τm )

λm Va

(37)

where Va is the average velocity of a monatomic gas at temperature T with molecular weight M: Va )

8RT πM

(38)

2RT M

(41)

If the incident particle collides with the surface normal to the x-axis, then the velocity components may be sampled as follows from fu, using the procedure in Appendix C of ref 24:

[]

[

Vp√-ln(U) Vx Vy ) Vp√-ln(U) cos(φ) Vz Vp√-ln(U) sin(φ)

]

(42)

where U is a uniform random number between 0 and 1, and φ is sampled uniformly from [0, 2π]. More-complicated surface interaction models are available in the literature,26,27 where the orientation of the reflected molecule is also a function of the incident velocity. However, the diffuse reflection model is adequate for real surfaces.28 The number of particles that have crossed the exit boundary (and, hence, left the computational space) are determined. In the case of an effusion source, this boundary is represented by the nozzle opening. For each particle leaving the nozzle, a new


particle is generated randomly at the evaporant surface, so that the total number of simulation particles remains constant. 5.1.3. Particle Indexing. In the third step, the particles are indexed according to the cell in which they reside at the current time instant. The number of particles Nc present in each cell is counted and stored. 5.1.4. Interparticle Collision. The fourth step, interparticle collision simulation, is a probabilistic process that sets DSMC apart from deterministic simulation methods such as molecular dynamics (MD).29 The probability of a collision between two molecules in a homogeneous gas is proportional to the product of their relative speed Vrel and total collision cross-section σT. For each cell, the mean value of the product σTVrel is calculated, and the maximum value (σTVrel)max recorded. We use Bird’s24 “no-time-counter” (NTC) technique to obtain the number of collision pairs. In the NTC method, j cFN(σTVrel)max∆t NcN 2Vc

Mc )

(43)

potential collision pairs are selected from the cell, where Nc is j c the timethe instantaneous number of particles in the cell, N averaged value for Nc, and Vc the cell volume. The collision is computed with probability σTVrel/(σTVrel)max, where the parameter (σTVrel)max may be set initially at a reasonably high value (for example, 3σTVp) and is updated if a larger value is encountered during the sampling. The new direction of the components of the post-collision relative velocity is determined randomly as follows:

[][ x Vrel

Vrel cos(θ) y Vrel ) Vrel sin(θ) cos(φ) Vrel sin(θ) sin(φ) Vz rel

]

(44)

Here, U is a uniform random variable from [0,1], θ ) 1 - 2U, and φ ) 2πU. 5.1.5. Flow Sampling. The final step is sampling the macroscopic flow properties. The velocity components and the number of particles in a particular cell are used to calculate macroscopic quantities at the geometric center of the cell. The macroscopic flow velocity (V), equilibrium temperature (T), and pressure (P) are given as follows: V)

T)

1 NcR

P)

1 Nc

Nc

∑v

(45)

i

i)1

Nc

∑V

2 i

)

i)1

NcFNkbT Vc

V2 R

(46) (47)

where kb is the Boltzmann constant. The nozzle effusion rate and the vapor flux distribution, which are the properties of interest in the current work, are determined as follows: (1) Effusion rate: The effusion rate is directly proportional to the number of particles exiting the nozzle per unit time. It is calculated as a time-averaged value: fe )

NoutFNM NAk∆t

(48)

where Nout is the total number of particles that have exited the nozzle in k simulation steps.

5985

Figure 18. Physical space for DSMC simulation of flow inside the evaporation source.

(2) Vapor flux distribution: Vapor flux distribution is determined from the frequency distribution of the angles at which the particles are leaving the nozzle plane. 5.2. DSMC Application to the Evaporation Source. To implement the DSMC method on the evaporation source, a physical computational space must be defined first. Using the symmetry of the nozzle, we use the volume shown in Figure 18 as the physical space. Therefore, DSMC simulation is performed only for one-quarter of the nozzle volume. Surface 1 in the figure is the melt surface and the entrance boundary for the particles. Surfaces 2, 3, 6, and 7 in the figure represent the real solid walls of the source; particle-wall interactions are assumed to be diffuse. Surfaces 4, 5, 8, and 9 in the figure represent the symmetry boundaries, whereas Surface 10 in the figure is the nozzle exit. The following assumptions are made to reduce the computational effort: (1) It is assumed that the effusive flow properties are dependent mainly on the conditions directly underneath the nozzle. Therefore, although the physical space is not strictly symmetric in the x-direction, the assumption will not affect the flow properties significantly. (2) All the surfaces of the source are assumed to be at a uniform temperature that is equal to the melt surface temperature. (3) The particles are assumed to behave similar to elastic hard spheres. (4) It is assumed that no interparticle collisions occur after the particles exit the nozzle (Surface 10 in Figure 18). This assumption is valid because the flow regime is molecular flow, because of very low chamber pressure. The algorithm is initialized as follows: cubic cells of dimension lc ) λm/3 are used to grid the physical space; therefore, the cell volume is Vc ) lc3. The simulation time step is ∆t ) τm/10. FN is chosen such that there are 20 simulation particles per cell. The particles are generated at random locations on the melt surface (Surface 1 in Figure 18) with the Maxwellian velocity distribution of eq 42. Now, for any given evaporant, source geometry, and melt surface temperature, the DSMC algorithm can be used to predict the effusion rate and the vapor flux distribution. As an example, Figures 19, 20, 21, and 22 show the DSMC simulation results for the copper source when the melt surface temperature is 1400 °C, and the geometry is that of the pilot-scale evaporation source (L ) 0.635, D ) 0.375, Lb ) 1.2, and Wb ) 0.375). The distance between the melt surface and the lid (Surface 6 in Figure 18) is assumed to be h ) 0.25. Figure 19 shows the evolution of effusion rate with the simulation time steps, ultimately reaching the steady-state value of 0.96 g/h (this value is obtained by averaging the values of the last 2000 time steps). The standard

5986


Figure 19. Evolution of copper effusion rate with DSMC simulation steps: T ) 1400 °C.

Figure 20. Copper vapor flux distribution predicted using DSMC method, T ) 1400 °C.

Figure 22. Copper vapor pressure drop along the nozzle centerline predicted using DSMC method: T ) 1400 °C.

is also plotted. It is observed that this function overpredicts the vapor flux for small angles and underpredicts for larger ones. A function of type

Figure 21. Copper vapor pressure contour lines along the source length and height predicted using the DSMC method: T ) 1400 °C.

deviation is determined to be Tm:45

k

(T e Tm)

-2

) 134.4 + 2.919 × 10 T

-5

(A14) -8 2

+ 7.193 × 10 T + 6.438 × 10 T (A15)

(A10) A6. Electrical Conductivity

(T > Tm) (A11)

The thermal conductivity of BN, kBN ) 23 W/(m K), is approximately constant in the temperature range of 1000-1650 °C.46 The values for the 20-1000 °C temperature range are obtained from the following relationship, which is fitted to the experimental data:38

(

kBN ) 2.2154 × 101 + 5.9937 × 101 exp -

)

T 2.7164 × 102 (A12)

For graphite, data from ref 42 for DFP-1 grade are used to fit a quadratic relationship: kG ) 120 - 9.2 × 10-2T + 2.8 × 10-5T2

(A13)

which is valid in the temperature range of 25-1600 °C.

Figure A3. Thermal conductivity plots.

-1

Plots for thermal conductivity versus temperature are shown in Figure A3.

kCu, solid ) 398.5 - 4.944 × 10-2T 9.653 × 10-5T2 Cu, molten

k ) 1.389 × 10 A

Electrical conductivity data for graphite and molybdenum are required to model the amounts of heat generated in the heater and the current leads, respectively. The electrical resistivity of molybdenum (σMo e ) is given by the following linear relationship obtained using the data in ref 47: -2 σMo e ) 3.9 + 2.9 × 10 T

(A16)

The temperature range for the data is 20-2500 °C. A7. Emissivity Normal emissivities are used for all the materials. For copper, the following linear relationship from Perry’s Chemical Engineering Handbook48 is applied: εCu(T) ) ε1 +

T - T1 (ε - ε1) T2 - T1 2

(A17)

5990


where T1 ) 1349.8 K, T2 ) 1549.8 K, ε1 ) 0.16, and ε2 ) 0.13. The emissivity of BN, εBN ) 0.81, is taken from the work of Walker and Casey.49 For graphite, the following expressions that have been derived by Junker20 are used: εG(T) ) af3(T) + bf2(T) + cf(T) + d

(A18)

where the temperature T is given in Kelvin and f(T) is defined in eq A6, with parameter values of a ) -1.3797 × 10-2, b ) 5.6707 × 10-4, c ) 1.1048 × 10-1, and d ) 7.1334 × 10-1. Additional parameters have the following values: µ1 ) 7.1835 and µ2 ) 0.9671. The emissivity of molybdenum is estimated from the following linear relationship obtained using the data in ref 50: εMo ) 6.5 × 10-2 + 1.3 × 10-4T

(A19)

Here, the temperature T is given in Celsius. For alumina, emissivity data from ref 51 are used to obtain the following relationship: εA ) 0.86 - 6 × 10-4T + 1.9 × 10-7T2

(for 20 °C < T < 1500 °C) (A20)

(for T > 1500 °C) εA ) 0.4 The emissivity plots are shown in Figure A4.

(A21)

Literature Cited (1) Green, M. A.; Emery, K.; King, D. L.; Hishikawa, Y.; Warta, W. Solar cell efficiency tables (version 28). Prog. PhotoVoltaics 2006, 14, 455–461. (2) Wieting, R. D. CIS product introduction: Progress and challenges. In AIP Conference Proceedings, National Center for Photovoltaics (NCPV) 15th Program Review Meeting, Munich, Germany, 1999; pp 3-8. (3) Jasenek, A.; Rau, U.; Weinert, K.; Ko¨tschau, I. M.; Hanna, G.; Voorwinden, G.; Powalla, M.; Schock, H. W.; Werner, J. H. Radiation resistance of Cu(In,Ga)Se2 solar cells under 1-MeV electron irradiation. Thin Solid Films 2001, 387, 228–230. (4) Tringe, J.; Nocerino, J.; Tallon, R.; Kemp, W.; Shafarman, W.; Marvin, D. Ionizing radiation effects in copper indium gallium diselenide thin-film solar cells. J. Appl. Phys. 2002, 91, 516–518. (5) Birkmire, R. W. Compound polycrystalline solar cells: Recent progress and Y2K perspective. Sol. Energy Mater. Sol. Cells 2001, 65, 17–28. (6) Hanket, G. M.; Singh, U. P.; Eser, E.; Shafarman, W. N.; Birkmire, R. W. Pilot-Scale Manufacture of Cu(InGa)Se2 Films on a Flexible Polymer Substrate. In Proceedings of the Twenty Ninth IEEE PhotoVoltaic Specialists Conference, New Orleans, LA, 2002; p 567. (7) Birkmire, R. W.; Eser, E. Polycrystalline thin film solar cells: present status and future potential. Annu. ReV. Mater. Sci. 1997, 27, 625– 653. (8) Kemell, M.; Ritala, M.; Leskela¨, M. Thin film deposition methods for CuInSe2 solar cells. Crit. ReV. Solid State Mater. Sci. 2005, 30, 1– 31. (9) Junker, S. T.; Birkmire, R. W.; Doyle, F. J., III. Mass and heat transfer modeling of a physical vapor deposition effusion source. AIChE J. 2005, 51, 878–894. (10) Hanket, G. M. Ph. D. Thesis, University of Delaware, Newark, DE, 1999. (11) Rocheleau, R. E. Ph. D. Thesis, University of Delaware, Newark, DE, 1981. (12) Maissel, L. I.; Glang, R. Handbook of Thin Film Technology; McGraw-Hill: New York, 1970. (13) Bunshah, R. F. Deposition Technologies for Films and Coatings: DeVelopment and Applications; Noyes Publications: Norwich, NY, 1982. (14) Kennard, E. H. Kinetic Theory of Gases; McGraw-Hill: New York, 1938. (15) Brown, G. P.; DiNardo, A.; Cheng, G. K.; Sherwood, T. K. The flow of gases in pipes at low pressures. J. Appl. Phys. 1946, 17, 802– 813. (16) Bird, R. B.; Stewart, W. E.; Lightfoot, E. N. Transport Phenomena; Wiley: New York, 1960.

(17) Rocheleau, R. E.; Baron, B. N.; Russell, T. W. F. Analysis of evaporation of cadmium sulfide for manufacture of solar cells. AIChE J. 1982, 28, 656–661. (18) Ohring, M. The Material Science of Thin Films; Academic Press: New York, 1992. (19) Junker, S. T.; Birkmire, R. W.; Doyle, F. J., III. Manufacture of Thin-Film Solar Cells: Modeling and Control of Cu(InGa)Se2 Physical Vapor Deposition onto a Moving Substrate. Ind. Eng. Chem. Res. 2004, 43 (2), 566–576. (20) Junker, S. T. Ph.D. Thesis, University of Delaware, Newark, DE, 2003. (21) Doody, J. Ph.D. Thesis, University of Delaware, Newark, DE, 2002. (22) Incropera, F. P.; DeWitt, D. P. Fundamentals of Heat and Mass Transfer, 4th Edition; Wiley: New York, 1981. (23) COMSOL Multiphysics 3 Manuals; COMSOL: Burlington, MA, 2006 (available via the Internet at www.comsol.com). (24) Bird, G. A. Molecular Gas Dynamics and the Direct Simulation of Gas Flows; Oxford University Press: Cambridge, U.K., 1994. (25) Oran, E. S.; Oh, C. K.; Cybyk, B. Z. Annu. ReV. Fluid Mech. 1998, 30, 403–441. (26) Cercignani, C.; Lampis, M. Kinetic models for gas-surface interactions. Transp. Theory Stat. Phys. 1971, 1, 101–114. (27) Lord, R. G. Some extensions to the Cercignani-Lampis gas scattering kernel. Phys. Fluids A 1991, 3, 706–710. (28) Prasanth, P. S.; Kakkassery, J. K. Direct Simulation Monte Carlo (DSMC): A numerical method for transition-regime flowssA review. J. Indian Inst. Sci. 2006, 86, 169–192. (29) Haile, J. M. Molecular Dynamics Simulation: Elementary Methods; Wiley: New York, 1992. (30) Birkmire, R.; Eser, E.; Fields, S.; Shafarman, W. Prog. PhotoVoltaics 2005, 13, 141–148. (31) Lide, D. R., Ed. CRC Handbook of Chemistry and Physics, 83rd Edition; CRC Press LLC: Boca Raton, FL, 2002; p 1-12. (32) Bondi, A. van der Waals Volumes and Radii. J. Phys. Chem. 1964, 68, 441. (33) Lide, D. R., Ed. CRC Handbook of Chemistry and Physics, 83rd Edition; CRC Press LLC: Boca Raton, FL, 2002; p 4-132. (34) Lide, D. R., Ed. CRC Handbook of Chemistry and Physics, 83rd Edition; CRC Press LLC: Boca Raton, FL, 2002; p 6-123. (35) Dushman, S. Scientific Foundations of Vacuum Technique, 2nd Edition; Wiley: New York, 1962. (36) Lide, D. R., Ed. CRC Handbook of Chemistry and Physics, 83rd Edition; CRC Press LLC: Boca Raton, FL, 2002; p 12-219. (37) Lide, D. R., Ed. CRC Handbook of Chemistry and Physics, 83rd Edition; CRC Press LLC: Boca Raton, FL, 2002; p 4-127. (38) Combat Boron Nitride Solid Grade AX05: Product Data Sheet, Saint-Gobain Ceramics, Amherst, NY (available via the Internet at www.bn.saint-gobain.com). (39) Sheppard, R. G., Mathes, D. M., Bray, D. J., Eds. Properties and Characteristics of Graphite for the Semiconductor Industry; Poco Graphite, Inc.: Decatur, TX, 2001. (40) Characteristics and Properties: Alumina Insulation Type ZAL45 and ZAL-45AA, ZIRCAR Ceramics, Inc. (available via the Internet at www.zircarceramics.com). (41) Eggers, D. F., Jr.; Gregory, N. W.; Halsey, G. D., Jr.; Rabinowitch, B. S. Physical Chemistry, 2nd Edition; Wiley: New York, 1964; p 752. (42) Touloukian, Y. S.; Buyco, E. H. Thermophysical Properties of Matter: Specific HeatsNonmetallic Solids; IFI/Plenum Data Corporation: New York, 1970; Vol. 5. (43) Brandes, E. A., Brook, G. B., Eds. Smithell’s Metals Reference Book, 7th Edition; Butterworth and Heinemann: Boston, MA, 1992; p 142. (44) Brandes, E. A., Brook, G. B., Eds. Smithell’s Metals Reference Book, 7th Edition; Butterworth and Heinemann: Boston, MA, 1992; p 144. (45) Brandes, E. A., Brook, G. B., Eds. Smithell’s Metals Reference Book, 7th Edition; Butterworth and Heinemann: Boston, MA, 1992; p 1411. (46) Touloukian, Y. S.; Powell, R. W.; Ho, C. Y.; Klemens, P. G. Thermophysical Properties of Matter: Thermal ConductivitysNonmetallic Solids; IFI/Plenum Data Corporation: New York, 1970; Vol. 2, p 656. (47) Brandes, E. A., Brook, G. B., Eds. Smithell’s Metals Reference Book, 7th Edition; Butterworth and Heinemann: Boston, MA, 1992; p 144. (48) Perry, H.; Green, W. Perry’s Chemical Engineering Handbook, 7th Edition; McGraw-Hill: New York, 1997; Table 5-6.

Ind. Eng. Chem. Res., Vol. 48, No. 13, 2009 (49) Walker, G. F. Casey, J. Measurement of total normal emittance of boron nitride from 1,200F to 1,900F with normal spectral emittance data at 1,400F, NASA Technical Report No. NASA TN D-1268, National Aeronautics and Space Administration, 1962. (50) Brandes, E. A., Brook, G. B., Eds. Smithell’s Metals Reference Book, 7th Edition; Butterworth and Heinemann: Boston, MA, 1992; p 176. (51) Glazman, M. S.; Landa, Y. A.; Litovskii, E.; Puchkelevich, N. A. Spectral emissivity of refractories. J. Refract. Ind. Ceram. 1989, 30, 312– 319.

5991

(52) Mukati, K.; Ogunnaike, B. A.; Eser, E.; Fields, S.; Birkmire, R. W. Scaleup of Cu(InGa)Se2 Thin Film Coevaporative Physical Vapor Deposition Process, 2. Evaporation Source Design. Ind. Eng. Chem. Res. 2009, 48, DOI: 10.1021/ie801596a.

ReceiVed for reView October 21, 2008 ReVised manuscript receiVed April 12, 2009 Accepted April 14, 2009 IE8015957

Scaleup of Cu(InGa)Se2 Thin Film Coevaporative Physical Vapor

Recommend Documents