GEC ICP Reactor Coupled with the Two-Term Boltzmann Equation

The GEC cell was introduced by NIST (National Institute of Standards and Technology) in order to provide a standardized platform for experimental and modeling studies of discharges in different laboratories. The plasma is sustained via inductive heating. The Reference Cell operates as an inductively-coupled plasma in this model.

This tutorial models the GEC ICP reactor by solving plasma fluid equations fully coupled with the homogeneous and time-independent electron Boltzmann equation in the classical two-term approximation. The approximated Boltzmann equation is solved for in each position of space and is coupled with the fluid equations by way of the electron mean energy. The rate constants of electron impact reactions and the electron transport parameters are obtained by suitable integration of the computed electron energy distribution function over electron scattering cross sections. Simulation results from this model are compared with results from the model GEC ICP Reactor, Argon Chemistry where the electron energy distribution function (EEDF) is assumed Maxwellian.

Figure 1: GEC ICP reactor geometry consisting of a 5 turn copper coil, plasma volume, dielectrics, and wafer with pedestal.

This application requires the Plasma Module and the AC/DC Module.

Model Definition

Fluid model

Inductively coupled discharges typically operate at low pressures (<10 Pa) and high charge density (>1017 m−3). High density plasma sources are popular because low pressure ion bombardment can provide a greater degree of anisotropy on the surface of the wafer.

The electron density and mean electron energy are computed by solving a pair of drift-diffusion equations for the electron density and mean electron energy. Convection of electrons due to fluid motion is neglected. For detailed information on electron transport, see Theory for the Drift Diffusion Interface in the Plasma Module User’s Guide.

The electron source Re and the energy loss due to inelastic collisions Rε are defined later. The electron diffusivity, energy mobility, and energy diffusivity are computed from the electron mobility using:

The source coefficients in the above equations are determined by the plasma chemistry using rate coefficients. Suppose that there are M reactions that contribute to the growth or decay of electron density and P inelastic electron-neutral collisions. In general, P >> M. In the case of rate coefficients, the electron source term is given by:

where xj is the mole fraction of the target species for reaction j, kj is the rate coefficient for reaction j (SI unit: m3/s), and Nn is the total neutral number density (SI unit: 1/m3). The electron energy loss is obtained by summing the collisional energy loss over all reactions:

where Δεj is the energy loss from reaction j (SI unit: V), and kk is the rate coefficient.

For nonelectron species, the following equation is solved for the mass fraction of each species. For detailed information on the transport of the non-electron species, see Theory for the Heavy Species Transport Interface in the Plasma Module User’s Guide.

The electrostatic field is computed using the following equation:

The space charge density ρ is automatically computed based on the plasma chemistry specified in the model using the formula:

For detailed information about electrostatics see Theory for the Electrostatics Interface in the Plasma Module User’s Guide.

For a nonmagnetized, nonpolarized plasma, the induction currents are computed in the frequency domain using the following equation:

The plasma conductivity is specified using the cold plasma approximation

where ne is the electron density, q is the electron charge, me is the electron mass, ω is the angular frequency μr and μi are the real and imaginary parts of the HF mobility. The effective collision frequency νeff and the factor

are obtained in coherence with the Boltzmann equation approach here used.

Boltzmann Equation in the Two-Term Approximation

The EEDF used in this model is obtained from the solution of the homogeneous and time-independent electron Boltzmann equation in the classical two-term approximation

where F0 is the isotropic part of an EEDF constant in time and space that satisfies the following normalization

In this model the Boltzmann equation in the two-term approximation is solved in the high-frequency limit including the effects of electron-electron collisions. The different terms are presented below:

The following definitions apply

Here:

•

γ = (2q/me)1/2 (SI unit: C1/2/kg1/2)

•

me is the electron mass (SI unit: kg)

•

ε = (v/γ)2 is energy (SI unit: V)

•

σε is the total elastic collision cross section (SI unit: m2)

•

σm is the total momentum collision cross section (SI unit: m2)

•

is the normalized total momentum collision cross section (SI unit: m2)

•

q is the electron charge (SI unit: C)

•

ε0 is the permittivity of free space (SI unit: F/m)

•

T is the temperature of the background gas (SI unit: K)

•

kb is the Boltzmann constant (SI unit: J/K)

•

ne is the electron density (SI unit: 1/m3)

•

Nn is the background gas density (SI unit: 1/m3)

•

Λ is the Coulomb logarithm, and

•

λ is a scalar-valued renormalization factor that ensures that the EEDF is normalized to 1as explained above. An ODE is implemented to solve for the value of λ.

•

M is the mass of the target species (SI unit: kg).

