APPENDIX 1
PHOENICSCVD
A customised CFD code for the
simulation of CVD processes
ABSTRACT
PHOENICSCVD is a specialpurpose code developed to simulate chemical vapour deposition (CVD) equipment, particularly in the semiconductor industry. The technical features of the code are described; these relate primarily to transport mechanisms, gasphase and surface chemistry, surfacetosurface radiation and plasma effects. Implementation is by means of input and data files, together with a userfriendly menu system in which much of the setup is automated . Extensive testing and validation of the code has been carried out.
1. INTRODUCTION
Chemical Vapour Deposition (CVD) is a technique that is widely used, particularly in the semiconductor industry. There is a continuous demand for scaleup of reactors to produce larger wafers; at the same time the feature sizes are reducing, bringing ever increasing demands on film quality and process performance. As a result, the life cycle of CVD equipment is short and new equipment has to be developed every few years. Conventionally, such requirements have been met by empirical trialanderror methods, straightforward scaleup of existing equipment and ad hoc modifications to existing processes. However, computer simulation offers the possibility of a more fundamental understanding of the essential physics and chemistry, leading to a systematic approach to reactor and process design.
PHOENICSCVD was developed to meet this challenge. It was created by a consortium of organisations within the ACCESSCVD Project, with financial support from the European Commission through its ESPRIT programme. The partners in the Project were CHAM (coordinator), Technical University of Delft, Fraunhofer Institut fuer Integrierte Schaltungen (IISB), Siemens Corporate Research and Development and ASM International.
The objectives of the ACCESSCVD Project were:
(i) the compilation of modelling equations for the transport, heat transfer and chemical phenomena important in CVD processes, together with the provision of a practical engineering tool to permit the inclusion of the increasingly important plasma effects,
(ii) incorporation of the models in a flow simulation code and their validation against experimental data,
(iii) the provision of an industryspecific userinterface to make the simulation package accessible to designers and operators of CVD equipment without the need for specialised CFD knowledge, and
(iv) the implementation of the package in an industrial environment by way of pilot applications, ensuring that it could meet the needs of daytoday use.
This document concentrates on the first of these objectives; the models that were developed for inclusion in the code are described in the following Section. For additional information the reader is referred to the special December 1995 issue of The PHOENICS Journal, which was exclusively devoted to the ACCESSCVD Project; further references to original source material can be found there.
2. MODEL AND SOFTWARE DEVELOPMENT
2.1 Transport Mechanisms
Various basic assumptions are made about the nature of the gas flow in the reactor. The first is that the gas mixture will behave as a continuum, i.e. that the mean free path for the molecules is much smaller than the characteristic dimensions of the reactor; this effectively imposes a lower limit on the operating pressure of, typically, about 30Pa. The assumption is made that ideal gas behaviour is a suitable approximation; the code has also been designed wi th laminar flow in mind, although this is not a requirement. Viscous heating due to dissipation is also neglected. All these limitations are entirely compatible with most CVD applications, although some very low pressure reactors should not be simulated using PHOENICSCVD.
The fundamental equations to be solved are those representing conservation of mass, momentum and energy, coupled with a similar equation for each of the individual gas species. The mass and momentum equations take the conventional CFD form; the energy and species equations have additional contributions to take account of chemical reactions and diffusion of species within the gas mixture.
2.1.1 Species transport equations
Species concentrations are expressed as dimensionless mass fractions, w _{i}; diffusion is expressed in terms of diffusive mass fluxes
(1)
where v_{i} is the i^{th} species velocity and v is the mass averaged mixture velocity. Diffusion arises from concentration gradients (ordinary diffusion) and also from temperature gradients (thermal or Soret diffusion); j_{i} therefore contains two separate contributions.
The conservation equation for the i^{th} species is then
(2)
where the terms on the right hand side represent, respectively, convection, diffusion and gas phase chemistry; m_{i} is molecular mass, is a stoichiometric coefficient and summation is over the reactions, with forward and backward reaction rates given, in molecular terms, by and .
2.1.2 Energy transport equation
The energy equation contains the standard convection and diffusion terms, augmented by the extra contributions referred to above:
(3)
where H_{i} is the species molar enthalpy and c_{p} is the mixture specific heat. The third term on the right hand side represents diffusion of energy associated with the interdiffusion of species and the final term the energy associated with the creation and destruction of species by the chemical reactions.
2.1.3 Ordinary diffusion
Ordinary diffusion fluxes in a binary mixture are given by Fick’s Law:
(4)
where D_{12} is the binary diffusion coefficient of the gas pair. The same approach can be applied to a multicomponent mixture in which a single carrier species dominates.
For nondilute mixtures this simple formulation can no longer be utilised; instead, the diffusive fluxes of each species is related to all the concentration gradients through the StefanMaxwell equations:
(5)
where M is the average molar mass of the mixture. For implementation in PHOENICSCVD it is more convenient to rearrange these equations to obtain expressions for each diffusion flux in terms of the others; this retains a conventional diffusion term in the species equations (useful for stability) with other terms being added as sources. A further advantage of the StefanMaxwell approach is that binary diffusion coefficients are used: these are simple functions of pressure and temperature; alternative methods require multicomponent diffusion coefficients, which are more complicated functions and correspondingly more expensive to evaluate.
As an alternative to the StefanMaxwell equations an approximate expression for multicomponent diffusion fluxes has been developed by Wilke; it retains the form of Fick’s Law,
(6)
but with an effective diffusion coefficient:
(7)
where x represents molar fraction. The Wilke approximation is accurate for dilute mixtures and offers computational savings when compared with StefanMaxwell.
The philosophy within PHOENICSCVD has been to offer a range of alternatives to the user where these can save significant computational time. Fick, Wilke and StefanMaxwell formulations are therefore all included in the code.
2.1.4 Thermal diffusion
Thermal diffusion causes small, light molecules to diffuse towards hot regions while larger, heavier molecules move the other way. This effect can be important in coldwall reactors which may have large temperature gradients. The mass fluxes are given by
(8)
where the thermal diffusion coefficient is a function of gas composition and temperature, but independent of pressure.
2.1.5 Thermodynamic properties of gas mixtures
For each gas species such properties as specific heat, heat of formation and standard entropy are needed in the energy equation given above (3); mixture specific heat is also required. These properties are functions of temperature and are usually defined by polynomial coefficients. PHOENICSCVD follows the conventions adopted in the Chemkin thermodynamics database: for each species two sets of seven coefficients are required, one for a low temperature range and one for a high range.
Mixture specific heat is simply defined by mass averaging.
2.1.6 Transport properties of gas mixtures
Solution of the conservation equations requires mixture properties such as density, viscosity and thermal conductivity to be known as functions of local conditions and composition.
Experimental data are available for only a few common gases and these properties are derived from kinetic theory for the individual gas species; LennardJones parameters are required and these are taken from a data file whose format is again compatible with the corresponding Chemkin database.
Density poses no problems and is derived straightforwardly from the ideal gas law. Viscosity of a single gas is given by
(9)
where N_{A} is Avogadro’s number and is a tabulated function depending on the LennardJones parameters for the species in question; the mixture viscosity is calculated from
(10)
with
. (11)
Thermal conductivity is given by
(12)
for monatomic gases; the Eucken correction is used for polyatomic gases:
. (13)
Mixture conductivity is derived from expressions having the same form as those used for viscosity.
2.1.7 Ordinary diffusion coefficients
Binary diffusion coefficients for the gas pairs are also obtained from kinetic theory, although a correction factor suggested by Wilke and Lee is also included; the result is
(14)
where is the average collision diameter for the pair and _{Dij}(T) is a tabulated function of the LennardJones parameters and temperature.
2.1.8 Thermal diffusion coefficients
Thermal diffusion modelling requires the use of multicomponent diffusion coefficients that are dependent on mixture composition. Derivation from kinetic theory is computationally expensive and increases as the cube of the number of gas species present; the problem lies in the evaluation of large determinants having elements which are themselves complicated functions of gas composition and temperature. As an alternative to this, a simpler model proposed by Clark Jones has also been inc luded in PHOENICSCVD; use of this method results in
(15)
where the thermal diffusion factors (a_{ij}) are functions of molecular masses and LennardJones parameters.
The user can select either the exact or the Clark Jones approach; as a further computational economy each formulation can be simplified by the use of the rigid elastic spheres approximation instead of the LennardJones potential.
2.1.8 Data input
PHOENICSCVD is equipped with two data files for transport and thermodynamic parameters. The user is able to add new species with complete freedom, following a standard format.
2.2 Chemistry Modelling
The fundamental process in Chemical Vapour Deposition is the chemistry that brings about the creation of the required film on the solid surface; inevitably this involves surface chemistry but generally there are also gasphase reactions. Many simulation studies have relied on very simple models for the chemistry, with rate constants estimated from experimental data; while this approach might yield accurate predictions for the process conditions under which th e data were collected it is of limited value in the main field of predictive modelling. A deeper understanding of chemistry is also essential for studying conformality, dopant incorporation and the important effects that intermediates can have on process characteristics such as selectivity and deposit morphology. More detailed chemistry models have appeared recently, but they have typically only been applied in conjunction with idealised reactor configurations.
PHOENICSCVD has been designed to facilitate the handling and solution of more sophisticated chemistry formulations, in addition to the simpler models, for both gasphase and surface reactions. Numerical modelling can then be used to gain fundamental insight into the underlying phenomena.
2.2.1 Gasphase chemistry
Although not contributing directly to the deposition process itself, gasphase chemistry plays an important role, particularly in atmospheric pressure reactors; the consumption of reactants and the formation of intermediaries may easily influence film quality.
The k^{th} reversible gasphase reaction is written as
(16)
where n _{ik} = 0.5(n _{ik} +n _{ik} ).
In modelling terms the gasphase chemistry contributes the final terms to the species and energy equations ((2) and (3)). The net reaction rate,  is obtained from
(17)
where k_{k} and k_{k} are the forward and reverse reaction rate constants respectively. These constants are strongly dependent on temperature and independent of pressure at sufficiently high pressures. At lower pressures though there is pressure dependence.
A number of alternative formulations for the forward reaction rate constant have been included in PHOENICSCVD. The simplest is the extended Arrhenius expression:
k_{k} = A Tb exp(E_{A}/RT) p^{c} . (18)
More complicated alternatives are the Lindemann form which blends two Arrhenius expressions, one for low and one for high pressure. Further complexity is also available in the Troe form, based on nine or ten parameters.
For reversible reactions the reverse reaction rate constant is calculated from equilibrium thermochemistry, utilising the equilibrium constant for the reaction; this is derived from polynomial fits for standard heat of formation, standard entropy and specific heat.
2.2.2 Surface chemistry
Surface chemistry can be written in the same way as the gasphase reactions described above:
(19)
where the A_{i} (i=1,N) are the gaseous reactants and products and the B_{j} (j=1,M) are the bulk, or adsorbed reactants and products; s _{il} and c _{jl} are the stoichiometric coefficients for the l^{th} surface reaction.
The growth rate is simply calculated from the reactions that deposit the bulk species:
(20)
where m_{s} is the molar mass of the bulk species and r _{s} is its density
Surface reactions did not appear explicitly in the species and energy equations ((2) and (3)) earlier. This is because they are included as boundary conditions rather than source terms.
2.2.2.1 Semiempirical surface chemistry models
In many chemistry models the surface chemistry is simplified into one or more irreversible reactions with only gaseous reactants and, apart from the deposited species, only gaseous products. Equation (19) then simplifies to
. (21)
Semiempirical methods have been developed for the definition of the surface reaction rates and two of these have been included as builtin options. The first involves the ‘reactive sticking coefficient’, defined as the fraction of gas molecules colliding with the surface that contribute to the deposited film. The reaction rate then becomes
(22)
where p_{A} is the partial pressure of the reactant and T_{s} the temperature at the surface; the sticking coefficient is g _{A} (£ 1) and the remainder of the expression comes from kinetic theory. Arrhenius expressions can be used for the sticking coefficient, specified via the chemistry data file.
A second builtin formulation is that by LangmuirHinshelwood; here the reaction rate is in the general form
. (23)
Again, the coefficients are specified in the data file.
2.2.2.2 Detailed surface chemistry models
The simplified surface chemistry models described above take no account of the adsorption of reactants and intermediaries at the surface, their decomposition and the desorption of reaction products. Although several more sophisticated models have been suggested in recent years there is still no standard form for them. In general, therefore, the inclusion of new models requires usercoding to be added to PHOENICSCVD; this facility is provided and a number of specific models have been incl uded in the code to help users.
The solution strategy for detailed surface chemistry models is rather different from that usually employed. At the surface the fractional coverages of the adsorbed species are calculated using a separate solver for the simultaneous calculation of stiff surface chemistry; the resultant fluxes of reactants, intermediaries and products are then imposed in the normal way as boundary conditions to the transport equations for the gas species.
2.2.3 Numerical complications
The solution of the transport equations for CVD processes is often hampered by numerical difficulties. These result from the chemistry source terms which lead to strong coupling between the mass fractions of different species at a single location; transport equations are usually solved on the basis that the strongest links are given by the convection and diffusion terms which link concentrations throughout the solution domain. Additionally, a single chemistry model is likely to contai n reactions with vastly differing time scales; this ‘stiffness’ can impose unacceptable limitations on the convergence that can be achieved.
Various techniques have been suggested to alleviate these problems, but they cannot usually be applied to two and threedimensional elliptic problems. In PHOENICSCVD convergence is improved by the linearisation of the chemistry source terms which increases the stability of the solution method. The stiffness is addressed by the provision of an automatic underrelaxation based on local creation and destruction rates. In extreme cases this may still not be sufficient and an alternative solver has been provided; this solves the species equations on a pointbypoint basis, using a NewtonRaphson solver for the simultaneous solution of all the species in a single computational cell. Such a technique is wellsuited for handling large sets of chemical reactions with widely varying time scales; however, because it is a ‘local’ method convergence may still be slow because the convective and diffusive contributions are less well dealt with. In some cases the use of both solvers in turn may be beneficial.
2.2.4 Data input
Specification of the chemistry to be adopted is made easy by the use of a data file. For each of the numbered reactions the user is able to specify the chemical reaction (in terms of species and stoichiometric coefficients), the nature of the model to be adopted (e.g. kinetic limited, diffusion limited, LangmuirHinshelwood) and the appropriate values of the parameters in the model.
2.3 Radiation Modelling
Radiation is a major element in heat transfer within CVD reactors, in addition to convection, conduction and diffusion. The nature of radiation is such that it does not fit easily into the framework of computational fluid dynamics, because of the links it introduces between remote parts of the solution domain. Radiation terms therefore appear as source terms calculated in a totally different manner from the other mechanisms of heat transfer. CVD gases are typ ically transparent to radiation and a surfaceto surface model is therefore sufficient. There are two basic formulations: viewfactor based or Monte Carlo. Both approaches were successfully adopted within the ACCESSCVD Project, although the Monte Carlo coding was restricted to axisymmetric geometries. Time constraints made it impossible for the Monte Carlo coding to be included in the first commercial release of PHOENICSCVD; it is intended to add it to a later code update.
In each approach the surfaces in the reactor are divided into zones within which the temperature and surface optical properties are assumed to be constant. The radiative heat source to the i^{th} zone is given by
(24)
where T_{j} is the surface temperature of the j^{th} zone, s is the StefanBoltzmann constant and R is the radiative exchange matrix. This heat source is used together with the other fluxes to calculate the zone surface temperature iteratively as the general solution proceeds. It is the radiative exchange matrix that is calculated by the viewfactor or Monte Carlo approach.
2.3.1 Viewfactor model
The viewfactor from the j^{th} to the i^{th} thermal zone is given by
(25)
where r is the distance between the points in each zone, n_{i} and n_{j} are the two surface normals and A_{j} is the surface area of the j^{th} zone; is the transmittance from j to i in the l^{th} spectral band. If all surfaces are opaque the viewfactors are purely geometrical properties of the reactor; the introduction of semitransparent solids means that the viewfa ctors must be recalculated if the degree of transparency changes, e.g. as a result of a changing surface temperature.
Calculation of the energy flux from j to i in the l^{th} spectral band is calculated using the Gebhard factors
(26)
which take account of direct and multiply reflected radiation; r^{(l)}_{m} and a^{(l)}_{m} are the reflectivity and absorptivity of the m^{th} zone in the l^{th} band. The net flux in the l^{th} band is then given by
(27)
where is the fraction of black body radiation in the l^{th} band. The total radiative heat source, as given above, includes the full radiation exchange matrix which is the sum of the matrices for each spectral band.
The surface temperature of each zone is calculated iteratively as the solution proceeds and the surface temperature used for the viewfactors should take account of this. For simplicity a userspecified value for each zone is used to avoid frequent recalculation of the radiation exchange matrix. Greater accuracy can be achieved if the specified surface temperatures are modified manually and the solution restarted; typically a single iteration of this type is sufficient.
The calculation of the viewfactors between zones is carried out by means of a raytracing algorithm which detects obstructions (opaque or semitransparent). The geometrical viewfactors are then derived using a modified double integral method which takes account of singular integrals when zones are in contact. For computational efficiency the spectral variation of surface properties is implemented in a banded manner; this is sufficient for most purposes.
2.3.2 Monte Carlo model
The Monte Carlo method is the same in principle; as before a radiative exchange matrix is generated from which heat sources can be determined. However, in this case the calculation is based directly on rays traced from each zone in turn until, after reacting with other surfaces, they reach their final destinations. A large number of rays are used, sampled from the appropriate angular and spectral distributions for the emitting (and any intermediate reflecting and transmitting) surface s, so that a statistically significant picture can be built up, showing how much energy from a surface will be received at each other zone. An advantage of the Monte Carlo approach is that it permits specular reflection to be simulated; this can be particularly useful for the highly reflective surfaces in reactor lamphouses. Furthermore, the spectral variation of properties can be exact, rather than in discrete bands.
2.3.3 Data input
A data file is used to specify the optical properties of any materials that are used in the simulation. The file contains, for each numbered material, a set of coefficients used in the calculation of the real and imaginary parts of the refractive index for each of sixty spectral bands; these cover the wavelength range from 10^{7}m to 10^{4}m, representing temperatures from 250K to 6000K. Simplified nonspectral modelling is also permitted and for this the data file c ontains constant values of emissivity and reflectivity. Additional materials can easily be added by users.
2.4 Plasma Modelling
In recent years the use of plasmaenhanced chemical vapour deposition (PECVD) has become more common; the ability to provide the energy required for CVD chemistry without the undesirable sideeffects of high thermal stress and poor controllability that are often a problem in conventional, thermallyenergised systems makes it an attractive option for a range of processes.
The modelling of plasmabased chemistry is inherently more complex than the modelling of comparable neutral gas processes. Apart from the need to model additional variables associated with the electrons, ions and electric field, there is the notably nonequilibrium state of the plasma: the energy is taken up much more readily by the light electrons than the heavier ions and neutrals which typically remain close to the ambient temperature.
Various methods for plasma modelling exist, usually based either on fundamental laws of physics or engineering approaches containing mainly empirical information. The first category are typically very computer intensive while the second yield only qualitative results. The perceived need within the ACCESSCVD Project was for a model that was simple enough to be compatible with the computational times for the rest of the simulation but still able to produce valuable quantitative results. Th e outcome was the Eddy DriftDiffusion Model (EDDM) which has been implemented in PHOENICSCVD. This model is explicitly geared towards capacitively coupled discharges in the PECVD regime; within those limits it is equivalent to the more fundamental fluiddynamical approach, from which it is derived, but without such demanding computational requirements.
2.4.1 Fluiddynamical plasma modelling
Unlike general kinetic models which describe a plasma in terms of distribution functions, fluiddynamical plasma models use macroscopic variables such as densities, bulk velocities and pressures. Such a model for an RF discharge would require equations for electron density and energy (or temperature), density and velocity for each type of ion and also electrical potential. With appropriate boundary conditions these partial differential equations would describe the evolution of a low p ressure discharge under an applied electric field. For periodic RF excitation it would be necessary to solve the transient equations until the asymptotic limit cycle was reached; the transient effects, typically lasting only a few milliseconds, are of little interest for CVD plasma processes.
This approach has been shown to be very reasonable, but extremely slow. One reason for this is that stability and accuracy require time steps in the transient calculation which are small compared with the driving RF frequency; as a large number of cycles are required before the limit cycle is reached this leads to unreasonable run times. This is another example of numerical stiffness (see chemistry section above): the electrons and the ions are described within the same framework, yet the dynamical phenomena act on vastly different time scales.
2.4.2 The Effective DriftDiffusion Model
The Effective DriftDiffusion Model overcomes the stiffness problem by making use of additional simplifications to the fluiddynamical equations; further, the limit cycle is calculated directly, with no attempt to model the transient evolution. The simplifications are based on the very circumstances which have caused the numerical stiffness in the first place: the pronounced scale separation of the dynamics. Specifically the following assumptions are introduced:
 the characteristic length scales of the plasma (Debye length, mean free path) are small compared with the scale of the reactor; the discharge can then be separated into the bulk, which is numerically resolved, and the thin boundary sheaths at electrodes and walls, which are represented by boundary conditions
 the characteristic frequencies of the electrons (particularly the dielectric relaxation frequency) are high compared with the applied frequency; the electrons are then quasistatic, or in Boltzmann equilibrium
 the dynamic frequencies of the ions (ion plasma frequency, ionization rate) are small compared with the RF frequency; the ions can then be regarded as stationary, reacting only to average field values.
