This User Guide accompanies the PHOENICS-CVD code (Version 2.1, PHOENICS 2021). Its purpose is to enable users to set up and run CVD simulations.
PHOENICS-CVD is an add-on to the general PHOENICS code and matches it in general structure; the experienced PHOENICS user will therefore only require the information about the additional CVD features of the code which this guide contains. Less experienced users may require help with respect to the use of PHOENICS itself; they should refer to the standard PHOENICS documentation (including the on-line information system, POLIS).
The features of PHOENICS-CVD will first be described in general terms, with a description of the different options available. Activation of the features will then be addressed, followed by a discussion of examples that provide the basis from which further cases can be developed.
The code is designed to model the behaviour of CVD reactors. This necessitates the modelling of fluid flow and heat transfer relating to a multi-component gas, including both gas-phase and surface chemical reactions. Additionally, the incorporation of plasma effects can be included.
The code can simulate up to 30 gas species undergoing multi-component diffusion in addition to the conventional convection and diffusion effects. Heat transfer by convection, diffusion and conduction is linked with surface-to-surface radiation, including the effects of semi-transparent solid materials. Chemical reactions, both gas and surface, are included by means of a data file containing information about the reactions. Up to 50 reactions of each type are permitted. Plasma-enhanced CVD modelling is achieved using an effective drift-diffusion model, providing a computationally efficient means of modifying chemical reaction rates to take account of the plasma.
Note that the restrictions on numbers of gas species and reactions can be modified by straightforward coding changes and recompilation; however, computational speed is likely to be a restricting factor.
The code also has a menu-driven user-interface, rendering it easy to use without the need for extensive CFD experience.
A range of options is provided for each of the features provided within the code. The theoretical basis of these different approaches is described in Appendix 1; this document is primarily related to implementation.
- Diffusion
Multi-component diffusion can be modelled in several ways. In increasing order of accuracy (and computational expense) these are Fick, Wilke and Stefan-Maxwell models; in each case the binary diffusion coefficients can be based on either the actual gas temperature or a specified reference temperature.
Thermal (Soret) diffusion is included using either the Clark Jones or exact formulation, incorporating the rigid spheres approximation or the use of Lennard-Jones parameters. Again, the increasing accuracy of the more sophisticated models is achieved at the expense of computation time.
- Chemistry
For gas-phase chemistry the following options are provided for the reaction rate calculation:
- extended Arrhenius dependence
- Lindemann form
- Troe form (9 parameters)
- Troe form (10 parameters);
in addition, provision has been made for user coding.
For surface chemistry the choice is between the Langmuir-Hinshelwood formulation and a reactive sticking coefficient expression. The complexity of surface chemistry is such that many models do not fit easily into standard formats; in this case the use of user-defined reaction rates is again permitted and the Kleijn and Roenigk & Jensen models (for silicon and tungsten deposition respectively), among others, have been included in the code to exemplify this option. Complex surface chemistry, involving adsorbed species, is available using user coding; several examples are provided as guidance.
Specification of the chemistry is by means of a data file in which each required reaction is included, with the model type and all necessary parameters; this will be described further in Section 4.4.
- Material properties
- Radiation
- Plasma
The properties of the gas mixture (density, kinematic viscosity, specific heat and thermal conductivity) depend on the local composition, in terms of mass fractions of the different species, and the local conditions, primarily temperature and pressure. For economy, two simpler alternatives are also provided; in these the properties are based on the carrier gas, using either the actual temperature or a reference value.
A surface-to-surface radiation model based on viewfactors is included. The radiation contributions are used in the conventional conjugate heat transfer solution of the energy equation. The model has the ability to incorporate semi-transparent solid materials although the gas itself is assumed to be transparent. Optical surface properties are derived from information in a data file (see Section 4.5 for further detail) and can be dependent on both temperature and wavelength; the use of spectrally varying properties can be of great importance in CVD applications.
Plasma modelling is achieved by means of the effective drift-diffusion model. This involves making assumptions that are generally appropriate in CVD applications; it is then possible to simulate the plasma by solving equations for the electron density, the electron temperature and the real and imaginary parts of a complex potential (see Appendix 1 for more information). If the time-average potential (bias or DC voltage) is not zero, an additional equation can be solved for the DC voltage. In this implementation the parameters affecting the plasma are specified in advance; the plasma therefore affects the gas flow (through the chemistry, which will depend on the state of the plasma) but there is no reciprocal interaction.
The plasma model is only for capacitively-coupled rather than inductively-coupled devices.
Standard features of PHOENICS are also available within the PHOENICS-CVD code. Here, only the items specific to PHOENICS-CVD will be described; general PHOENICS information can be found by using POLIS, the PHOENICS On-Line Information System.
The primary means of problem specification is the 'q1' input file; it is this that is interpreted by the main PHOENICS processing module (EARTH). The 'q1' file is written in a special form using PIL (PHOENICS Input Language); this can be produced directly by the user or created automatically at the end of a menu session. Although it is possible to use the menu alone without having any understanding of the 'q1' file it is strongly recommended that users develop at least a working knowledge of 'q1' file structure. In this way it is possible to verify that problem details have been correctly set-up within the menu and also to modify settings without the need to re-enter the menu each time. Additionally there are some more advanced features of the code that require the user to make 'q1' settings outside the menu framework.
In addition to the 'q1' input file use is made of data files specific to the individual code features.
This section will cover in turn: problem set-up using the menu, 'q1' file settings for CVD problems and data specification in the data files. Examples of 'q1' files are provided in the CVD Input Library.
- Use of the menu
- Preliminary activities
Before starting to use the menu it is advisable to make sure that all the information that will be needed is available. The geometrical information is particularly important and thought should be given at this stage to the way in which objects in the solution domain will be defined. Additionally, the various data files should be checked to ensure that they contain the appropriate data. These files are described in detail in Sections 4.3-4.6. The built-in material property file (phoenics/d_earth/props) should contain the required solids.
- General menu points
The menu is made up of two linked parts; the first is specific to CVD applications and the second is based on the general PHOENICS menu system. This reflects the fact that many of the settings required for a CVD simulation are also required for other types of modelling. The most efficient procedure is to define the gas species immediately (because they will be needed in the general menu for initialisation and the setting of boundary conditions); the plasma solution should also be activated if required. The remaining CVD settings are best left until the geometry has been defined.
- Starting a new CVD case
- Starting the CVD Menu
- The CVD Menu Panels
The main CVD panel looks like this (the date may differ):
The items on this panel are as follows:
Reference pressure sets PRESS0, the reference pressure. This is always added to the solved pressure, P1, before use in obtaining pressure-dependent properties. It can also be set from the Main Menu, Properties panel.
The first four ‘Settings’ sections should be visited before the case geometry has been created. Some of the settings made here are required for inlet and outlet boundary conditions.
It is of course possible, though less efficient, to create the geometry, then activate the variables here and re-visit the boundary conditions which need to be set.
The radiation section should be visited after the geometry creation is complete, as radiative transfer zones must be created on all internal and external surfaces.
These sections are now described in turn.
Gas species ‘Settings’ leads to this panel:
From here the species database can be inspected, and the required species selected from a scrolling list by clicking ‘Select’. The list is read from the SPECIDAT data file.
Selections are made using standard Windows techniques – <ctrl>+click for multiple selection, <Shift>+click to select a range.
The solution of the selected species is automatically activated on exit from the CVD part of the menu.
Clicking ‘View’ displays a list of the currently selected species. They can be de-selected here.
All but one of the selected species are solved for. The remaining one is deduced from continuity. It is usual to select the most prevalent species as the deduced, or carrier species.
The carrier is selected from the list of selected species by clicking ‘Modify’:
Diffusion and gas properties ‘Settings’ leads to this panel:
The thermal (Soret) diffusion can be activated, and the thermal diffusion option selected:
The options for thermal diffusion are, in order of increasing complexity:
- Clark-Jones form using rigid spheres approximation
- Clark-Jones form using Lennard-Jones parameters
- Exact form using rigid spheres approximation
- Exact form using Lennard-Jones parameters
For economy an update frequency is included: this is the frequency with which the thermal diffusion source terms will be updated.
The options for the multi-component diffusion treatment are:
- Fick
- Wilke
- Stefan-Maxwell
The gas mixture properties and binary diffusion coefficients can be evaluated at the actual local gas temperature, or at a fixed reference temperature.
Chemistry ‘Settings’ leads to this panel:
From here the reaction database can be inspected, and the required reactions selected from a scrolling list. The reaction database is read from the CHEMIDAT data file.
The G in column 1 indicates a Gas phase reaction, an S indicates a surface reaction.
Surface reactions take place on selected PLATE objects, and on selected faces of BLOCKAGE objects.
Under-relaxation based on the chemistry sources can be useful in avoiding convergence problems caused by strong chemical reactions; this is activated here if required, and the linear relaxation factor set.
The use of a stiff solver is possible for cases where the chemistry is much more important than the flow effects.
Surface-to-Surface Radiation Model
Radiation ‘Settings’ (for a new case, or a case with the Surface-to-Surface radiation model turned off) leads to this panel:
When the model is turned on, the panel changes to this:
When the radiation model is first turned on, radiation surfaces are created on all internal surfaces of BLOCKAGE and PLATE objects, and on the domain boundaries.
The default boundary conditions are taken from the values set at the top of the panel. If heat sources have already been set on the parent objects, these will be used instead of the default thermal condition.
The default material for optical properties is selected from a list, read from the OPTICDAT data file.
The default external thermal condition can be fixed flux or fixed temperature.
If the radiation model is turned off, the radiation zones are only deleted on exit from the CVD part of the menu.
All the zones can be deleted and new zones generated by clicking on ‘Generate’. This may be necessary if the geometry has been changed, and the existing zones no longer fit the surfaces of the objects.
The radiative surfaces appear in the Object Management Dialog and on screen as RAD_SURF objects. The default settings made can then be changed as required.
The Plasma ‘Settings’ panel looks like:
The settings for the equations allow the required parameters to be set.
The boundary condition patches are created from selected surfaces of BLOCKAGE and PLATE objects.
- Object Settings
In the CVD-VR system, a new object type RAD_SURF has been introduced. These objects are automatically created to cover all internal surfaces and exposed parts of the domain boundaries when the radiation model is turned on, or when the ‘Generate’ button is pressed on the Radiation – Settings panel.
The Attributes dialog for a RAD_SURF object looks like:
For ‘external’ surfaces (domain edge and on non-participating solids), the boundary condition can be Fixed temperature or Fixed Flux. On ‘internal’ surfaces (on participating solids) only a zero-applied-heat-flux condition is allowed.
The surfaces can be subdivided to make a finer radiation ‘grid’.
It is also possible to create RAD_SURF objects manually, just as any other object. In that case care must be taken to set the ‘Internal’/’External’ flag correctly, and the ‘Side’ flag such that the source is inside participating objects and outside non-participating ones.
Surface Chemistry and Plasma boundary conditions
Surface chemistry patches are created on the outer surfaces of BLOCKAGE or PLATE objects.
Clicking the CVD Settings button on a BLOCKAGE leads to:
From here, surface chemistry can be activated on any exposed face of the blockage. Note that:
West – East is in the X direction
South – North is in the Y direction
Low – High is in the Z direction
In a similar fashion, the plasma conditions for AC and DC can be selected. The default is ‘Insulated’. The other options are ‘Electrode’ and ‘Earthed’. The voltage at all ‘Electrode’ surfaces is set in the main CVD Menu Plasma settings pages. In addition to ‘Electrode’ and ‘Earthed’, DC has a fourth option, ‘User’, which allows the user to set the linking coefficient (‘Coefficient’) and voltage (‘Value’) or current directly. If the coefficient is set to FIXFLU, then the Value represents the current density in A/m2. Otherwise Value has units of Volts and Coefficient of A/m2/V.
For internal or external PLATE objects, the settings are directly on the Attributes dialog:
The ‘Plasma boundary conditions’ only appear if the Plasma model has been turned on from the CVD-Menu Plasma settings page.
- Menu Operation Example
The simplest way to understand menu operation is to work through one or more examples; in this case we shall use an ASM single wafer reactor being used for silicon deposition (a diagram of the geometry is provided as Figure 1).
| INLET FLOW |**********************|
C| PLUG PROFILE |**********************|
E|||||||||||||||||||||||**********************|
N|VVVVVVVVVVVVVVVVVVVVV|***** BLOCKAGE BL2****|W
T| |**********************|A
R| | |**********************|L
E| V | ___________________|L
| |**| |
L| |F*| |
I| | |I*| |
N| V |N*| |
E| |__| |
| |
| | |W
C| V |A
E| --> --> | |L
N| V |L
T| |
R| |
E|___________________________ |
|****** WAFER **************| | |
L|***************************|_______ V |
I|***********************************| |
N|***********************************| |W
E|***********************************| |A
|************************************| | |L
|************ BLOCKAGE BL1***********| V |L
|************************************| |
|************************************| |
|************************************| |
|************************************|PRESSURE|
|************************************| OUTLET |
To start a new CVD case in VR-Editor, click on File, Start new case and select CVD from the list of SPPs.
The data files CHEMIDAT, OPTICDAT, THERMDAT, SPECIDAT and TRANSDAT do not have to be copied to the working directory. If they are there, they will be used. If they are not, the copies in /phoenics/d_earth/d_opt/cvd/inplib will be used. This allows users to modify the files locally without the danger of damaging the default files.
The CVD menu is entered from the Main Menu, Models panel.
Existing CVD Q1 files will be recognised as CVD, and the relevant menu panels will be shown. In many cases, the individual surface chemistry and plasma boundary condition patches will be absorbed into the equivalent settings for the blockage or plate objects to which they refer. Radiative surface patches will be transformed into radiative surface objects.
Figure 1. Reactor geometry
Start the VR-Editor then click on ‘File – Start New Case – CVD’
Open the Main Menu (
on the hand-set
or 
on the
toolbar), then click on ‘Models, CVD settings’.Set the Reference pressure to 133Pa (equivalent to 1 torr).
Gas Species
Click on Gas species ‘Settings’ to define the gas composition. No species are currently selected so click on ‘Select’ to choose from the species data file (SPECIDAT).
Select the following species (keep the <ctrl> key pressed to enable multiple selection):
N2, H2, SIH4, SIH2, SI2H6, SI2H4, SI3H8.
As each is selected it will be highlighted. Make a note of the identification numbers: they will be used throughout to identify the species (with a preceding 'S' in the general part of the menu). At this point you may click View to check your selection.
One species must be selected as the carrier gas. The mass fraction of this species will be deduced from the others rather than being solved directly; it is also used in some simplifying assumptions if these are selected later. The carrier should be the dominant species (in mass or molar fraction terms). Click ‘Modify’ and specify N2 as the carrier gas. Click ‘Previous panel’.
Diffusion
Click on Diffusion and gas properties ‘Settings’ to select the models for multi-component diffusion and gas properties. Leave thermal diffusion inactive.
The default Multi-component diffusion setting is the Wilke approximation which is satisfactory. The binary diffusion coefficient option can be left as ‘Actual gas temperature’.
The default for the gas mixture properties is to use the gas composition at the prevailing temperature but two simpler options are also provided. Return to the main CVD panel.
Chemistry
Click Settings for Chemistry.. The first task is to choose the reactions: click ‘Select’ to see what is available. The list comes from the CHEMIDAT file: G refers to gas phase reactions and S to surface reactions. The Kleijn model for silicon deposition is represented by reactions 6, 7, 9, 10-16; select these and return to the menu.
As surface reactions are included we must define the surfaces on which they occur: this will be done later when the geometry is being defined.
We will now set a relaxation patch based on the chemistry (in addition to the flow based relaxation already activated): click in the tick-box by ‘Activate relaxation patch’ to activate the patch. The default factor of 0.5 is satisfactory so return to the main CVD panel.
Geometry and boundary conditions
At this point we need to create the geometry, so click Previous panel to return to the Top menu, then OK to leave the menu.
Turn the mesh toggle
on, and click in
the domain to bring up the Mesh Settings dialog. Switch to cylindrical-polar co-ordinates,
and set the domain size to:0.01 radians in X,
0.21m in Y, and
0.34m in Z.
Leave the grid set to AUTO in all directions. Close the Mesh Settings dialog.
The next task is to set up objects to which initial or boundary conditions will be attached; in the ASM geometry we need to identify solid regions (FIN & WAFER), blank regions (BL1 & BL2), inlet, outlet, wall and wafer as shown in the diagram.
To create the objects, click on Object (
on the hand-set
or 
on the
toolbar) to bring up the Object Management Dialog. On the menu bar of that dialog click
‘Object, New, New Object’.Set name to BL1. Set Size and Position as 0.01, 0.155, 0.1, 0.0, 0.0, 0.24
Click ‘General’ then ‘Attributes’. Leave the material as the default ‘Solid with smooth wall friction’. The volume covered by this object will be blanked out of the solution. Set a surface temperature of 290K on the North and Low faces. Click OK to close the Attributes dialog and OK to close the Object Specification dialog.
Click ‘Object, New, New Object’. Set name to BL2. Set Size and Position to 0.01, 0.11, 0.1, 0.0, 0.1, 0.0. Set a surface temperature of 290K on the South and High faces. Click OK to close the Attributes dialog and OK to close the Specification dialog.
Click ‘Object, New, New Object’. Set name to FIN. Set Size and Position to 0.01, 0.005, 0.04, 0.0, 0.1, 0.1. In General, Attributes, change the material to 111 Steel. Click on ‘Other materials’, select ‘Solids’, then ‘111 Steel at 27 deg C’. Click OK to close the Attributes dialog and OK to close the Specification dialog.
Click ‘Object, New, New Object’. Set name to WAFER. Set Size and Position to 0.01, 0.12, 0.005, 0.0, 0.0, 0.235. In General, Attributes, change the material to ‘111 Steel’. Activate Surface chemistry on the Low face. Click ‘CVD Settings – Surface Chemistry / Plasma’, then toggle ‘No’ to ‘Yes’ for the Low face. Click ‘OK’ to close the dialog. Set the heat source to Fixed temperature of 900K. Click on ‘Adiabatic’ next to Heat source. From the list of energy sources select ‘Fixed Temperature’, then enter 900 in the ‘Value’ input box. Click OK to close the Attributes dialog and OK to close the Specification dialog.
Click ‘Object, New, New Object’. Set name to INLET. Set Size and Position to 0.01, 0.1, 0.0, 0.0, 0.0, 0.0. In General, set the Type to INLET. In Attributes, set the inlet density to user-set, 1.573E-3 kg/m3, the temperature to 290K, the Z velocity to 1.344 m/s, and the inlet value of S140 (SIH4) to 0.113. Click OK to close the Attributes dialog and OK to close the Specification dialog.
Click ‘Object, New, New Object’. Set name to OUTLET. Set Size and Position to 0.01, 0.055, 0.0, 0.0, 0.155, 0.34. In General, set the Type to OUTLET. Click OK to close the Attributes dialog and OK to close the Specification dialog.
Click ‘Object, New, New Object’. Set name to WALL. Set Size and Position to 0.01, 0.0, 0.24, 0.0, 0.21, 0.1. In General set the Type to PLATE. In Attributes, set a surface temperature of 290K. Click OK to close the Attributes dialog and OK to close the Specification dialog.
Title
Click on Menu. On the top page, enter a suitable title for the case.
Initialisation
Solution is made more efficient if realistic initial values can be specified for the solved variables. In Initialisation set S116 (N2) to 0.887, S140 (SIH4) to 0.113 and TEM1 to 290.
The initial values for the species equations must be in mass fractions, not molar fractions.
Variables
On the Models panel, click on Settings by Solution control/Extra variables. This confirms that the correct velocity components, pressure, temperature and species are being solved; it also shows the stored variables, including the carrier species, the deposition rate variable, DEPO and PRPS (to define material properties); click Previous panel If additional variables were required this would provide the route for setting them.
Numerical Controls
In Numerics, set 400 sweeps. The default relaxation settings will suffice. (Note: plasma variables, if active, will usually require a very small setting of FALSDT. This is because they are unaffected by the gas flow that normally determined the appropriate relaxation level.)
Radiation
Go to ‘Models, CVD, Radiation Settings’. The optical property temperature determines the temperature at which the properties will be calculated: Set the default reference temperature for optical properties to 290K.
The optical property index defaults to 198, which represents black body. In this case the outer walls are steel, which has an index of 111: Change the default material for optical properties to ‘111 Steel (spectral set equal to iron)’.
The default setting for an 'outside' zone is Fixed heat flux of 0W (adiabatic). Change this to ‘Fixed temperature’ and set the default external temperature to 290K.
Activate radiation. 12 radiation surfaces should be created, each using the set default values. It will be necessary to modify some of these later. Leave the CVD menu and the Main Menu.
In the Object Management Dialog, select all the blockage, plate, inlet and outlet objects, and hide them. The radiation surfaces will now be visible. It may help to go to Settings, Editor parameters and increase the X scale factor to 10.
Select the 3 surfaces which cover the FIN object (B10-12), and from the right-click context menu open the dialog box. In Attributes, change the optical property temperature to 300. Close the dialog, and click YES to propagate the change.
Select the 2 surfaces covering the WAFER object (B13-14), and in the same way change the optical property temperature to 900.
Monitoring controls
The print-out controls can be left at their default settings but it is sensible to match the monitoring (probe) location to the problem. While the code is running the values of selected variables in the monitoring cell will be displayed. It is worth selecting a 'sensible' cell for the probe which will provide useful information on the progress of the solution.
Move the probe to 0.005,0.05,0.16. (The Position buttons on the hand-set move the probe when no object is selected).
Solution
Run the solver by clicking on Run - Solver.
The solution should look something like this in the VR-Viewer post processor:
This example has demonstrated the general use of the menu. There are many more options that have not been used and you might like to investigate them further. It is advisable to save the 'q1' file just created before doing so: click File – Save as a case and enter a descriptive case name. All the input and output files will be saved using the given name as a base.
- Preliminary activities
- The 'q1' input file
In this section the 'q1' settings required to activate the different model options are described; a working knowledge of the PHOENICS Input Language is assumed. Sample 'q1' files are included in the CVD Input Library and should be referred to in conjunction with the descriptions below. The first corresponds to the same problem as used for the menu description above but it was created without using the menu; experienced PHOENICS users might like to compare it with the 'q1' file created by the menu and convince themselves that the two, although different in format, are equivalent in content.
The CVD options are activated by setting NAMGRD=CVD; this line should therefore always be included.
- Variables
The primary variables are the mass fractions of the gas species. These should be identified as scalar variables between C1 and C30 inclusive (corresponding to integers 16 to 45) and it is essential that no other variables should occupy any of these allocations. One (and only one) of the mass fractions should be stored rather than solved; its value will then be deduced by the requirement that the mass fractions must sum to unity. It is recommended that the dominant species (in terms of mass or molar fraction) should be the one to be stored; the stored species will also be the one used internally if material properties are based on the carrier gas rather than the local mixture. All solved species should be solved whole-field by use of the SOLUTN setting.
Certain other variables must also be stored:
(i) the appropriate flow variables - pressure (P1), velocities (U1,V1,W1 depending on problem dimensionality) and temperature (TEM1)
(ii) material identifier (PRPS) and deposition rate (DEPO) in Angstrom/minute.
Additional variable names are also recognised if storage is activated for them in the 'q1' file; the appropriate values will be inserted automatically within the code and can then be printed in the 'result' file or viewed using the VR-Viewer post-processor. These are FR01-FRnn, MD01-MDnn, SP01-SPnn, SN01-SNnn, SPHT and KOND; they refer to molar fractions for each species, multi-component diffusion coefficients for each species, positive and negative gas-phase chemistry sources for each species, mixture specific heat and mixture thermal conductivity respectively. Here 'nn' is the total number of gas species. (Note that the molar fractions, diffusion coefficients and sources will be in an internal order - this will be different from the order used in the 'q1' file only if the first species named in the 'q1' file, C1, is not the one that is stored rather than solved.) Care must be taken when using the optional storage to ensure that only species mass fractions occupy the positions 16 to 45 in the variable list. Note that SPHT will contain an effective specific heat, defined as (h-h0)/T, where h is enthalpy, T temperature and subscript 0 refers to 250K.
The different gas species involved in the simulation must be identified using the NAME command. The available species are listed in the data file SPECIDAT which can be modified to suit the user's requirements (subject to a limit of 300 species and the proviso that the corresponding species must all appear in the transport and thermodynamics data files (TRANSDAT and THERMDAT). The name of a species is made up of 'S' followed by the identifying number in SPECIDAT (so S80 is hydrogen if the default SPECIDAT is used).
In the VR Editor
The solution for pressure, velocity, temperature and all the selected species is activaed automatically by the Editor, together with storage for DEPO and EMIS. Storage of other derived variables, such as the molar fractions, can be activated from Main menu – Models – Solution control/Extra variables’.
- Diffusion settings
- Material property settings
- Chemistry settings
- Radiation settings
- Plasma settings To activate the plasma model the ‘q1’ needs to contain the following lines:
- Additional features
- Summary of 'q1' settings
- Object-related settings
Settings are required in the 'q1' file to determine which of the various models will be used for multi-component and thermal diffusion. These make extensive use of the SPEDAT command, the use of which can be understood from the 'q1' files provided and the examples below.
In the first case the choice is between Fick, Wilke and Stefan-Maxwell. Multi-component diffusion is activated by the PRNDTL command for the mass fraction variables; the appropriate setting is -GRND8. Model selection is specified by the variable MCDOPT using the SPEDAT command - values of 1, 2 and 3 respectively should be used. For Stefan-Maxwell diffusion it is often necessary to implement linear relaxation, using SPEDAT to define the relaxation factor, SMRLX. The binary diffusion coefficient can be calculated using actual temperature or a reference value, selected by values of 4 and 2 respectively for BINOPT (again specified using SPEDAT); if a reference temperature is required it must be allocated to TMP1A.
Thermal diffusion options also make use of SPEDAT to define THMDIF, THMOPT and THMFRQ. The first of these should be set to T to activate thermal diffusion; the second determines the model to be used:
1 - Clark Jones form using rigid spheres approximation
2 - Clark Jones form using Lennard-Jones parameters
3 - exact form using rigid spheres approximation
4 - exact form using Lennard-Jones parameters.
These are in increasing order of computation time; as a computational economy it is possible to update the thermal diffusion source terms less frequently than every sweep by setting THMFRQ to a value greater than 1.
Both Stefan-Maxwell and thermal diffusion require the settings USOURC=T and UDIFNE=T; these settings should anyway be made to include the energy transport caused by diffusion fluxes, as should DIFCUT=0.0. (Note that it is possible, though not recommended, to leave out these settings if neither Stefan-Maxwell nor thermal diffusion is active; in that case the energy flux associated with diffusion fluxes will also be omitted.)
As an example the settings
SPEDAT(SET,CVD,MCDOPT,I,3)
SPEDAT(SET,CVD,SMRLX,R,0.5)
SPEDAT(SET,CVD,BINOPT,I,4)
SPEDAT(SET,CVD,THMDIF,L,T)
SPEDAT(SET,CVD,THMOPT,I,3)
SPEDAT(SET,CVD,THMFRQ,I,3)
(in addition to the appropriate PRNDTL, USOURC, UDIFNE and DIFCUT settings) activate thermal diffusion based on the exact form using the rigid spheres approximation with source updates every three sweeps, binary diffusion coefficient calculation based on actual temperature and Stefan-Maxwell diffusion with a relaxation factor of 0.5.
In the VR Editor
The Editor makes all the above settings based on the choices made in the Diffusion panel of the CVD Menu.
Material property settings are implemented by storing the variable PRPS and the specification of the selected property values in the 'q1' file. Density (RHO1), specific heat (CP1), laminar viscosity (ENUL), and thermal diffusion coefficient (PRNDTL(TEM1)) settings will then be taken from the property information table. The commands FIINIT(PRPS)=70 and CSG10='q1' enable the user to specify this property information in the 'q1' file (provided that property setting 70 has not been added to the 'props' file in the PHOENICS d_earth directory). The following indented lines should then be included in the 'q1' file:
MATFLG=T; IMAT=1
70 GRND8 GRND8 GRND8 GRND8 1.000
0 0
0 0
0 0
0 0
The method to be used for the calculation of mixture properties is given by a SPEDAT command for the variable MCPROP: 1 refers to the carrier gas at a reference temperature, 2 to the carrier gas at actual temperature and 3 to the actual gas composition at the actual temperature. Again TMP1A determines the reference temperature to be used, if required. One line of zeros should be included for each GRNDn entry in the property definition line. Additional materials can also be specified here if IMAT is suitably increased.
In the VR Editor
The Editor automatically inserts the material property lines for material 70, and makes it the domain material. The mixture property settings are made based on the choices made in the Diffusion panel of the CVD menu.
Gas-phase chemistry is activated by means of a PATCH which will normally cover the whole solution domain (solid regions will be excluded automatically by the code). The patch name should start with 'CHEM' and the type should be VOLUME. To allow for source linearisation both coefficient and value must be set to GRND1; COVAL statements are required for all species and the temperature TEM1.
Surface chemistry is activated by means of a PATCH having a name that starts with 'SURF' (note that there must be at least one additional character); the patch type should be the appropriate area (e.g. SOUTH) and the patch should only extend over the gas cells adjacent to the solid surface. COVAL statements should be provided for all species (including those that are not involved in the surface reactions), TEM1 and P1; this ensures that the Stefan velocity (flux to the surface) is correctly included. The coefficient should be set to FIXFLU and the value to GRND1 in all cases except P1; there the coefficient should be set to the 'batch factor'. For a conventional surface chemistry patch this will be 1.0. A higher number can be used if the computational surface represents several physical surfaces; by this means a batch reactor can be approximately simulated without modelling each wafer individually (useful if a parametric study is being undertaken for which run times would otherwise be unacceptable). If the coefficient is given an inappropriate setting (e.g. FIXFLU) a batch factor of 1.0 will be assumed.
Automatic under-relaxation based on the chemistry sources can be useful in avoiding convergence problems caused by strong chemical reactions; this is introduced by means of a PATCH called RELT having PHASEM type. COVAL statements should be provided for all species involved in the chemistry, with a coefficient of GRND1 and a value of SAME. A linear relaxation factor, CHMRLX, can then be specified using SPEDAT: values are between 0 and 1, with low values signifying tight relaxation. Note that conventional relaxation may also be needed, particularly for any species that do not participate in the chemistry. Linear relaxation is recommended for TEM1 and a small value (say 0.3) may be necessary to ensure stability.
The reactions to be included in the simulation are specified using their identification numbers in the CHEMIDAT data file; additionally, the total number of gas-phase and surface reactions, NGREAC and NSREAC respectively, must be given. In all cases SPEDAT is used. The individual reaction numbers are specified using GREAC and SREAC in the two cases. As an example,
SPEDAT(SET,CVD,NGREAC,I,2)
SPEDAT(SET,CVD,GREAC(1),I,6)
SPEDAT(SET,CVD,GREAC(2),I,7)
SPEDAT(SET,CVD,NSREAC,I,1)
SPEDAT(SET,CVD,SREAC(1),I,12)
declares two gas-phase reactions (6 and 7) and one surface reaction (12).
An output file 'data' is created by the code, containing information about the species present and the reactions that have been selected. This should always be checked for consistency; any mass imbalance in a selected reaction will in any case cause the code to issue a warning message.
In the VR Editor
The whole-domain gas-phase chemistry ‘CHEM’ patch is created automatically by the Editor. The surface-reaction ‘SURF’ patches are automatically generated from each BLOCKAGE and PLATE object based on the settings for the faces. The relaxation patch is created if the relaxation option is selected on the Chemistry page of the CVD menu.
Radiation is implemented by setting S2SR=T in the 'q1' file; additionally it is necessary to provide storage for variable EMIS which is used internally to define surface optical properties. To ensure that the correct surface-to-surface radiation formulation is used, the 'q1' file should also contain the line:
SPEDAT(SET, CVD, RADCVD, L,T).
The model is based on a surface-to-surface viewfactor approach, i.e. it is assumed that the gas is transparent. Thermal zones must be defined to cover all solid surfaces and viewfactors are then calculated between all pairs of thermal zones (taking account of any semi-transparent solid regions). From these a radiation exchange matrix is calculated that fully defines the energy transfer between surfaces in terms of the surface temperatures, which are calculated unless they are defined as fixed. Each thermal zone requires PATCH and COVAL commands to define the location, type and boundary conditions for the zone. The patch names must start with @RI; the patch type must always be an area, e.g. SOUTH. For zones on the surfaces of solids for which conjugate heat transfer is being employed (referred to below as internal surfaces) the patch should be on the solid side of the surface; otherwise (for example, at the domain edge) the patch should be on the gas side of the surface.
For each thermal zone there must be a COVAL command for TEM1; the coefficient and value determine the type of thermal boundary condition:
Coefficient Value Boundary condition
GRND1/2 GRND1 Zero applied heat flux (internal surface)
GRND1/2 VAL Fixed heat flux of VAL (external surface)
0.0 VAL Fixed temperature of VAL (external surface)
Zero applied heat flux is the only permitted setting for an internal surface. External surfaces (i.e. surfaces of non-participating solids) normally have a specified surface temperature or heat flux. GRND1/2 indicates that either GRND1 or GRND2 can be used: GRND1 is usually preferred although in some circumstances the less stable GRND2 can be more efficient.
In addition to its radiation patch, a fixed temperature thermal zone should be accompanied by an appropriate laminar wall boundary condition for temperature in the conventional manner (see examples); the menu will do this automatically.
A COVAL command is also necessary for the variable EMIS at each thermal zone. The coefficient determines the solid material according to the numbering system in the data file OPTICDAT. The value, if positive, specifies a reference temperature which will be used to calculate optical properties; a negative value is used to indicate the thickness of a semi-transparent coating. If a semi-transparent coating is indicated, the EMIS coefficient will determine the optical properties in the layer; the optical properties of the background material will be determined by the PRPS value there: it is therefore necessary to ensure that consistency is maintained between PRPS values and the material labels in the optical data file. The temperatures specified for the surface property calculation can only be a guess before the run has been carried out; for the highest accuracy it would be necessary to carry out additional runs for which the thermal zone temperatures are taken from the result file for the preceding run (these values are listed at the end of the result file). Typically one or two continuation runs would be sufficient.
In order that the accuracy of the viewfactor calculation can be seen by the user, two energy conservation checks are provided. The first relates to the viewfactors themselves and the second to the radiation exchange matrix; in each case the calculated value is printed out for each thermal zone (to the screen and also to the output file 'rad.dat' (see below)), together with the correct value. Errors of more than 1-2% should be regarded with suspicion: first check that the radiation patches have been correctly specified, then try increasing the number of test rays used (see below) and then increase the number of thermal zones until satisfactory agreement is achieved. In any case, each thermal zone should be defined so that the surface temperature does not vary significantly over it.
An output file 'rad.dat' is created by the code; in addition to information about the thermal zones that have been defined it contains the array of viewfactors between the thermal zones, the accuracy checks referred to above and the radiation exchange matrix elements. This file should be saved if further simulations using the same geometry and thermal zones are to be carried out. It can then be renamed as RADI.DAT; if a file of that name is present in the current directory when a run is started the information it contains will be used and no viewfactor calculation will be carried out. This can be a significant time-saving feature but requires care. Use of the wrong RADI.DAT will usually cause an error because the number of thermal zones will be incorrect; however, if the number of thermal zones in RADI.DAT is the same as that in the case being run the simulation will continue based on incorrect information. If a second run is carried out using the same thermal zone definition but with differing surface properties it is possible (provided that no semi-transparent solids are present) to use the previous viewfactor calculation as the basis for the derivation of the new radiation exchange matrix; to achieve this it is necessary to modify 'rad.dat' to remove everything below the viewfactor array listing before renaming it as 'RADI.DAT'.
Several additional settings are available within the surface-to-surface radiation model; these make use of the SPEDAT command. XMIR, YMIR and ZMIR are logicals which, if set to T, indicate symmetry in the corresponding co-ordinate directions (only available in 3-D); AXIBFC enables a 2-D BFC case to be specified as axisymmetric. By default the spectral variation of surface optical properties from OPTICDAT is used; the constant values in the file can be chosen instead by setting NOSPCT to T. In 3-D viewfactor calculations a coarse calculation is used in the construction of the viewfactors if an internal accuracy check is satisfied; this can result in view factor symmetries being only approximately satisfied. The logical FINE3D can be used to force the fine calculation under all circumstances; this should restore the symmetries, albeit at a cost in computation time. The viewfactors are usually renormalised before the radiation exchange matrix is derived as a means of ensuring good energy conservation; again, this may remove certain symmetries and lead to small anomalies in the results. VFNORM can be set to F to remove the renormalisation (energy conservation will then be imposed by modifications to the radiation exchange matrix instead). The accuracy of the viewfactor calculation can improved by increasing the number of test rays used in the determination of each zone-to-zone, or cell-to-cell, viewfactor; the number of rays is specified using NUMRAY. Note though that there is an internal limit of 5, and that NUMRAY has no effect in 3-D cases for which only one ray is used for each cell as a computational economy.
View factors in cyclic geometries can be calculated but only if BFC=T. If a polar grid is used it is necessary to carry out a preliminary calculation for which BFC=T is set in the 'q1' file after the geometry has been defined; the 'rad.dat' file can then be used as input for the real case (BFC=F) as already described.
In the VR Editor
The ‘@RI’ patches required to calculate the viewfactors are automatically generated from the RAD_SURF objects which are themselves generated when the radiation model is turned on, or whenever the ‘Generate’ button is clicked on the Radiation page of the CVD menu.
The other additional settings described can be activated from the ‘Other settings’ page of the CVD menu Radiation page.
SPEDAT(SET,CVD,PLASMA,L,T): Global switch for plasma model
SPEDAT(SET,CVD,LACPLS,L,T): Switch for AC plasma
SPEDAT(SET,CVD,LDCPLS,L,T): Switch for DC plasma
If the AC plasma model is activated it is necessary to solve for an additional four variables: T0, NE, PHI1 and PHI2; these are respectively the electron temperature, the electron density and the real and imaginary parts of the complex potential (note that electron density is quoted in units of 1016 electrons/m3 to avoid excessively large numbers). As for the other scalar variables, whole-field solution should be activated.
If the DC plasma model is activated (LDCPLS=T) it is necessary to solve for an additional variable EPOT. This is the DC or Bias voltage. As for the other scalar variables, whole-field solution should be activated.
Additional variable names are also recognised if storage is activated for them in the 'q1' file: PHIT and GAMM are the amplitude and phase angle of the complex potential and GION is the ionisation rate.
The diffusive nature of the equations in the effective drift-diffusion model requires the modification of the standard PHOENICS transport equations to remove the other terms; this is achieved by setting TERMS(...,N,N,Y,P,Y,N) for each plasma variable. The diffusion coefficient is set by PRNDTL(...)=-GRND7 in each case.
The source terms in the T0 and NE equations are set up using two patches which cover the whole domain (solid regions being dealt with internally); the patch names are 'COOL' and 'IONIZE'. Each patch has type VOLUME and COVAL statements are required for T0 and NE; in each case both coefficient and value should be set to GRND7. Other parameters in the plasma equations are specified using arrays T0PAR, NEPAR and PHIPAR in SPEDAT statements; the meanings are as follows:
NEPAR(1):Bohm velocity (m2/s)
NEPAR(2):Ionization rate - upper bound on Arrhenius term (/s)
NEPAR(3):Ionization rate - multiplying factor in Arrhenius term (/s)
NEPAR(4):Ionization rate - activation energy in Arrhenius term (K)
NEPAR(5):Ionization rate - field independent rate (/(m3/s))
T0PAR(1):Ohmic heating coefficient ((m2K)/(V2s))
T0PAR(2):Diffusion coefficient (m2/s)
T0PAR(3):Cooling source term - time constant (/s)
T0PAR(4):Cooling source term - equilibrium temperature (K)
PHIPAR(1):Sheath parameter (m4/1016 electrons)
PHIPAR(2):Electron mobility (m2/(Vs))
PHIPAR(3):Phase angle of externally applied AC potential (radians)
PHIPAR(4):Magnitude of externally applied AC potential (V).
PHIPAR(5):DC electrode parameter, Hall coefficient = n*Ze/dx[1/m3*C]
PHIPAR(6):DC bias Voltage (V)
In addition to the terms in the plasma equations it is also necessary to specify appropriate boundary conditions. For T0 no action is required; this ensures zero energy flux across boundaries. For the AC potentials, patches will be required for electrodes and earthed surfaces; all these patches should have the appropriate area type and COVAL statements for PHI1 and PHI2 with a coefficient of FIXFLU and value of GRND7. Electrode patches should have a name starting with 'ELE'; the externally applied potential referred to above will be applied there. All other PHI1 and PHI2 patches of area type will be at zero potential (corresponding to earthed surfaces). The names of these patches can start with ‘ERT’.
The DC potential EPOT can use the same ‘ELE’ and ‘ERT’ patches used by PHI1 and PHI2. At electrode patches, coefficient and value should both be GRND7. At Earthed patches the coefficient is still GRND7 but the value is zero. If a current density is to be specified, the coefficient should be FIXFLU and the value the current density (A/m2). In this case the name of the patch is immaterial.
NE boundary conditions should be specified at any surface through which electron flux is permitted. The patch type should again be an area, so the potential patches can be used for this purpose; the coefficient for NE should be GRND7 and the value 0.0.
Realistic initialisation for plasma variables can be very beneficial; in particular, an inappropriate value for T0 can have severe consequences because it is used in an exponential term in the electron density equation.
The plasma variables influence the chemistry and thus the flow; however, the flow only influences the plasma via the parameters in the plasma equations, which are dependent on gas composition. In the current release these parameters are taken to be constant and the plasma equations are then independent of the other variables. It is therefore possible to solve for the plasma first, in isolation, and then restart with solution for gas species and flow variables once a reasonable plasma solution has been obtained. This can be economical because no time is spent on the flow and chemistry solution while the plasma variables (and thus the reaction rates) are still unrealistic.
The use of linear relaxation for plasma variables can cause convergence difficulties when block correction is implemented. False timestep is then recommended; the appropriate value can be small, typically O(10-6).
In the VR Editor
The solution for the plasma variables, and the whole-domain ‘COOL’ and ‘IONIZE’ patches are automatically created when the plasma model is turned on. The appropriate boundary condition patches are created from BLOCKAGE and PLATE objects based on the settings made. The parameters in the plasma equations are set from the Plasma settings pages of the CVD menu. The plasma variables can be initialised from the Main menu – Initialisation page.
To facilitate CVD modelling two additional features have been provided; the first relates to inlets and the second to shower plates.
Inlets are normally defined in terms of mass flow (per unit area), with incoming velocity, temperature and species mass fractions also being specified. As an alternative it is possible to give the flow in terms of sccm (standard cubic centimetres per minute) over an area that must also be given; the incoming gas composition must then be defined in terms of molar, rather than mass, fractions. To select this option the patch name must begin with 'VIN'. COVAL statements are required for P1 and all other variables defining the incoming stream; the coefficients are FIXFLU for P1 and ONLYMS for other variables, as usual for an inlet. The value for P1 is the inflow (in sccm), the value for the velocity normal to the inlet face is the corresponding area and the values for the gas species are molar fractions. Thus
PATCH(VIN1,LOW,........)
COVAL(VIN1,P1,FIXFLU,100.0)
COVAL(VIN1,W1,ONLYMS,0.1)
COVAL(VIN1,S80,ONLYMS,0.2)
specifies an inflow of 100sccm over an inlet area of 0.1m2; the gas composition is 20% S80 (hydrogen) by volume. Note that this feature cannot be used in BFC cases.
Shower plates are common in CVD reactors and a special patch has been defined to model them. The region in question should first be defined to have the appropriate solid properties (using PRPS); an additional patch should then be defined in the same place having a name beginning with 'SHWR' and an area type that indicates the plane of the plate. COVAL statements are required for all solved scalar variables in addition to the velocity component normal to the plate. All values should be GRND7 and the coefficients should be ONLYMS, with the exception of those for P1 and the normal velocity which have a special meaning. The coefficient for the normal velocity should be the open area fraction (porosity) of the plate and that for P1 should be a loss parameter. This is defined as dp/(dynamic viscosity x mean velocity), where dp is the pressure drop across the plate; for circular holes this will be equal to 32l/a2, where l is the thickness of the plate and a is the diameter of the holes. As an example,
PATCH(SHWR1,HIGH,........)
COVAL(SHWR1,P1,32000.0,GRND7)
COVAL(SHWR1,W1,0.5,GRND7)
COVAL(SHWR1,S80,ONLYMS,GRND7)
COVAL(SHWR1,TEM1,ONLYMS,GRND7)
CONPOR(PLAT1,-1,........)
COVAL(PLAT1,PRPS,0.0,111)
represents a steel plate with 50% open area and a loss parameter of 32000 (possibly representing a plate 1mm thick with holes of 1mm diameter). Note that a new PRPS number should really be specifically defined so that the plate would have the properties that take account of its porous nature; the two patches should cover the same range of cells.
The use of a stiff solver is possible for cases where the chemistry is much more important than the flow effects; this can be activated by setting USOLVE=T. Using this approach involves solving for all species, including the carrier; the code will automatically make this adjustment but the user should ensure that all appropriate COVAL statements are provided for the species that is nominally stored (such statements would not be needed under normal circumstances).
In the VR Editor
The CVD sccm inlet using ‘VIN’ as a patch name has not been implemented in the CVD menu. The USER_DEFINED object type can be used to create such a patch if required.
The shower plate has not been implemented as a special object type. Again, the USER_DEFINED object can be used to create the patch commands listed above.
The stiff solver can be turned on from the Chemistry settings page of the CVD menu.
This section contains a review of the settings described in the preceding sections.
Group 7:
Variables to be solved or stored:
Flow variables, including: TEM1
Mass fractions C1-C30 (one stored, rest solved)
Radiation variable (if used) EMIS
Other required stored variables: PRPS, DEPO
Optional stored variables: FR01-FRnn, MD01-MDnn,
SP01-SPnn, SN01-SNnn,
SPHT, CNDT
Plasma variables (if used): T0, NE, PHI1, PHI2, EPOT
Use whole-field solver (default) for all solved scalar variables
Species identification:
NAME(C*) = Snnn (nnn is the identifier in SPECIDAT)
Group 8:
USOURC = T, UDIFNE = T, DIFCUT = 0.0
TERMS(...,N,N,Y,P,Y,N) for each plasma variable
USOLVE = T if stiff solver used
Group 9:
PRNDTL(Cn) = -GRND8
PRNDTL(...) = -GRND7 for each plasma variable
Property information as required, using CSG10 = 'q1'
Group 13:
Gas phase chemistry (only seen in Q1EAR and RESULT):
PATCH(CHEM*,VOLUME, ........ )
COVAL for solved species and TEM1
Surface chemistry (only seen in Q1EAR and RESULT):
PATCH(SURF*,area type, ........ )
COVAL for all solved species, P1 and TEM1
Under-relaxation based on chemistry:
PATCH(RELT,PHASEM, ........ )
COVAL for each solved species
Radiation (only seen in Q1EAR and RESULT):
PATCH(@RI***, area type, ........ )
COVAL for EMIS and TEM1
Plasma sources:
PATCH(COOL,VOLUME, ........)
COVAL for T0 and NE
PATCH(IONIZE,VOLUME, ........)
COVAL for T0 and NE
Plasma boundary conditions (only seen in Q1EAR and RESULT):
PATCH(ELE*, area type, ........)
COVAL for PHI1, PHI2, EPOT and NE
PATCH(***, area type, ........)
COVAL for PHI1, PHI2, EPOT and NE
Group 17:
Linear relaxation for TEM1
Group 19:
S2SR = T (if radiation activated)
SPEDAT(SET, CVD, RADCVD, L, T).
SPEDAT(SET,CVD,parameter,type,value) for (some of) the following parameters:
Diffusion:
MCDOPT integer
SMRLX real
BINOPT integer
THMDIF logical
THMOPT integer
THMFRQ integer
Gas properties:
MCPROP integer
Chemistry:
CHMRLX real
NGREAC integer
GREAC(1:n) integer(s)
NSREAC integer
SREAC(1:n) integer(s)
ANJAC logical
Radiation:
XMIR logical
YMIR logical
ZMIR logical
AXIBFC logical
NOSPCT logical
FINE3D logical
VFNORM logical
NUMRAY integer
Plasma:
NEPAR(1:5)real
T0PAR(1:4)real
PHIPAR(1:6)real
- Variables
- dir is E,W,N,S,L,H for East, West, North, South, High or Low. If the line is missing, or the setting is NO there is no surface chemistry on that face.
- Batch_factor is the batch enhancement factor for that face.
- dir is as above
- status can be one of EARTHED, ELECTRODE or INSULATED.
- dir is as above
- status can be one of EARTHED, ELECTRODE, INSULATED or USER.
- Batch_factor is the batch enhancement factor for that face.
- Species data file (SPECIDAT)
- Chemistry data file (CHEMIDAT)
- Data file entry format
The chemistry data file, CHEMIDAT, contains all the information required to calculate the reaction source terms in the species mass fraction equations. The reactions are identified by a character (G for gas-phase or S for surface) and an integer (the reaction number). Each reaction is allocated a number of lines with these identifiers at the start of the first line. Comments can be included between entries using ! as the first character on such lines.
A sample from the file indicates the nature of the entries:
!-----G5-------------
G 10 3 SI2H6<=>SI2H4+H2 2.100E+10 0.000E+00 1.930E+05 5.000E-01
3 SI2H6 -1 SI2H4 1 H2 1
0
0 0
!----------------------------------------------------------------------------------
! Kleijn model: surface reactions
! the parameters for the sticking coefficients of S1 and S3 were taken from
! C.R. Kleijn, "Transport Phenomena in CVD reactors", PhD thesis, TU Delft, 1991
!----------------------------------------------------------------------------------
!-----S1-------------
S 11 5 SIH4<=>SI+2H2 0.000E+00 0.000E+00 0.000E+00 0.000E+00
3 SIH4 -1 SI 1 H2 2
0
-101 7 8.000E-02-7.500E+03 1.200E-06 4.990E+04 1.700E+04-1.538E+05
6.100E-01
SI 2.300E+03 2.809E-02 2.252E-05
The format of the first line is A1,I3,I2,1X,A32,1P4E10.3. The entries are, in order: the character identifier, the reaction number, the number of additional lines in the entry, the reaction in symbols (for information only) and the four Arrhenius parameters.
The format of the second line is 1X,I2,4(1X,A16,I3). The entries are, in order: the number of species in the reaction (maximum 8), the name and stoichiometric coefficient for each species (up to 4). If more than 4 species are involved a continuation line should be inserted with format 3X,4(1X,A16,I3). There is no need to fill either line if the number of species does not require it.
The format of the third line is 1X,2I4. This line will be used for third-body enhancement factors but is currently not used: zeros should be inserted.
The format of the fourth line is 1X,2I4,1P6E10.3. The entries are, in order: IRT (defining the reaction type - see below), the number of additional parameters (maximum 90) and the values of those parameters (up to 6). If more than 6 additional parameters are used a continuation line should be inserted with format 9X,1P6E10.3; again there is no need to fill the line.
For surface reactions there is another line; its format is 1X,A16,1P3E10.3. The entries are the name of the deposited species, followed by its density and molecular mass; the final entry is the surface concentration of vacant sites (mol/m2) which is only required if adsorbed species are being used. Users intending to use adsorbed species in their surface chemistry modelling should consult CHAM for advice before proceeding.
- Reaction rate calculations
The formulation to be used for a reaction is specified in the data file by means of the integer IRT (on line 4). A negative value indicates the use of user-coding; this should be inserted in functions GAF and GAR, for forward and reverse gas-phase chemistry, or GFLUX, for surface chemistry (both these routines are in the 'cvdchm' file). Non-negative values select built-in options as described in Appendix 1; only the specification of the parameters is described here. All the options are selected by the value of IRT, modulo 100; that leaves the 'hundreds column' free to be used for passing additional information about the reaction
(a) Gas-phase reactions
The built-in reaction are:
0 - extended Arrhenius (4 parameters)
1 - Lindemann form (6 parameters)
2 - Troe form (9 parameters)
3 - Troe form (10 parameters).
The Arrhenius parameters appear on line 1 of the data file entry in the order Ak, bk, Ek, ck, with ck only used in the extended Arrhenius option; for all other choices the first three parameters are Ak,inf, bk,inf, Ek,inf. The line of additional parameters contains as many as necessary of, in order,
Ak,0, bk,0, Ek,0, ak, T1,k, T3,k, T2,k;
the Lindemann form requires only the first three, the 9 parameter Troe form requires the first six, and the 10 parameter Troe form requires all the parameters.
The 'hundreds column' of IRT is used to determine whether gas or plasma temperature should be used in calculating reaction rates, and whether the reaction is reversible; a value of 5 indicates plasma temperature, 0 represents a reversible reaction and 1 an irreversible one.
As examples, IRT=3 activates the 9 parameter Troe form using gas temperature, IRT=501 the Lindemann form using plasma temperature and IRT=-2 user-coding.
(b) Surface reactions
The built-in options are:
0 - sticking coefficient formulation
1 - Langmuir-Hinshelwood formulation.
The Arrhenius parameters in the sticking coefficient appear on line 1 of the data file entry in the order Ak, bk, Ek (the fourth entry is not used, but should be present); if Ak=0, a sticking coefficient of unity is used.
User-coding has been included for the Kleijn model for SiH4 and Si2H6, the Roenigk and Jensen model for SiH4 and the Kleijn et al model for tungsten production from WF6 and H2; these correspond to 'units column' values of 1 to 4 respectively. A plasma model has been added for demonstration purposes only; it uses a value of 5.
The 'hundreds column' of IRT is used to determine the rate-determining mechanism for the surface deposition:
1 - kinetics limited
2 - diffusion limited
3 - diffusion/kinetics limited
4 - no deposition.
As examples, IRT=100 activates the sticking coefficient form with the reaction limited by kinetics and IRT=-301 the Kleijn model for SiH4 deposition, taking account of both kinetics and species diffusion to the surface in determining the reaction rate.
The use of diffusion limitation should be used with care as it can result in significant underestimation of reaction rates. It can be necessary when the surface reaction formulation specifies a sub-linear dependence on reactants which would become unrealistic at low concentration levels.
- Data file entry format
- Optical properties data file (OPTICDAT)
- Transport and thermodynamics data files (TRANSDAT and
THERMDAT)
These data files are in the same format as that adopted by the CHEMKIN package. A sample from THERMDAT follows:
Additional parameters are written for the BLOCKAGE and PLATE objects to capture the surface chemistry and plasma boundary condition settings made. The radiation surfaces are captured as RAD_SURF objects.
Blockage
The extra lines written for a BLOCKAGE object for surface chemistry are:
> OBJ, SUF_REAC_dir, YES
> OBJ, BATCH_F_dir, batch_factorWhere:
The extra settings made for plasma boundary conditions are:
> OBJ, PLASMA-AC_dir, status
Where:
If the line is missing, INSULATED is assumed.
> OBJ, PLASMA-DC_dir, status
Where:
If the line is missing, INSULATED is assumed.
If USER, there will be another line
> OBJ, PLASMA-V_dir, Coefficient, Value
Plate
The extra lines written for an external PLATE object for surface chemistry are:
> OBJ, SUF_REAC, YES
> OBJ, BATCH_F, batch_factorWhere:
If the line is missing, NO surface chemistry is assumed.
The extra settings made for plasma boundary conditions are:
> OBJ, PLASMA-AC(DC), status
Where status can be one of EARTHED, ELECTRODE or INSULATED (or USER for PLASMA-DC). If the line is missing, INSULATED is assumed. For internal plates, _L is appended to SUF_REAC, BATCH_FAC and PLASMA for Low-side conditions, and _H for High-side conditions.
Radiation Surface
The attributes for a RAD_SURF object are:
Object type:
> OBJ, TYPE, RAD_SURFManually or automatically generated:
> OBJ, SETTING, AUTOIf the object has been generated automatically, the user is prevented from changing the object side and whether the surface is ‘internal’ or ‘external’. ‘Internal’ surfaces lie on participating solids. All other surfaces are ‘external’.
External linkage coefficient:
> OBJ, TEMP_COEF, temp_coefTemp_coef can be GRND1 or GRND2
External heat flux (for ‘external’ object only):
> OBJ, HEATFLUX, QsurfQsurf is the external heatflux in W/m3
Or external temperature (for ‘external’ object only):
> OBJ, SURF_TEMP, TsurfAlternatively, Tsurf the external surface temperature in K.
Material for optical properties:
> OBJ, OPTIC_INDX, imatimat is the index of the material in the OPTICDAT file
Temperature for evaluating optical properties:
> OBJ, OPTIC_TEMP, ToptTopt is the temperature in K at which the optical properties of the surface will be evaluated
Or thickness of optical coating:
> OBJ, COAT_THICK, coatAlternatively, coat is the thickness of a semi-transparent coating in m.
Object side:
> OBJ, SIDE, HIGHIf the line is absent LOW is assumed. This is only required for RAD_SURF objects not lying on the domain boundaries.
The species data file is used to list the species that are available, up to a limit of 300. A sample demonstrates the format:
140 G SIH4
141 G SIH3
142 G SIH2
The identification numbers must run consecutively, starting at 1. The following letter indicates a gas species (G) or a surface adsorbed species (A); the species name which follows must be exactly the same as that used in other data files. Adsorbed species can be used if appropriate surface chemistry is selected.
The optical data file, OPTICDAT, defines the surface properties for use in the surface-to-surface radiation model. The bulk of the information it contains for each material relates to the temperature dependence of the refractive index and absorption coefficient in each of the 60 spectral intervals. If the spectral variation is not activated only the first line of each entry is used (although the others must be present).
The first entry is the material label; this should be between 100 and 199 for opaque solids, and between 200 and 299 for semi-transparent materials. The next two entries are the constant values of emissivity and reflectivity respectively. The rest of the line relates to the Monte Carlo model which is not part of this release.
A sample entry follows:
# Al=ALUMINIUM Data from E.D.Palik
100 0.1 0.9 0 0.000 0.000 0.000 0.000
1 99.000 0.000 0.000 0.000 0.000E+00 0.208E+03 0.000E+00 0.000E+00 1
1 99.000 0.000 0.000 0.000 0.000E+00 0.223E+03 0.000E+00 0.000E+00 2
1 99.000 0.000 0.000 0.000 0.000E+00 0.204E+03 0.000E+00 0.000E+00 3
1 99.000 0.000 0.000 0.000 0.000E+00 0.205E+03 0.000E+00 0.000E+00 4
1 99.002 0.000 0.000 0.000 0.000E+00 0.212E+03 0.000E+00 0.000E+00 5
1 98.991 0.000 0.000 0.000 0.000E+00 0.205E+03 0.000E+00 0.000E+00 6
. . . . . . . . . .
. . . . . . . . . .
4 0.120 0.000 0.000 0.000 0.000E+00 0.230E+01 0.000E+00 0.000E+00 55
4 0.095 0.000 0.000 0.000 0.000E+00 0.199E+01 0.000E+00 0.000E+00 56
4 0.075 0.000 0.000 0.000 0.000E+00 0.172E+01 0.000E+00 0.000E+00 57
4 0.059 0.000 0.000 0.000 0.000E+00 0.145E+01 0.000E+00 0.000E+00 58
4 0.046 0.000 0.000 0.000 0.000E+00 0.121E+01 0.000E+00 0.000E+00 59
4 0.035 0.000 0.000 0.000 0.000E+00 0.968E+00 0.000E+00 0.000E+00 60
For each spectral interval (numbered on the right hand end of the line) there are 8 real parameters. These are used to define the complex refractive index of the material, n + ik; the first 4 parameters relate to the real part and the second 4 to the imaginary part, as follows:
n = n0 + n1T* + n2T*2 + n3T*3
k = exp(k0T*) (k1 + k2T* + k3T*2)
where T* = (T-300)/1000. The spectral intervals cover the range of wavenumbers from 102 to 105 cm-1. The viewfactor approach splits the spectral range into a smaller number of bands within which the properties can be regarded as approximately constant; the numbers in the left hand column define appropriate bands for each material,
Users can add to the optical data file in order to simulate the behaviour of materials that are not already included, bearing in mind the compatibility with PRPS values. The version of the data file provided with the code demonstrates this: substances such as aluminium, which are included in the built-in PHOENICS property table, are given the same material number in the optical properties data file; others, such as tungsten, which are not, have been allocated higher numbers and their material properties should be defined in the 'q1' file or by adding to the 'props' file.
SIH4 121386SI 1H 4 G 0300.00 4000.00 1000.00 1
0.06893873E+02 0.04030500E-01-0.04183314E-05-0.02291394E-08 0.04384766E-12 2
0.11070374E+04-0.01749116E+03 0.02475166E+02 0.09003721E-01 0.02185394E-04 3
-0.02681423E-07-0.06621080E-11 0.02925488E+05 0.07751014E+02 4
The first line contains the species name, a number which is not used and then the elemental composition of the species; the numbers at the end of the line are the low, high and transition temperatures respectively. The next three lines contain polynomial coefficients used to determine enthalpy, entropy and specific heat as a function of temperature within the two ranges thus defined. The first 7 parameters apply to the high temperature range and the remaining parameters to the low temperature range. An exception is made for adsorbed species in the Arora & Pollard surface chemistry model. In that case the high data are used instead for adsorption data and the high and transition temperatures are set equal.
An extract from TRANSDAT follows:
SIH4 2 207.6 4.084 0.000 0.000 1.000 ! mec
SIH3 2 170.3 3.943 0.000 0.000 1.000 ! mec
SIH2 2 133.1 3.803 0.000 0.000 1.000 ! mec
Here only the numbers in the third and fourth columns are used: they are the Lennard-Jones parameters e/k and s respectively for the species named in the first column.
As with the other data files, the user is expected to add entries for missing materials in the standard format. Note that the names of any species must be consistent in all data files in which they appear.
Three input 'q1' files are presented in Appendix 2. These have been constructed for the simulation of an ASM single-wafer reactor for silicon deposition, a tungsten reactor and a parallel plate plasma reactor, also for silicon deposition. All geometries are axisymmetric and laminar, steady-state conditions are assumed. Although these examples are based on real cases it should be noted that some of the parameters are selected for demonstration purposes. The 'q1' files are also provided in the code installation as library cases.
In the first case the inlet flow is 88.7% nitrogen and 11.3% silane by mass, at a temperature of 290K. The susceptor is held at a temperature of 900K and all outside walls are cooled to a temperature of 290K. Five gas-phase and five surface reactions are activated.
The tungsten reactor has a showerhead inlet, the details of which are not modelled; instead a uniform inlet is specified. The inlet flow is 6% hydrogen and 94% tungsten hexafluoride by mass, corresponding to approximately 90% hydrogen by volume. Gas-phase chemistry is less important in this case and only a single surface reaction is used in the simulation. The showerhead structure is a very significant part of the problem because of its thermal effect: it is heated radiatively by the hot tungsten susceptor (733K) and in turn passes heat to the gas by convection and to the outer wall by radiation.
The plasma enhanced reactor is for silicon deposition from pure silane. Three gas-phase and two surface reactions are included, all dependent on the state of the plasma. The gas-phase reaction rates all exhibit dependence on ionisation rate while the surface reaction rates all have an Arrhenius formulation, based on surface temperature, multiplied by the ion flux to the surface. (The reactions defined in CHEMIDAT for this simulation are for demonstration purposes only.)