The role of thermal energy accommodation and atomic recombination probabilities in low pressure oxygen plasmas

Surface interaction probabilities are critical parameters that determine the behaviour of low pressure plasmas and so are crucial input parameters for plasma simulations that play a key role in determining their accuracy. However, these parameters are difficult to estimate without in situ measurements. In this work, the role of two prominent surface interaction probabilities, the atomic oxygen recombination coefficient γO and the thermal energy accommodation coefficient αE in determining the plasma properties of low pressure inductively coupled oxygen plasmas are investigated using two-dimensional fluid-kinetic simulations. These plasmas are the type used for semiconductor processing. It was found that αE plays a crucial role in determining the neutral gas temperature and neutral gas density. Through this dependency, the value of αE also determines a range of other plasma properties such as the atomic oxygen density, the plasma potential, the electron temperature, and ion bombardment energy and neutral-to-ion flux ratio at the wafer holder. The main role of γO is in determining the atomic oxygen density and flux to the wafer holder along with the neutral-to-ion flux ratio. It was found that the plasma properties are most sensitive to each coefficient when the value of the coefficient is small causing the losses of atomic oxygen and thermal energy to be surface interaction limited rather than transport limited.


Introduction
Low temperature plasmas produced in low pressure gases are widely used for industrial applications, prominent examples include the etching and deposition of nano-scale structures in the semi-conductor industry [1][2][3], electric propulsion of spacecraft [4,5] and as potential negative ion sources for neutral beam heating of fusion plasmas [6,7]. For such applications plasma simulations are often used to aid in reactor and process design and to better understand fundamental processes [6,[8][9][10][11][12][13][14][15]. These simulations require a wide range of fundamental input data such as electron impact cross sections, heavy particle rate coefficients and surface interaction probabilities [16]. The accuracy of any plasma simulation is inherently linked to any uncertainty in these fundamental input data [17]. Previous investigations have studied the sensitivity of plasma simulations to the uncertainty in gas phase reaction cross sections and rate coefficients [17] and the consequences of including or neglecting various collisional processes [18,19]. Such investigations yield valuable information to assess how the accuracy of plasma simulations may be improved and the number of reactions can be minimised in order to simplify chemical reaction mechanisms.
There have also been a number of experimental and computational studies that have discussed the role of surface recombination probabilities for atomic neutral species [20][21][22][23] and surface quenching probabilities for excited neutrals [24][25][26] in defining the plasma dynamics. In this context surface recombination refers to, for example, an oxygen atom recombining with another oxygen atom on a surface, and returning to the plasma as an oxygen molecule. It has been shown that changes in the value of these coefficients can lead to large variations in not only the reactive neutral species densities, to which they are directly linked, but also to the densities of negative ions and the overall electron heating dynamics. As a result the predictive capability of plasma simulations where reactive atomic and radical species or excited neutrals are present is strongly affected by the accuracy of these surface interaction probabilities. The situation is complicated further by the fact that these parameters are prone to vary with specific plasma conditions and reliable predictive formulae for them do not yet exist.
Experimental investigations have been carried out to measure surface recombination probabilities for various atomic or excited species under different plasma and surface conditions [20,[27][28][29][30]. However, the often complex dependencies of these probabilities on both the plasma and surface condition means that deriving a general scaling law that is applicable to any atom/molecule in contact with any surface exposed to any type of plasma is, at best, challenging. As a result these probabilities still represent a major uncertainty in plasma modelling and thus constrain the parameter space over which plasma simulations can be expected to yield accurate results.
Another surface interaction probability that has been less studied in the context of non-equilibrium plasmas is the thermal energy accommodation coefficient, often denoted as α E . This coefficient represents the fraction of atoms or molecules which attain thermal equilibrium with a surface after a collision with it and is often defined in the form [31]: where T i and T r are the temperatures of particles incident on the surface and reflected from the surface respectively. T s represents the surface temperature, or equivalently the temperature of particles reflected from the surface that attain thermal equilibrium with the surface. The value of the thermal energy accommodation coefficient has been determined for a range of gas-surface combinations experimentally or through numerical simulations [32][33][34][35][36][37][38][39][40][41]. Many of these investigations have focused on noble gases or air interacting with metallic surfaces [32,34,35,38].
However, fewer studies have dealt with reactive gases, or interactions with non-metallic surfaces. This is significant for the application of these coefficients in plasma modelling as reactive gases are used in the majority of industrial applications and non-metallic surfaces such as inert metal-oxides commonly make up a large part of the overall surface in contact with the plasma.
Given that the value of α E determines the rate at which gas phase atoms or molecules exchange energy with surfaces it will play an important role in determining the neutral gas temperature in non-thermal plasma sources and have a significant influence on the dynamics of the plasma as a whole. Several physical phenomena in non-thermal plasmas may result in variations of the value of α E from those measured or calculated for 'clean' surfaces, or cause variations from process to process in the same reactor. For example, it has been proposed, based on the results of molecular dynamics simulations, that the value of α E depends on the presence of adsorbed atoms or molecules on surfaces [38]. In non-thermal plasmas the coverage of surfaces by adsorbed atoms or molecules is known to be affected by several plasma parameters such as the flux of atomic neutral species, ion bombardment and surface temperature [42][43][44]. All of these parameters are susceptible to change from process to process in a given plasma reactor either intentionally as a mechanism to improve process outcomes, or unintentionally due to process drifts [12]. As such, it is reasonable to expect some change in the value of α E under varying plasma conditions, however, it is difficult to quantify how significant this change may be.
Given the uncertainty in the values of thermal accommodation and atomic species recombination coefficients it is necessary to understand, quantitatively, how the variation of these parameters affects the properties of a given plasma. While much is known in a qualitative sense about how these coefficients affect gas temperatures and atomic species densities, it is important to consider the ranges over which these parameters have the most significant effect on the overall plasma properties, in particular those quantities relevant for industrial applications. Quantities of specific relevance in this regard are reactive neutral and ion fluxes to surfaces and the ratio of the two quantities [45]. In addition, the ion bombardment energy at the substrate is known to be an important factor in many industrial plasma processes [13,45,46]. In this context a variety of different schemes exist to control these quantities in low pressure plasmas, such as substrate biasing [47][48][49][50] and voltage waveform tailoring [11,26,47,49,[51][52][53][54][55][56]. However, the known importance of recombination and accommodation coefficients in determining atomic species densities and gas temperatures suggests that they have the potential to influence all of these properties of a given plasma which may affect the range over which they can be controlled in industrial applications.
This study seeks to quantify the influence of varying the atomic oxygen recombination coefficient γ O and thermal energy accommodation coefficient α E on the plasma properties in a low pressure inductively coupled oxygen plasma through two-dimensional numerical simulations inclusive of charged and neutral particle dynamics. It has been found that the thermal energy accommodation coefficient has a strong influence on the plasma properties through its role in determining the neutral gas temperature and consequently the neutral gas density. The role of the atomic oxygen recombination coefficient is primarily in determining the atomic oxygen density and its flux at the wafer holder as well as the reactive neutral-to-ion flux ratio to the wafer holder. The greatest sensitivity of the plasma to α E and γ O occurs at low values of each coefficient when the losses of atomic oxygen and thermal energy are limited by the surface interaction probability as opposed to the rate of diffusion to surfaces. The numerical simulations used to carry out this investigation are described in section 2. The results of this investigation are discussed in detail in section 3 and our conclusions are presented in section 4.

Description of the model
In this work a two-dimensional, modular, fluid-kinetic simulation code, the hybrid plasma equipment model (HPEM) has been used to model a planar inductively coupled plasma reactor. The geometry of the simulation domain is shown in figure 1 and is based on the experimental reactor presented in [57,58]. The simulation domain consists of a cylindrically symmetric chamber 27.5 cm in radius with a gap of 10 cm between the wafer holder and the alumina window separating the gas from a 4-turn inductive coil which is used to drive the plasma. The HPEM code has previously been described in detail elsewhere [9] so only a brief overview of the relevant details is given here. The modular structure of the HPEM allows for the relevant physics of a given system to be addressed in different ways using various modules. In this work the modules used are the fluid kinetics Poisson module (FKPM), the electron Monte Carlo simulation (eMCS), the electromagnetics module (EMM) and the plasma chemistry Monte Carlo module (PCMCM).
The FKPM forms the basis of the simulation and solves the continuity, momentum and energy equations for all heavy particle species separately in combination with Poisson's equation for the electric potential. Electron momentum is solved for using the drift-diffusion approximation. The accommodation coefficient α E and the value of γ O are assumed the same for all surfaces. The EMM solves for the electric and magnetic fields associated with the current in the inductive coil which is adjusted to deliver the specified input power. The eMCS simulates the interaction of electron pseudo-particles with electric fields obtained from the EMM and FKPM using a Monte Carlo algorithm which is described in detail in [59]. This algorithm yields the electron energy distribution function (EEDF) at all points in the simulation domain along with the corresponding electron impact rate coefficients and electron transport properties which are used as input for the FKPM in the following iteration.
After the simulation has reached convergence the PCMCM is used to simulate ion and neutral energy distribution functions to the wafer holder, which is marked in the simulation geometry shown in figure 1. The PCMCM launches ion and neutral pseudo-particles from locations in the plasma weighted by their source functions obtained from the FKPM and integrates the trajectories of these pseudo-particles through the time and space dependent electric fields obtained from the FKPM. All particles reaching the wafer holder are recorded and as such the IEDFs and ion/neutral fluxes presented in section 3 are averages over the entire surface of the wafer holder.
The thermal energy accommodation coefficient α E is used to define a 'jump' in the temperature of a given species in the computational cell adjacent to the reactor wall. Following Lofthouse et al [60] the temperature gradient adjacent to a surface is given by: where T g is the temperature of the species, γ is the ratio of specific heats (i.e specific heat at constant pressure and specific heat at constant volume), P r is the Prandtl number and λ is the collisional mean free path for momentum transfer.
Requiring that the heat flux from the plasma to the wall is constant this condition is implemented by defining an effective thermal conductivity in the computational cell adjacent to the wall such that: where κ is the unmodified thermal conductivity and Δx is the computational mesh spacing. The reaction scheme used in the simulation includes electrons, ground state molecular oxygen (O 2 (X S - The full list of reactions is given in table 1. This reaction scheme has been designed to incorporate the dominant reactions between the main charged and neutral species expected to be present under the conditions studied, based on previous  No. Reaction 1.81×10 −18 T 0 3.8 exp(700/T g ) [86] works incorporating larger reaction schemes [19,61]. The exothermicity of each reaction for use in the gas temperature calculation is calculated based on the enthalpy of formation of each species and for electron impact dissociation from the measurements of Cosby [62]. The metastable molecules ) are assumed to be quenched at surfaces with probabilities of 5×10 −5 and 2×10 −2 respectively. These coefficients are generally poorly known and as such the values used here represent estimates in the same range as those used in previous works [24-26, 61, 63].
In order to gain insight into how the change of atomic oxygen recombination γ O and thermal energy accommodation coefficients α E affect the discharge dynamics the values of the two coefficients are varied over ranges which may be expected for the conditions investigated. Here, the same value of γ O has been applied at the surface for both O( 3 P) and O( 1 D). Due to a lack of experimental data, and in order to simplify analysis of the results presented, no quenching of O( 1 D) at the surface causing it to return as O( 3 P) is considered. Values of γ O reported in the literature range from ≈ 0.01 [28] to ≈ 0.5 [29] depending on the surface material, in the pressure range considered here. Thus, the value of γ O has been varied over this range in order to cover the maximum variation of the plasma parameters resulting from differences in its value for different wall materials, or variations from process to process where the wall material is susceptible to change through deposition or etching.
As discussed previously, relatively few studies have been devoted to the measurement of α E in low temperature plasmas. Given that α E is known to depend on a range of different properties of both the gas and the surface, it is difficult to predict in advance which values should be used for any given system. In this context, α E is varied over a wide range from 0.1 to 0.95. It should be noted that this range is not necessarily claimed to be realistic for an oxygen plasma. It is chosen here to emphasise the maximum range of variation of the other plasma parameters possible under variation of α E in order to assess its importance. However, it is expected that, depending on the plasma-surface combination in question, there exist plasma processes where the effective value of α E may lie at almost any point in the range considered here. Figure 2 gives examples of the spatial profiles of (a) the ground state atomic oxygen density and (b) the average gas

No.
Reaction 3.90×10 −17 [85] a Units: Rate coefficients in cm 3 /s; T 0 =T g /300; T g has units K, T e has units eV. b f(ò) indicates that the rate coefficients are obtained from the EEDF calculated via the eMCS. Electron-impact excitation of O 2 into rotationally, vibrationally and electronically excited states and of O( 3 P) into electronically excited states, according to the cross sections given in [64,65,70] is accounted for in the eMCS to properly simulate electron energy losses. c Assumed to be the same as for O( 3 P). d Superelastic cross section obtained by detailed balance from the corresponding excitation process. e Cross section estimated by shifting and scaling the corresponding cross section for the ground state by the excitation threshold of the metastable. f Cross section estimated by shifting and scaling the corresponding cross section for O 2 (a 1 Δ g ). g In reality this reaction proceeds through the production of O 3 , which is neglected in this model for simplicity. h Assumed the same as for O 2 (a D g 1 ) state. i Lower rate assumed for collisions of + O 2 with metastable molecules based on observations of [78] for the case of O + . temperature as a function of radial and axial dimensions. The average gas temperature is defined as a mole fraction weighted sum of the temperatures of the individual neutral species in the plasma. The parameters used here represent the base case where the pressure p = 10 mTorr, and the power P = 500 W. The gas flow rate is set at 20 sccm which leads to an average gas residence time of around 1 s. The temperature of all surfaces in contact with the plasma is held constant at 325 K. Here, the accommodation coefficient α E =0.5 and the atomic oxygen recombination coefficient γ O =0.05. Under these conditions the dominant power deposition process is through electron heating directly below the coil. This heating allows for electrons to gain sufficient energy to dissociate molecular oxygen, via reactions R18-R23 in table 1, producing two oxygen atoms with a total of 0.5-3 eV of thermal energy depending upon the particular dissociation pathway [62]. This thermal energy is produced through a Franck-Condon process where electron impact produces a dissociative excited state with excess potential energy being dissipated into the translational modes of the fragments. The rate of this process is maximum just under the coil at a radial location of around 9 cm from the centre of the reactor. However, the mean free path of oxygen atoms under the low pressure conditions considered here is sufficiently long that these atoms easily diffuse from their point of creation, resulting in a relatively flat atomic oxygen density profile as a function of both radial and axial dimensions. Similarly, the dominant process transferring electron energy into gas heating is through the dissociation process which produces atomic oxygen with excess thermal energy. Overall, the gas temperature takes on a profile typical of conditions where thermal diffusion to the walls is the dominant energy loss process. The density of atomic oxygen increases towards the walls of the reactor under these conditions because the thermal energy accommodation coefficient α E is much larger than atomic oxygen recombination coefficient γ O . As a result, the gas temperature decreases towards the walls, which act as sinks for thermal energy. Losses of atomic oxygen at the walls are proportionally lower as γ O is small. The result is that the atomic oxygen density profile does not decrease towards the walls and instead follows approximately the total neutral gas density which increases near the walls due to the decrease in gas temperature. This is facilitated by the essentially constant value of the gas pressure, which varies by only a few percent over the whole reactor. In general, the radial and axial profiles of many plasma parameters change as γ O and α E are varied. However, for the low pressures considered here these changes in spatial profile are usually small. For example, the relative variation in the value of each plasma parameter at the centre of the reactor with the two surface coefficients is similar to the change in their volume-averaged values. As a result, in the remainder of this article we discuss the values of each parameter at the centre of the reactor and species fluxes to the centre of the wafer holder to illustrate major trends as a function of both surface coefficients.
The ground state atomic oxygen density and gas temperature at the radial and axial centre of the reactor are shown as a function of the atomic oxygen recombination coefficient figure 3(a). The atomic oxygen density decreases from 8.3×10 13 cm −3 for γ O =0.01 to around 5.9×10 12 cm −3 when γ O =0.5. The gas temperature remains comparatively constant over the same range, increasing from 860 K to 990 K. The decrease of the atomic oxygen density with increasing γ O is a result of the increasing conversion of atomic oxygen to molecular oxygen at the surface as γ O increases.
Under the low pressure conditions considered here, where the mean free path between gas phase collisions is long, the diffusion of oxygen atoms to the reactor walls is rapid and as a result wall losses are the dominant destruction pathway of atomic oxygen since the rate of gas-phase reactions is low. As such, even relatively small changes in γ O can lead to significant changes in the atomic oxygen density as shown here. This is particularly true at low values of γ O where the density of atomic oxygen changes rapidly with varying γ O . As γ O is increased the loss of atomic oxygen to the walls transitions from being limited by the rate of surface recombination at the walls at low γ O to being limited by the rate of diffusion towards the walls at high γ O . This trend is approximately in accordance with the formula developed by Chantry [87] for the surface loss rate of active species in low pressure gases and the conclusions of other previous work [29] which suggest that atomic species densities scale with 1/ γ O under low pressure conditions. The small increase in the gas temperature with increasing γ O is a result of the increasing density of molecular oxygen, both in the ground and metastable states, as the contribution of atomic oxygen to the overall gas density decreases. This results in a higher rate of dissociation per unit volume and thus greater thermal energy deposition through the production of oxygen atoms with high translation energies by the Franck Condon effect, as discussed earlier.
The ground state atomic oxygen density and gas temperature at the radial and axial centre of the simulation as a function of the thermal energy accommodation coefficient α E for γ O =0.05 are shown in figure 3(b). The atomic oxygen density is also dependent on the value of α E , increasing from 1.5×10 13 cm −3 to around 5.1×10 13 cm −3 as α E is increased from 0.1 to 0.95. Under this range the degree of dissociation (i.e. the density of atomic oxygen divided by the total gas density) decreases slightly as α E is increased from 0.1 to 0.3 and remains approximately constant for the remainder of the range. The initial decrease in dissociation degree with increasing α E is a result of the decreasing electron temperature as shown in figure 4(b). This leads to a decreased rate coefficient for electron impact dissociation of molecular oxygen, thus decreasing the relative density of atomic oxygen. Overall, despite the decrease in dissociation degree with increasing α E between 0.1 and 0.3, the atomic oxygen density increases with α E as a result of higher total gas densities due to the lower gas temperatures.
As expected, the gas temperature is more significantly affected by changes in the thermal energy accommodation coefficient α E than changes in γ O , scaling approximately as 1/α E from 745 K at α E =0.95 to 3065 K at α E =0.1. This relationship is similar to that derived by Chantry [87] for the scaling of active species loss rates with their surface loss coefficients and reflects the fact that the dominant loss of thermal energy from the neutral species in the plasma is via interactions with the wall. As the value of α E decreases the loss of thermal energy from the neutral species transitions from a diffusion limited regime to a surface interaction limited regime. This dependence is analogous to the dependence of the atomic oxygen density on γ O .
The change in the overall gas density with varying α E has consequences for other plasma parameters such as the plasma potential and electron temperature at the radial and axial centre of the simulation, as shown in figure 4(b). Both quantities decrease in a similar manner to the gas temperature with increasing α E . The mechanism behind these trends is an increase in the electron collision frequency caused by increasing gas density at higher values of α E . Higher electron collision frequencies lead to decreased electron losses to the walls via diffusion. This in turn means that the required positive ion losses to the walls, and consequently the ionization rate, are also decreased. As a result, the electron temperature decreases in order to fulfil the requirement for a lower ionization rate. This in turn leads to a decrease in the plasma potential which is typically proportional to the electron temperature. The same trends are produced by global models where the gas temperature is fixed as a constant and the gas pressure is increased, which causes a similar change in the overall gas density as considered here.
The values of the plasma potential and the electron temperature are affected to a smaller degree by the value of the atomic oxygen recombination coefficient γ O for constant α E , shown in figure 4(a). The slight decrease in plasma potential and electron temperature with increasing γ O is due to a small increase in the electron collision frequency with increasing γ O as the atomic oxygen density decreases and molecular oxygen becomes the dominant component of the neutral gas. This is a result of the higher total electron impact cross section for molecular oxygen, in the ground or metastable states, than that of atomic oxygen due to the opportunity to excite vibrational and rotational modes.
The densities of charged species at the radial and axial centre of the reactor as a function of γ O and α E are shown in figures 5(a) and (b) respectively. As γ O increases the electron density decreases slightly, while the O − density increases. As a result, the electronegativity (i.e. the density of negative ions ) increase. This results in increased negative ion formation through electron attachment which also acts as a loss channel for electrons. However, the total electron density will be sustained by requiring a constant power deposition.
By As the thermal energy accommodation coefficient α E increases from 0.1 to 0.2 the electron density decreases. With a constant power deposition, the electron density is determined by the volume integral of n e R E (T e )N, where R E is the rate coefficient for electron power loss and N is the gas density. With the increase in gas density and moderate change in R E (T e ), the electron density decreases. Over the same range of α E the density of O − increases as a result of the increased electron attachment rate as the density of molecular oxygen, in both the ground and metastable states, increases. As α E increases above 0.2 the electron collision frequency continues to increase causing electrons to be more effectively confined to the centre of the reactor which results in an increase in the central electron density.
Conversely, for α E 0.2 the density of O − remains approximately constant as the rate of electron attachment to oxygen molecules increases at a similar rate to the rate of electron detachment from O − with increasing gas density. The density of + O 2 ions at the centre of the reactor increases with α E due to the increasing ion-neutral collision frequency as the gas density is increased. This results in enhanced confinement of ions at the centre of the reactor. The density of O + ions at α E 0.2 increases slowly with increasing α E , behaving similarly to + O 2 over the same range. However, from α E =0.1 to α E =0.2 the density of O + ions decreases due to a decrease in the dissociation degree, and a decrease in the ionization rate coefficient of atomic oxygen which occurs due to the decreasing electron temperature. The ionization rate of atomic oxygen, producing O + , is more sensitive to changes in electron temperature than the ionization rate of molecular oxygen, producing + O 2 , as a result of the higher electron impact ionization threshold for atomic oxygen. This means that the density of O + is more significantly affected by the change in electron temperature over this parameter range.
For surface processing applications the most relevant plasma parameters are species fluxes to surfaces where wafers are placed and treated, rather than plasma bulk properties in themselves. Figure 6 shows the ion energy distribution functions (IEDFs) for + O 2 and + O , as calculated by the PCMCM, for several combinations of α E and γ O . Figure 6(a) shows the IEDFs for γ O =0.01 and 0.5 (α E =0.5). The IEDFs shown are integrated over all angles of ion impact. At γ O =0.01 the IEDFs of + O 2 and O + are similar in both the total flux and average energy with the peak in ion flux occurring at an energy approximately equal to the plasma potential. As γ O is increased to 0.5 the total ion flux remains approximately the same, however, the flux of + O 2 becomes greater than that of O + . This is a result of lower atomic oxygen densities in the gas phase leading to decreased production of O + , whose dominant formation pathway is through direct electron impact ionization of atomic oxygen. The small decrease in the energy at the peak in ion flux is as a result of the decrease in plasma potential with increasing γ O shown in figure 4(a). The IEDFs for α E =0.1 exhibit a peak in ion flux at higher energy than those for α E =0.9. This is a result of the higher plasma potential at lower α E , as shown in figure 4(b). The IEDFs at α E =0.1 also exhibit a smaller spread in ion energies. The reason for this is that ions experience fewer collisions with neutral species on transit to the surface due to the lower neutral gas density at low values of α E . In addition, the relative contribution of O + to the total ion flux decreases as α E is increased. This is a result of the lower dissociation degree and electron impact ionization rate for atomic oxygen at higher values of α E discussed previously, leading to decreased O + production and thus lower fluxes of O + .
The total fluxes of + O 2 , O + and ground state atomic oxygen at the wafer holder as a function of γ O for α E =0.5 are shown in figure 6(a). Here, the fluxes of + O 2 and O + are almost equal for γ O =0.01 but diverge as γ O is increased with the flux of + O 2 increasing and the flux of O + decreasing. The trend in the neutral atomic oxygen flux to the wafer holder is similar to that for O + with increasing γ O . These trends are consistent with the variation in the gas phase atomic oxygen density with increasing γ O discussed earlier.  The total fluxes of + O 2 , O + and ground state atomic oxygen at the wafer holder as a function of α E for γ O =0.05 are shown in figure 7(b). As α E increases from 0.1 to 0.3 the fluxes of both ions to the wafer holder decrease and remain approximately constant with increasing α E thereafter. The more significant decrease in the flux of O + compared to that of + O 2 is as a result of the decrease in the rate coefficients for electron impact dissociation and ionization with decreasing electron temperature as α E increases leading to a lower degree of dissociation and O + density. The flux of + O 2 to the wafer holder does not correlate directly with the density of + O 2 at the reactor centre, shown in figure 5(b). This is a result of the changing shape of the + O 2 density profile as the ion-neutral collision frequency increases with increasing α E , meaning that the relationship between the ion density at the centre of the reactor and the ion flux to the wall changes as a function of α E . The atomic oxygen flux to the wafer holder increases with increasing α E consistent with the higher atomic oxygen densities at the reactor centre at larger values of α E due to the higher total gas density.
The reactive neutral-to-ion flux ratio at the wafer holder as a function of γ O for α E =0.5 is shown in figure 8(a) and the reactive neutral-to-ion flux ratio at the wafer holder as a function of α E for γ O =0.05 is shown in figure 8(b). Here, the ion flux is defined as the sum of the fluxes of both + O 2 and O + . In both cases the trend in neutral-to-ion flux ratio approximately follows the respective trend of the ground state atomic oxygen flux. The neutral-to-ion flux ratio is shown to vary over a range of 94 to 8 with increasing γ O and between 10 and 50 with increasing α E .
It has previously been predicted by Brichonet al via molecular dynamics simulations of chlorine plasma etching of silicon that the damaged layer thickness on an etched surface, a critical parameter for atomic layer etching, is strongly dependent on the reactive neutral-to-ion flux ratio in the range considered here (i.e. 10-200) [45]. Specifically, it was found that higher neutral-to-ion flux ratios were favourable to minimise damaged layer thickness. Since the chemistry used in this work is different, the results of the study by Brichonet al cannot be directly applied here, but it is reasonable to expect that similar dependencies are possible when using plasmas produced in oxygen, as considered here, or in other gases, for atomic layer processes. In this context, it is expected that reactor walls with very low atomic species recombination coefficients and high thermal energy accommodation coefficients should be favourable for minimising damaged layer thickness in such applications. These conditions have the additional advantage of offering low ion bombardment energies, which are advantageous for certain processes, due to the dependence of the plasma potential on α E . A potential disadvantage of operating in this parameter regime is that the neutral flux, and neutral-to-ion flux ratio are very sensitive to changes in γ O , meaning that the outcomes of plasma processes may be sensitive to small changes in wall condition.

Conclusions
The dependence of the plasma bulk properties and fluxes to the wafer holder on the atomic oxygen surface recombination coefficient γ O and the thermal energy accommodation coefficient α E in a low pressure inductively coupled oxygen plasma have been investigated using two-dimensional fluidkinetic plasma simulations. It has been shown that the thermal energy accommodation coefficient has significant influence on most of the plasma properties investigated through its role in determining the neutral gas temperature and consequently the neutral gas density. The primary role of the atomic oxygen recombination coefficient is in determining the atomic oxygen density and its flux at the wafer holder in combination with the reactive neutral-to-ion flux ratio at the wafer holder. The greatest sensitivity of the plasma parameters to α E and γ O occurs when the value of each coefficient is low. This is as a result of the transition from diffusion limited losses of atomic oxygen and thermal energy to surface interaction limited losses as the value of each coefficient is decreased. The trends presented offer insight into the roles played by both investigated coefficients in low pressure plasmas used in industry. In particular, the thermal energy accommodation coefficient has been shown to be an important input parameter in plasma simulations which warrants further study, both experimentally and computationally.