With these assumptions the plasma can be modelled in three parts: plasma transport, RF modulation and energy balance.
Plasma transport:
The quasineutrality of the plasma, resulting from the speed of the dielectric relaxation, means that the electron charge density equals the net ion charge density:
(28)
where e and n are the electron charge and number density, and qa and r a the charge and density of the ions of species a .
Under such conditions a residual, ambipolar field, E, offsets the diffusion pressure of the electrons:
(29)
where m and D are the electron mobility and diffusivity.
This in turn drives the ions and is responsible for the plasma transport:
. (30)
Radio Frequency Modulation:
Superimposed on the DC transport is the RF modulation. The electrical current density is given by Ohm’s Law,
(31)
and is subject to Kirchhoff’s Law, as a result of which there are no volume sources (ionisation in the bulk is inefficient enough to be neglected):
Ñ .d j = 0. (32)
As the RF modulation is comparatively weak the potential can be replaced by the leading terms in a Fourier expansion
(33)
where w _{RF} is the driving frequency
The real and imaginary parts of the complex potential therefore satisfy Laplace’s equation.
Energy balance:
The energy balance describes the production , transport and dissipation of the periodaveraged electron energy T_{0}:
(34)
where k is the thermal conductivity and U is the averaged energy loss term.
The model is completed by appropriate boundary conditions taking account of the sheath region. The first of these deals with ion velocity, where Bohm’s criterion is assumed:
. (35)
The sheath region is assumed to be depleted of electrons and unable to support drift or diffusion. Instead the current between bulk and electrodes is carried by the displacement current, with the sheath acting as a capacitor. With C_{sh} as the sheath capacitance and V as the electrode voltage this results in
. (36)
The electron depletion in the sheath also means that the plasma bulk is thermally insulated from the walls, i.e.
(37)
The above set of equations and boundary conditions define the EDDM.
2.4.3 Comparison with a fluiddynamical model
Although the EDDM is derived from a fluiddynamical model the method of implementation is very different and it is important that its behaviour should be as expected. Testing was carried out on a 1D parallel plate reactor containing an electropositive (argon like) gas at a pressure of 250mTorr and with an applied RF voltage of 150V.
The electron and ion density for the fluiddynamical model were compared with the plasma density for the EDDM (which assumes identical ion and electron densities). The results were qualitatively identical and quantitatively reasonable in the bulk region; of course, the EDDM makes no predictions in the sheath region.
Comparison of the averaged electron temperature showed that in each case the electron temperature increased towards the electrodes, in accordance with observation; this is due to the enhanced Ohmic heating in this region. In this case the quantitative agreement was very satisfactory.
The benefit of the EDDM is apparent from the fact that it required approximately two hundred times less computational time to converge.
2.4.4. Implementation in PHOENICSCVD
The EDDM as incorporated in PHOENICSCVD has been simplified to contain only a single ion species. This restriction simplifies the model and results in the solution of only four plasma equations in addition to the conventional transport equations. For most applications
there is a dominant carrier species which can reasonably be used to determine the parameters needed to define the plasma behaviour.
2.5 Use of PHOENICSCVD
To facilitate the use of PHOENICSCVD by CVD engineers with little, if any, CFD experience it has been equipped with a menudriven user interface which makes problem setup straightforward. The menu system is based on the standard PHOENICS menu system, with additional sections for the specification of options which are specific to CVD applications.
Special attention has been devoted to the requirements of the surfaceto surface radiation model. This requires that all solid surfaces are covered by thermal zones, a task that is far from straightforward for the user. The menu system has therefore been designed to provide a starting point from which modifications and refinements can easily be made.
For more experienced users the alternative input using a PHOENICS Q1 file is still available, either as an alternative or a supplement to menu driven operation.
Additional data required by PHOENICSCVD have already been mentioned in earlier sections; these are provided by the user in data files describing transport and thermodynamic properties, chemistry formulations and optical properties.
Extra features incorporated in PHOENICSCVD are (i) the provision for the user to specify inflowing gas streams in conventional engineering units (sccm) rather than the standard PHOENICS mass flowrates, (ii) the introduction of special coding to allow the inclusion of shower plates as commonly used in CVD reactors, and (iii) the ability to model wafer batches approximately using a smaller number of wafers.
3. FURTHER SOURCES OF INFORMATION
Throughout the ACCESSCVD Project opportunities were taken to promote awareness of the work being done. This primarily took the form of conference and journal publications:
Author (s): 
Chr Werner 
Title: 
Advanced models for reactor simulation in silicon IC manufacturing. 
Presented at: 
Invited paper at the Topical Conference on the Synthesis and Processing of Electronic Materials (San Francisco, California, November 1994). 