The source term, S represents energy loss due to inelastic collisions. Because the energy loss due to an inelastic collision is quantized, the source term is nonlocal in the energy space. The source term can be decomposed into four parts where the following definitions apply:

where xk is the mole fraction of the target species for reaction k, σk is the collision cross section for reaction k, Δεk is the energy loss from collision k, and δ is the delta function at ε = 0.

The mean electron energy is defined by the integral

The external excitation in the Boltzmann equation comes from an electric field. If the Boltzmann equation is to be solved for a given mean electron energy (as in the present model) a Lagrange multiplier is introduced to solve for the reduced electric field, such that the following equation, presented in the weak form, is satisfied:

where tilde denotes test function. For detailed information about electrostatics see Theory for the Boltzmann Equation, Two-Term Approximation Interface in the Plasma Module User’s Guide.

Coupling Between the Fluid Model and the Boltzmann Equation

The fluid model equations and the Boltzmann equation are solved fully coupled. The Boltzmann equation is solved at every point in which the fluid model is solved. This is the space where the mean electron energy, mole fraction of the different heavy species, and electron density (the gas density is assumed constant and as such it does not intervene in the computation of the ionization degree and the reduced angular frequency) are computed. The EEDFs obtained from the solution of the Boltzmann equation are used to compute macroscopic rate constants and transport coefficients intervening in the fluid model equations thus closing the loop.

Electron Rate Coefficients and Transport Parameters

Rate coefficients for electron impact reactions are computed using cross section data and the computed EEDF by the following integral

The reduced electron mobility associated with the DC transport is computed using

and the real and imaginary part of the mobility, used to computed the plasma conductivity, are computed with

boundary conditions

Electrons are lost to the wall due to random motion within a few mean free paths of the wall and gained due to secondary emission effects, resulting in the following boundary condition for the electron flux:

and the electron energy flux:

For the heavy species, ions are lost to the wall due to surface reactions and the fact that the electric field is directed toward the wall:

The walls of the reactor are grounded.

plasma chemistry

Because the physics occurring in an inductively coupled plasma is rather complex, it is always best to start a modeling project with a simple chemical mechanism. Argon is one of the simplest mechanisms to implement at low pressures. The electronically excited states can be lumped into a single species, which results in a chemical mechanism consisting of only 3 species and 7 reactions (electron impact cross-sections are obtained from Ref. 3):

Table 1: table of collisions and reactions modeled.
reaction	formula	type
1	e+Ar=>e+Ar	Elastic	0
2	e+Ar=>e+Ars	Excitation	11.5
3	e+Ars=>e+Ar	Superelastic	-11.5
4	e+Ar=>2e+Ar+	Ionization	15.8
5	e+Ars=>2e+Ar+	Ionization	4.24
6	Ars+Ars=>e+Ar+Ar+	Penning ionization	-
7	Ars+Ar=>Ar+Ar	Metastable quenching	-

Stepwise ionization (reaction 5) can play an important role in sustaining low pressure argon discharges. Excited argon atoms are consumed via superelastic collisions with electrons, quenching with neutral argon atoms, ionization or Penning ionization where two metastable argon atoms react to form a neutral argon atom, an argon ion and an electron. In addition to volumetric reactions, the following surface reactions are implemented:

Table 2: table of surface reactions.
reaction	formula	Sticking Coefficient
1	Ars=>Ar	1
2	Ar+=>Ar	1

When a metastable argon atom makes contact with the wall, it reverts to the ground state argon atom with some probability (the sticking coefficient).

electrical excitation

From an electrical point of view, the GEC reactor behaves as a transformer. A current is applied to the driving coil (the primary) and this induces a current in the plasma (the secondary). The plasma then induces an opposing current back in the coil, increasing its resistance. The current flowing in the plasma depends on the current applied to the coil and the reaction kinetics. The total plasma current can vary from no current (plasma not sustained) to the same current as the primary, which corresponds to perfect coupling between the coil and the plasma.

In this example a fixed power of 1500 W is applied to the coil. Some of this power is dissipated in the coil, some is deposited into the plasma.

Results and Discussion

In this section, result are presented from the model. For reference, results of the fluid model when assuming a Maxwellian EEDF are also presented.

The macroscopic quantities that change the most are the electron temperature and the number density of the argon excited state. When comparing the simulations results of the models using an assumed and computed EEDFs one observes that the model with the computed EEDF: have an electron temperature 1 eV higher in most of the plasma bulk; the maximum value for the number density of the argon excited state is only about 15% larger, however the spatial distribution is considerable different; the electron density is about 15% lower and the profiles are similar; the power absorbed by the plasma maintains almost the same absolute value but there is a deeper penetration of the absorbed power into the plasma bulk. This is because the lower plasma density implies a larger skin depth.

