COMSOL 6.4 - In-Plane Microwave Plasma

Wave-heated discharges can be very simple, such as when a plane wave is guided into a reactor using a waveguide, or very complicated, such as in the case with ECR (electron cyclotron resonance) reactors. In this simple example, a wave is launched down a waveguide where it intersects a flowing gas at low pressure, resulting in formation of an argon plasma. Microwave plasmas typically have high number density without requiring significant power absorption. The plasma potential is also quite low compared to capacitive or DC discharges. Therefore, microwave plasmas share many of the characteristics of inductive discharges.

Figure 1: Diagram of geometry modeled. A TE or TM mode wave enters from the top port and intersects the gas flow leading to the formation of a plasma.

The model requires the Plasma Module and RF Module.

Model Definition

Electron transport is modeled by solving the continuity equation, the momentum equation under the drift-diffusion approximation, and the mean electron energy equation (for detailed information on electron transport, see Theory for the Drift Diffusion Interface in the Plasma Module User’s Guide)

The source coefficients in the above equations are determined by the plasma chemistry. The electron rate expression is defined as

where νe,j is the stoichiometric coefficient, and the reaction rate is defined as

where kjf is the forward rate constant and kjr is the reversed rate constant. Both the Electron Impact Reaction feature and Reaction feature can contribute to the electron rate expression. However, when using the Reaction feature it is important to note that the associated electron energy gain or loss is not included in the source term of the electron mean energy equation.

The rate constants can be computed from electron impact cross-section data

where γ = (2q/me)1/2 (SI unit: C1/2/kg1/2), me is the electron mass (SI unit: kg), ε is the electron energy (SI unit: V), σ is the electron impact collision cross section (SI unit: m2), and f is the electron energy distribution function.

When Townsend coefficients are used, the reaction rate is defined as

where αj/Nn is the reduced Townsend coefficient for reaction j (SI unit: m2) and Γe is the electron flux as defined above (SI unit: 1/(m2·s)). Townsend coefficients can increase the stability of the numerical scheme when the electron flux is field driven as is the case with DC discharges.

The total electron energy loss or gained is calculated by summing the collisional energy changes from all reactions defined with the Electron Impact Reaction feature as

where Δεj is the energy loss from reaction j (SI unit: V) and F is the Faraday constant (SI unit: C/mol). For excitation and ionization collisions Δεj corresponds to the energy of the excited state being excited/deexcited or ionized, for attachment Δεj is set to zero, and for elastic collisions

where me and mk are the electron and heavy species mass in kg, Te is the electron temperature in eV, and Tgas is the gas temperature in K.