Author (s): 
K J Kuijlaars, C R Kleijn and H E A van den Akker 
Title: 
A detailed model for low pressure CVD of tungsten. 
Presented at: 
International Conference on Metal Coatings and Thin Films (San Diego, California, April 1995): also accepted for publication in Thin Solid Films. 


Author (s): 
K J Kuijlaars, C R Kleijn and H E A van den Akker 
Title: 
Modelling of transport phenomena and detailed chemistry in chemical vapour deposition equipment. 
Presented at: 
1995 Meeting of the Electrochemical Society (Reno, Nevada, May 1995). 


Author (s): 
R P Brinkmann, R Fuerst, Chr Werner and M Hierlemann 
Title: 
A reduced fluiddynamical discharge model for applications in technologyoriented CAD. 
Presented at: 
First Symposium on Process Control, Diagnostics and Modeling in Semiconductor Processing, 1995 Meeting of the Electrochemical Society (Reno, Nevada, May 1995). 


Author (s): 
A Kersch 
Title: 
Development of RTPcontrol with Monte Carlo simulation. 
Presented at: 
Invited paper at the First Symposium on Process Control, Diagnostics and Modeling in Semiconductor Processing, 1995 Meeting of the Electrochemical Society (Reno, Nevada, May 1995). 


Author (s): 
T Schafbauer and A Kersch 
Title: 
Temperature control in RTP using reduction of equipment models. 
Presented at: 
First Symposium on Process Control, Diagnostics and Modeling in Semiconductor Processing, 1995 Meeting of the Electrochemical Society (Reno, Nevada, May 1995). 