Figure 10 presents computed EEDFs for several axial positions between the wafer and the dielectric window. When comparing with Maxwellian EEDFs with the same mean energy (not shown in the figure), the computed EEDFs tend to be less populated in the regions of high energy above 16 eV and low energy below 5 eV, and tend to be more populated in the mid energies between 5 and 16 eV. The electron-electron collisions cause the EEDF to tend toward a Maxwellian distribution function, populating the high energy tail that otherwise would present a much steeper cutoff around the ionization energy region of 16 eV.

When answering the question “which EEDF should I use?” one can conclude that if the goal of a work is to study the electric properties of the ICP discharge the use of a Maxwellian EEDF (or a generalized EEDF with a less populated tail) could be enough since the electron density and the power absorbed by the plasma are fairly similar. However, if the goal is to study a complex chemistry that strongly depends on the electron temperature the model with computed EEDF might be necessary.

Figure 2: Electron density obtained with the computed EEDF.

Figure 3: Electron density obtained with a Maxwellian EEDF.

Figure 4: Electron temperature obtained with the computed EEDF.

Figure 5: Electron temperature obtained with a Maxwellian EEDF.

Figure 6: Number density of excited argon atoms obtained with the computed EEDF.

Figure 7: Number density of excited argon atoms obtained with a Maxwellian EEDF.

Figure 8: Power deposition into the plasma obtained with the computed EEDF.

Figure 9: Power deposition into the plasma obtained with a Maxwellian EEDF.

Figure 10: Computed EEDFs at a radial position of 3 cm and for several axial positions. For reference, it is presented a Maxwellian EEDF with the mean energy of 5.4 eV, which is the value found at z=3.5 cm for the model with the computed EEDF.

References

1. G.J.M. Hagelaar and L.C. Pitchford, “Solving the Boltzmann Equation to Obtain Electron Transport Coefficients and Rate Coefficients for Fluid Models,” Plasma Sources Sci. Technol., vol. 14, pp. 722–733, 2005.

2. D.P. Lymberopolous and D.J. Economou, “Two-Dimensional Self-Consistent Radio Frequency Plasma Simulations Relevant to the Gaseous Electronics Conference RF Reference Cell,” J. Res. Natl. Inst. Stand. Technol., vol. 100, p. 473, 1995.

3. Phelps database, www.lxcat.net, retrieved 2017.

Plasma_Module/Space-Dependent_EEDF_Modeling/argon_gec_icp_boltzmann

Modeling Instructions

Root

The following instructions show how to create a model for a ICP reactor that solves plasma fluid type equations fully coupled with the Boltzmann equation in the two-term approximation.

Go to the and open the argon_gec_icp model.

Select to solve the Boltzmann equation in the two-term approximation, add the contribution from electron-electron collisions and the oscillating electric field.

Set the number of elements in the extra dimension to 50, refined the mesh in the extra dimension in the low energies region, and choose to compute the maximum energy automatically.

From the menu, choose .

From the Application Libraries root, browse to the folder Plasma_Module/Inductively_Coupled_Plasmas and double-click the file argon_gec_icp.mph.

Plasma (plas)

In the window, expand the node, then click .

In the window for , locate the section.

From the list, choose .

Select the check box.

In the ω/N text field, type 2.11E-14[m^3/s].

In the N text field, type 50.

In the R text field, type 30.

Select the check box.

Choose to define the transport parameters using the . This way, the electron mobility is obtained from suitable integration of the EEDF over the cross section for momentum transfer. The remaining transport parameters are computed using the electron mobility.

Plasma Model 1

In the window, expand the node, then click .

In the window for , locate the section.

From the list, choose .

Multiphysics

Choose to define the plasma conductivity using an effective collision frequency for momentum transfer computed from the EEDF.

Plasma Conductivity Coupling 1 (pcc1)

In the window, expand the node, then click .

In the window for , locate the section.

From the list, choose .

Since the computation time of this type of problems are very long the mesh is coarsened.

Mesh 1

In the window, expand the node.

Size 1

In the window, expand the node, then click .

In the window for , locate the section.

From the list, choose .

Click

Mesh 1

In the window, collapse the node.

Multiphysics

In the window, collapse the node.

Plasma (plas)

In the window, collapse the node.

Label the study and group the plots so that it is easy to identify that they refer the model that uses a Maxwellian EEDF.

Maxwellian EEDF

In the window, click .

In the window for , type Maxwellian EEDF in the text field.

Add a study to solve for the EEDF only. The solution of this study is going to be used as the initial condition for the study that solves the fluid type equations coupled with the Boltzmann equation.

Add Study

In the toolbar, click