For heavy species, the following equation is solved for the mass fraction of each species (for detailed information on the transport of the nonelectron 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.

In a microwave reactor the high frequency electric field is computed in the frequency domain using the following equation:

The relationship between the plasma current density and the electric field becomes more complicated in the presence of a DC magnetic field. The following equation defines this relationship:

Here, σ is the plasma conductivity tensor, which is a function of the electron density, collision frequency, and magnetic flux density. Using the definitions:

where q is the electron charge, me is the electron mass, ne is the collision frequency, and ω is the angular frequency of the electromagnetic field. In this example, the inverse of the plasma conductivity is diagonal because there is no external DC magnetic field:

The gas flow is modeled assuming a constant velocity in the x direction.

plasma chemistry

The chemical mechanism for the plasma consists of only 3 species and 7 reactions (electron impact cross section are obtained from Ref. 2):

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	-

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

electrical excitation

The Port boundary condition is used to excite the plasma. A total absorbed power of 30 W is fed into the port.

In a second study, the electrical excitation is changed to the TM mode, where the electric field has only an in-plane component. The total absorbed power is the same as the TE mode case.

Results and Discussion

The electron density is plotted in Figure 2 and peaks slightly downstream of the crossing point. The electron density is also slightly asymmetric in the y-plane due to the fact that the electromagnetic waves are absorbed asymmetrically. The electron “temperature” is plotted in Figure 3. The electron temperature is relatively low everywhere, in part due to the high operating pressure (1 Torr). The electron “temperature” peaks directly underneath the waveguide where the wave is absorbed. The norm of the electric field can be seen in Figure 5. The electric field is high inside the waveguide and there are no losses. Once the wave is exposed to the plasma, the energy is absorbed by the electrons, raising the electron temperature enough to generate new electrons through ionization. The ionization rate is high enough to sustain the plasma. The contour of the critical plasma density is also plotted in Figure 5. The electromagnetic wave cannot penetrate into regions exceeding the critical plasma density. Since the electron “temperature” is relatively low, one would expect the plasma potential to be low. The plasma potential is plotted in Figure 4 and is only around 10 volts.

In the TE mode, electrons do not experience any change in the high-frequency electric field during the microwave time scale. This means that the phase coherence between the electrons and electromagnetic waves is only destroyed through collisions with the background gas. The loss of phase coherence between the electrons and high-frequency fields is what results in energy gain for the electrons. Therefore, the momentum collision frequency is simply given by:

where νe is the collision frequency between the electrons and neutrals.

When switching to the TM mode, the electron density, electron “temperature” and plasma potential are quite similar to the TE mode case. This can be seen in Figure 7, Figure 8 and Figure 9. The electric field is very different however, Figure 9. The electric field cannot penetrate past the contour of critical electron density, and has its greatest magnitude in this location. The power deposition, Figure 10 is also highly localized to the contour of critical electron density. The TM mode causes in-plane motion of the electrons on the microwave time scale, so in regions where the high-frequency electric field is significant (the contour where the electron density is equal to the critical density), the time-averaged electric field experienced by the electrons may be nonzero. This destroys the phase coherence between the electrons and the fields, causing the electrons to gain energy. This is an example of a nonlocal kinetic effect, which is difficult to approximate with a fluid model. However, since this effect is similar to collisions with a background gas, the nonlocal effects can be approximated by adding an effective collision frequency to the momentum collision frequency:

where νeff is the effective collision frequency to account for nonlocal effects. In this example, since the Doppler broadening parameter is set to 20, this corresponds to an effective collision frequency of:

Figure 2: Plot of the electron density.

Figure 3: Plot of the electron temperature in the reactor.

Figure 4: Plot of the plasma potential in the reactor.

Figure 5: Plot of the electric field norm. The white contour represents the critical plasma density, where the electron density is equal to 7.6E16[1/m3].

Figure 6: Plot of the power deposition into the plasma. The white contour represents the critical plasma density, where the electron density is equal to 7.6E16 1/m3.

Figure 7: Electron density for the TM mode case.

Figure 8: Plot of the electron temperature for the TM mode case.

Figure 9: Close up of the high frequency electric field norm for the TM mode case.

Figure 10: Close up of the power deposition into the plasma for the TM mode case.

References

1. M.A. Lieberman and A.J. Lichtenberg, Principles of Plasma Discharges and Materials Processing, John Wiley & Sons, 2005.

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

Plasma_Module/Wave-Heated_Discharges/inplane_microwave_plasma

Modeling Instructions

From the menu, choose .

New

In the window, click

Model Wizard

In the window, click

In the tree, select > .

Click .

Click

In the tree, select > .

Click

Global Definitions

Parameters 1

In the window, under click .

In the window for , locate the section.

In the table, enter the following settings:


Name	Expression	Value	Description
P0	30[W]	30 W	Absorbed power

Geometry 1

Rectangle 1 (r1)

In the toolbar, click

In the window for , locate the section.

In the text field, type 0.05.

In the text field, type 0.1.

Click

Rectangle 2 (r2)

In the toolbar, click

In the window for , locate the section.

In the text field, type 0.25.

In the text field, type 0.05.

Locate the section. In the text field, type -0.1.

In the text field, type 0.1.

Rectangle 3 (r3)

In the toolbar, click

In the window for , locate the section.

In the text field, type 0.05.

In the text field, type 0.1.

Locate the section. In the text field, type 0.15.

Click

Click the

button in the toolbar.

Definitions

Walls

In the toolbar, click

In the window for , locate the section.

From the list, choose .

Select Boundaries 1–3, 6, 8, 11, and 13 only.

In the text field, type Walls.

Plasma (plas)

In the window, under click .

In the window for , locate the section.

Click

Select Domain 1 only.

Cross Section Import 1

In the toolbar, click

and choose .

In the window for , locate the section.

Click

Browse to the model’s Application Libraries folder and double-click the file Ar_xsecs.txt.

Click

Reaction 1

In the toolbar, click

and choose .

In the window for , locate the section.

In the text field, type Ars+Ars=>e+Ar+Ar+.

Locate the section. In the kf text field, type 3.734E8.

Reaction 2

In the toolbar, click

and choose .

In the window for , locate the section.

In the text field, type Ars+Ar=>Ar+Ar.

Locate the section. In the kf text field, type 1807.

Species: Ar

In the window, click .

In the window for , locate the section.

Select the checkbox.

Species: Ar+

In the window, click .

In the window for , locate the section.

Select the checkbox.

Surface Reaction 1

In the toolbar, click

and choose .

In the window for , locate the section.

In the text field, type Ars=>Ar.

Locate the section. From the list, choose .

Surface Reaction 2

In the toolbar, click

and choose .

In the window for , locate the section.

In the text field, type Ar+=>Ar.

Locate the section. From the list, choose .

In the window, click .

In the window for , locate the section.

Select the checkbox.

Locate the section. Select the checkbox.

Plasma Model 1

In the window, click .

In the window for , locate the section.

Specify the u vector as


10	x
0	y

In the T text field, type 350.

In the pA text field, type 1[torr].

Locate the section. In the μeNn text field, type 4E24.

Initial Values 1

In the window, click .

In the window for , locate the section.

In the ne,0 text field, type 1E17.

Wall 1

In the toolbar, click

and choose .

In the window for , locate the section.

From the list, choose .

Ground 1

In the toolbar, click

and choose .

In the window for , locate the section.

From the list, choose .

Electron Outlet 1

In the toolbar, click

and choose .

Select Boundary 14 only.

Outflow 1

In the toolbar, click

and choose .

Select Boundary 14 only.

In the window for , locate the section.

Select the checkbox.

Electromagnetic Waves, Frequency Domain (emw)

In the window, under click .

In the window for , locate the section.

From the list, choose .

Port 1

In the toolbar, click

and choose .

Select Boundary 9 only.

In the window for , locate the section.

From the list, choose .

Select the checkbox.

In the Pdep text field, type P0.

Materials

Material 1 (mat1)

In the window, under right-click and choose .

Select Domains 2 and 3 only.

In the window for , locate the section.

In the table, enter the following settings:


Property	Variable	Value	Unit	Property group
Relative permittivity	epsilonr_iso ; epsilonrii = epsilonr_iso, epsilonrij = 0	5	1	Basic
Relative permeability	mur_iso ; murii = mur_iso, murij = 0	1	1	Basic
Electric conductivity	sigma_iso ; sigmaii = sigma_iso, sigmaij = 0	0	S/m	Basic