Author (s): 
R P Brinkmann 
Title: 
A unified boundary sheath/presheath model for DC and RF driven low pressure discharges. 
Presented at: 
IEEE International Conference of Plasma Science (Madison, Wisconsin, June 1995) 
Author (s): 
R P Brinkmann and R Fuerst 
Title: 
Aspects of an efficient plasma model for the simulation of RF driven gas discharges in the PECVD regime. 
Presented at: 
International Conference on Phenomena in Ionized Gases (Hoboken, New Jersey, July 1995). 


Author (s): 
R P Brinkmann and R Fuerst 
Title: 
An efficient twozone model of capacitively driven low pressure gas discharges. 
Presented at: 
12th International Symposium on Plasma Chemistry (Minneapolis, Minnesota, August 1995). 


Author (s): 
A Kersch, HJ Timme, T Schafbauer and A Ajmera 
Title: 
Equipment simulation of openloop rapid thermal processing. 
Presented at: 
RTP '95 (Amsterdam, August 1995). 


Author (s): 
A Kersch and M Schaefer 
Title: 
Modeling the wafer temperature in a LPCVD reactor. 
Presented at: 
SISDEP '95 (Erlangen, September 1995). 


Author (s): 
M Hierlemann, C Werner, H Simka, K F Jensen and M Utz 
Title: 
Kinetic modeling of the gas phase decomposition of germane by computational chemistry techniques. 
Presented at: 
EURO CVD 10 (Venice, September 1995). 


Author (s): 
R P Brinkmann, R Fuerst and Chr Werner 
Title: 
Fast twodimensional modeling of RFdriven low pressure discharges in the PECVD regime. 
Presented at: 
11th International Conference on Gas Discharges and their Applications, GD ‘95 (Tokyo, September 1995). 


Author (s): 
K J Kuijlaars, C R Kleijn, H E A van den Akker, T G M Oosterlaken and G C A M Jansen 
Title: 
Modelling of blanket and selective tungsten LPCVD. 
Presented at: 
Advanced Metalization and Interconnect Systems for ULSI Applications (Portland, Oregon, October 1995). 


Author (s): 
Chr Werner, M Ilg and K Uram 
Title: 
3D equipment simulation for chemical vapor deposition. 
Presented at: 
AVS Meeting (Minneapolis, Minnesota, October 1995). 


Author (s): 
C R Kleijn 
Title: 
Numerical modelling of Chemical Vapour Deposition as a tool for process optimisation and reactor design 
Presented at: 
Paper published in Ultra Clean Society journal of Technology (Tokyo, Japan, January 1996). 
At the end of the Project a special twoday Seminar and Workshop was held in Wimbledon (September 2829, 1995). During the Seminar presentations were made by the partners, describing all aspects of the Project: model development, validation and pilot implementation. The presentations were later published in a special edition of The PHOENICS Journal (December 1995, Volume 8, No 4), devoted exclusively to the ACCESSCVD Project.