- Email: [email protected]
Modelling the electromagnetic field and plasma discharge in a microwave plasma diamond deposition reactor
ELSEVIER
Diamond
and Related
Materials
4 (1995) 1145-1154
Modelling the electromagnetic field and plasma discharge in a microwave plasma diamond deposition reactor W. Tan, T.A. Grotjohn Department of Electrical
Engineering, Michigan State University, East Lansing, MI 48824, USA
Received 7 September 1994; accepted in final form 29 March 1995
Abstract A numerical model which includes an electromagnetic field model and a fluid plasma model has been developed for a microwave cavity plasma reactor used for diamond thin film deposition. The microwave plasma reactor modelled is a cylindrical, single mode excited cavity with an input power probe for coupling the microwave energy into the cavity. The time-varying electromagnetic fields inside the resonant cavity are obtained by applying the finite difference time domain (FDTD) method to solve Maxwell’s equations. The microwave electric field interactions with the plasma discharge are described using a finite difference solution of the electron momentum transport equation. The characteristics of the discharge are simulated using a fluid plasma model which solves the electron and ion continuity equations, electron energy balance equation, and the Poisson equation. The spatial electric field patterns, power absorption patterns, and quality factor of the cavity loaded with a hydrogen discharge are investigated in the moderate pressure range (40-60 Torr). The physical behaviour of the diamond deposition discharge, such as plasma density, electron temperature, and plasma potential, are also simulated and analysed for various input conditions. The simulated results are compared with experimental data. Keywords:
Microwave discharge; Plasma modelling; Diamond deposition reactor; Electromagnetic
1. Introduction Microwave excited plasmas have demonstrated good potential in moderate pressure (l-150 Torr) plasma processing applications including diamond film deposition Cl]. One way to create a microwave excited plasma is utilizing a single mode microwave resonance cavity which couples microwave energy into the discharge via joule heating [2]. The resonant electromagnetic field energy is imparted to the electron gas, where the energy is then transferred to ions and neutrals by elastic and inelastic collisions. In order to understand the behaviour of the plasma inside a cavity reactor source, the microwave fields, gas transport and energy transport, as well as their relations, need to be investigated carefully. The coupling among these is complex and analytic models are hard to apply to handle this problem. Thus, a systematic numerical model is required to solve the plasma behaviour and the exciting microwave fields selfconsistently. One type of numerical model for plasma discharge simulation is the fluid (continuum) plasma model, which solves the moments of the Boltzmann equation, including Elsevier Science S.A. SSDI 0925-9635(95)00291-X
modelling
the continuity equation, momentum transport equation and energy transport equation. This technique has been widely used to investigate the behaviour of the discharges inside parallel plate r.f. plasma reactors [ 3,4]. The electric fields used for plasma excitation inside the reactor are solved by Poisson’s equation with the boundary conditions assigned on the electrodes. In terms of microwave discharge simulations for diamond deposition, the microwave electric field distribution has been assumed uniform [ 51. In this paper, a numerical model is developed which self-consistently solves the electromagnetic fields as well as the excitation mechanism and characteristics of the diamond deposition discharges inside a cylindrical microwave cavity plasma reactor. The electromagnetic fields are solved by using a finite difference time domain (FDTD) [6] method. Previous work successfully adapted the FDTD method to solve the electromagnetic field distribution inside a compact ECR source and also compared the FDTD method solution directly to analytical solutions on simple cavity geometries [7]. Moreover, the FDTD method has been applied to study the electromagnetic fields inside a microwave cavity
1146
W. Tan, T.A. Grotjohn/Dinmond
und Related Mat&u/s
plasma reactor used for diamond thin film deposition by assuming a uniform plasma load inside the cavity [S]. In this paper this uniformity assumption is removed, and the characteristics of the discharge are simulated by a fluid plasma model which solves the electron and ion continuity equations, the electron energy balance equation, and the Poisson equation. A final electromagnetic field and plasma discharge solution is obtained by iteratively solving the electromagnetic field model and the fluid plasma model. The details of coupling the FDTD electromagnetic field model and the fluid plasma model are provide in the following sections. In particular, a TMO13 mode, 17-78 cm i.d. microwave cavity ,plasma reactor loaded with a hydrogen discharge is simulated using this numerical model. The TM,,, mode is an azimuthally symmetric mode. The hydrogen plasma example was chosen because diamond film deposition processes often consist of high percentages of hydrogen in the discharge. And, as shown by Koemtzopoulos et al. [S], adding small percentages (e.g. 1%) of methane to a hydrogen discharge has only a minimal effect on the electron energy distribution function. The method used to excite electromagnetic waves inside the cavity and the boundary condition assignments for both the electromagnetic field model and fluid plasma model will be discussed. The methods to evaluate the cavity quality factor Q will be presented. The spatial electric field patterns, power absorption patterns and quality factor of the cavity loaded with a hydrogen discharge are investigated in the moderate pressure range of 40-60 Torr. The physical behaviour of the hydrogen discharge, such as plasma density, electron temperature and plasma potential, are also simulated and analysed for various input conditions. The simulated results are compared with experimental data.
2. Model description The numerical model developed in this study includes an FDTD electromagnetic field model and a fluid plasma model. These numerical models and the coupling of these two models are described in this section. 2.1. FDTD electrtzmagnetic$eId model
4 (1995) 1145-1154
density. The FDTD method is formulated by discretizing Eq. (1) with a centred difference approximation in both the time and space domains. In this study, Eq. (1) is discretized in two-dimensional (P and z) cylindrical coordinates [S]. One of the major motivations to solve Maxwell’s equations is to investigate the power dissipation inside the microwave plasma source due to the presence of plasma discharges. The power dissipation density, Pa&t), with a power absorbing load (such as a discharge) present is Pa&t)
= J(r,WW)
(2)
The current density, J, which is induced by the microwave fields, can be determined by solving the momentum transport equation ofelectrons. For example, if ion motion is ignored, the momentum transport equation of electrons can be written as m, $ v(r,t) = - eE(r,t) - mev&r,t)Y(r,t)
(3)
which is also called the Langevin equation [9], and the current density as J (r,t) = - @&Wr,r)
(4)
where v is the average electron velocity, e is the electron charge, m, is the electron mass, veff is the effective collision frequency [!I] and n, is the electron density. Eqs. (3) and (4) are solved by the finite difference method [S]. From the above equations, the spatially dependent current density J is determined by the spatially dependent discharge characteristics, such as electron density n, and effective collision frequency veff. Moreover, the effective collision frequency which describes the momentum transfer due to electron-neutral elastic collisions is a function of the electron energy distribution (electron temperature). For a complete microwave plasma excitation solution, the electron density and electron temperature must be solved self-consistently by using a fluid plasma model. 2.2. Fluid plasma mode!
The electromagnetic fields inside the cavity can be described by the time-dependent Maxwell’s equations, which are written as
0) where E and H are the electric and magnetic fields, p is the permeability, E is the permittivity and J is the current
The behaviour of the charged particles in a weakly ionized gas can be described by the particle, momentum and energy balance equations for electrons and ions, which are obtained from the moments of the Boltzmann equation. Moreover, these equations are often combined with Poisson’s equation to provide a self-consistent space charge field. Fur the moderate pressures (I150 Torr) of this study, the momentum balance equations of electrons and ions are given by Eqs. (8) and (9) according to Ref. [lo). In the steady state, these equa-
W. Tan, T.A. Grotjohn/Diamond and Related Materials 4 (1995j 1145-1154
tions are written as [ 111: (5) V.J, = nennkion - apine
(6)
V.Ji = nennkion - a,nin,
(7)
- D,Vn,
(8)
J,= -n,p,E
Ji = nil”iE - Di V ni V.q, = - 4 J,.E -2
-
(9)
nenn(Ei&on + Eextkxt + Edd
2
i ksl&in,
qe= ; kgzJe-
c
; k,D,n,
(10)
>
UT,
(11)
Eq. (5) is the Poisson equation, where Y is the potential (E = - Vlu). Eqs. (6) and (7) are the electron and ion continuity equations, respectively. Eqs. (8) and (9) are in essence the momentum balances for electrons and ions, respectively. Eq. (10) is the electron energy balance equation, with the total electron energy flux given by Eq. (1 1), where the thermal energy is assumed to be greater than the kinetic energy. In the above equations, n, and ni are the electron and ion densities, respectively; J, and Ji are the electron and ion fluxes, respectively; T, is the electron temperature; qe is the electron energy flux; kionr kext and kdis are the inelastic rate constants, with threshold energies &ion,E,,, and Edis for ionization, excitation and dissociation, respectively; and CI, is the recombination rate constant. Note that elastic losses have been included in Eq. (lo), where \I,,, is the electron-neutral momentum transfer frequency and m, is the species mass. Moreover, D,,i and ~~,i are the electron and ion diffusivities and mobilities, respectively. The fluid plasma model describes the characteristics of discharges by solving Eqs. (5)-( 11). Finite difference methods and staggered mesh techniques [12] are used to discretize these equations in cylindrical coordinates. Since the resonant electromagnetic modes of the microwave plasma sources used for diamond thin film deposition in this study are normally TMoln mode, the electric field spatial distributions are C$symmetric in nature and the E+ component is zero. Therefore the discharge behaviour can be assumed to be 4 symmetric and the discretization of the equations reduced to a two-dimensional problem. Thus, the simulation region remains in the r-z plane only. Note that the heating term qJ*E in Eq. (10) primarily represents the microwave power absorbed by the plasma, which is determined by the microwave
1147
electric field and the microwave induced current density using the FDTD electromagnetic field model described earlier. Therefore, for simplicity, this heating term given in Eq. (2) as Pabs can be viewed as an input parameter for the electron energy balance equation. The flow chart of the microwave cavity plasma simulation is shown in Fig. 1. The solution of the discretized fluid plasma equations, Eqs. (5)-( 1l), on a grid of N nodes involves a total of 4N unknowns for the fluid model since at each point the unknowns are Y, n,, ni and T,. To solve the system of 4N non-linear type discretized equations, both direct method (Newton’s method) and iterative method techniques are used. In particular, the Poisson equation, electron continuity equation and ion continuity equation are tightly coupled and solved by Newton’s method. The unknowns solved at each grid point are the electron density, ion density and plasma potential, and the electron and ion flux (steady state d.c. flux) can be determined by Eqs. (8) and (9) using the solution of these unknowns. Then, the electron density and flux are coupled into the discretized electron energy balance equation, Eq. (lo), to solve the electron temperature. The calculated electron temper-
Solve electron energy balance equation to determine electron temperature(steady
state)
Solve Poisson equation and continuity equations to determine particle densities, flux and plasma potential (steady state) by Newton’s method
I
Solve Maxwell’s equation by FDTD model to determine the power absorption
‘ig. 1. Flow chart for microwave
cavity plasma
reactor
simulation.
1148
W. Tun, T.A. GrotjohnJDiamond and Related Materials 4 (1995) 1145-1154
ature is then coupled back to the continuity equations and Poisson equation to modify the reaction rates and parameters and update the electron density, ion density and plasma potential. The final stable solution of density and temperature is achieved by iteratively solving the continuity/Poisson equations and the energy balance equation. As described above, an input term in the discretized electron energy balance equation is the absorbed microwave power density Pabs, which comes from the FDTD electromagnetic field model. The FDTD model which solves the electromagnetic fields inside the discharges provides output information, such as the microwave power absorption of the discharges, by simultaneously solving the microwave electron momentum transport equation, Eq. (3). The discharge characteristics, such as the plasma density and electron temperature, are determined in the steady state under this power absorption condition. The discharge characteristics information is then coupled to the FDTD model to modify the plasma conductivity and calculate a new discharge power absorption. Therefore, in this iterative manner, the power absorbed by the discharges will converge to a stable value, and the electromagnetic fields inside the reactor source and the plasma characteristics are solved selfconsistently.
3. Implementation for a microwave cavity plasma reactor
which the microwave energy is coupled into the cavity through a coaxial input probe. The specific resonant mode is determined by the geometrical size of the cylindrical cavity, and it can be adjusted by the movable sliding short. In this study, the height of the cavity is adjusted for TM0i3 mode excitation with a 2.45 GHz input microwave frequency. The coaxial input probe and movable short position are adjusted for maximum power transfer to the plasma load, as indicated by minimum reflected power from the cavity. A hydrogendominated discharge is excited and sustained inside the quartz dome by the input microwave power. The volume and the temperature of the discharge are dependent on the input power and pressure. The substrate and substrate holder are placed on the top of a quartz tube which is used to adjust the height of the substrate. A metal screen is placed at the bottom of the cavity to prevent the microwave field from propagating into the downward region. The simulation structure for the FDTD model and the plasma fluid model is shown in Fig. 3. The sliding short, the base plate, the substrate holder, the coaxial probe (both the inner and outer conductor) and the cavity sidewalls form the boundary of the FDTD simulation region. The boundary conditions for the electric fields on these surfaces are that only the normal components of the electric fields exist and the tangential electric fields on these surfaces are zero. This is based on assuming that these boundaries are formed by perfect conductors. At the input end of the coaxial probe is an
In this study, the numerical model is applied to simulate a microwave cavity plasma reactor used for diamond thin film deposition, as shown in Fig. 2. The cavity is a cylindrical, 17.78 cm i.d. microwave cavity in
Wave Excitation /Plane (TEM Wave)
Perfect Conductor
kFig. 2. Configuration of a microwave cavity plasma reactor diamond film deposition. The cavity diameter is 17.78 cm.
used for
Fig. 3. Simulation structure reactor with quartz tube.
cross-section
Quartz Tube (q~3.75) for microwave
cavity plasma
W. Tan, T.A. Grotjohn/Diamond
and Related Materials 4 (1995) 1145-1154
open boundary where the region in which the field has to be computed is unbounded. At this open boundary a truncation method is used to prevent any artificial reflection of outgoing waves. This non-reflecting boundary condition has been investigated extensively in solving electromagnetic wave scattering problems [ 131. Basically, this boundary condition allows electromagnetic waves arriving at the open boundary of the coaxial probe to leave the simulation region without any artificial reflections at the termination of the grid structure. The non-reflecting boundary condition which is applied in this model can be written as: (E:+‘(i,k,)
1149
+ ~C~%41’21- WlNW(w~1’21 + Ad)
AzrAr2rr) are approximately equal. For the FDTD electromagnetic field model, the number of grids used in the simulation is 30 x 83 in the r- and z-directions respectively, and the time step is 0.5 ps. The grid spacing for the FDTD model is uniform in both r- and zdirections. Since the grid structure spacing for the FDTD model is different from that of the fluid plasma model, a linear interpolation technique is used to interpolate the absorbed power to each grid point of the plasma model, and the plasma density and electron temperature to each grid point of the FDTD model. For the partially ionized and partially dissociated H2 plasma discharge simulation, the major particle interaction processes are the electron-H, molecule inelastic collision, electron-H, molecule elastic collision and electron-hydrogen ion recombination. The electron-H, inelastic collisions include the H2 molecule ionization, excitation and dissociation processes. The rate coefficients in Eqs. (6), (7) and ( 10) for these collision processes can be expressed using the Arrhenius relationship [ 111 as
x [Ez+‘(i,(k, - 1)) - Ez(i,k,)]
kion= &n exp ( - aion/KBT,)
= E:(i,(k,-
1))
+ KWW’21
- WlK~~/W”“1
x [E:“(i,(k,-
1)) - E:(i,k,)]
+ AZ)) (12)
and E$+‘(i,k,,) = E$(i,(ko - 1))
(13)
where i denotes the grid location in the r-direction and k, is the grid terminating point in the z-direction. Additionally, the quartz dome and the quartz tube, where the dielectric constant is equal to 3.75&,, are included in this model. The boundaries of the fluid plasma simulation are the quartz dome and the substrate holder. The diameter of the substrate holder is assumed to be the same as that of the quartz tube. The boundary conditions of the fluid plasma model are, at the substrate: n,=ni=O T,= T,
(14)
Y=O
kext= A,,, exp ( - GJGT,) his = &is exp ( - EsdKBTe)
where &ion, eext and ldis are the threshold energies for H,
molecule ionization, excitation and dissociation and Ai,,, A,,, and Adis are the pre-exponential factors, which are obtained by approximating the rate constant data at low electron temperatures to these relationships [11,15]. For simplicity, only the reactions with higher rate coefficients are considered in this study for H, discharges. The types of inelastic collisions and their corresponding rate parameters are summarized in Table 1. The collision frequency for electron-H, molecule momentum transfer is relatively independent of the electron temperature, so it can be written as [9]
and at the quartz wall: n,=n,=O
v,(H2)=
1.44 x 1012 x
Pressure (Torr) (17)
T,(K) Table 1 H, reaction
The boundary conditions at the quartz wall and the substrate wall assume that recombination occurs when electrons and ions reach the walls and that thermal equilibrium exists at the walls. The substrate and the substrate holder can be either grounded or floating. In the simulations of this paper the substrate potential is assumed grounded. The potential is assumed to be at the floating potential Yr on the quartz disk [ 141. The grid structure used in the plasma simulation is in the Y-and z-directions. In the z-direction a uniform grid spacing is used, and in the r-direction the grids are constructed such that all the unit cell areas (or volumes,
(16)
rates
Reaction
Expression
Ionization
e+H*-+e+H:+e H; + H,+H:
+ H
Threshold energy
Pre-exponential factor or rate coefficient (m3 s-j)
15.4 eV OeV
1.0 X 10-l” H+ is the diminant ion species
Excitation
e+H,-+Hf+e
12.0 eV
6.5 x lO_”
Dissociation
e+H,+e+H+H
10.0 eV
1.0 x lo-‘4
e + ion-neutral
OeV
1.0 x lo-‘4
Recombination
W. Tan, T.A. GrotjohnjDiamond and Related Materials 4 (1995) 1145-1154
1150
where T, is the neutral temperature, which can be represented by the translational temperature of H2 gas. Since v, is relatively independent of the electron temperature, the effective collision frequency veff is set equal to V It should be noted that the only neutral species ctnsidered in the ion and electron simulations was the Hz species and that the dominant ionic species in the plasma is Hz [5,11-J. Hence these simulations assume a low hydrogen dissociation percentage.
showed no significant variation along the #-direction as expected for the TM,,s mode. Once the electromagnetic field strength and mode is determined inside the cavity, the cavity quality factor Q of the plasma discharge loaded cavity reactor can be found by calculating the electromagnetic energy U stored in the cavity and the power absorbed in the cavity Pa,,*. The energy stored in the cavity U is expressed as U
1E12dV
=E
s 4. Experimental results
(18)
and the Q factor is given by
In order to verify the cavity resonant mode, the radial component of electric field strength at the cavity outer wall was measured by calibrated micro-coax electrical probes inserted through the cavity wall versus z and 4 [ 11. The centre conductor of the probe is inserted 2 mm beyond the inside surface of the cavity side walls, and the other end is connected to a microwave power meter. Since the probe power reading is proportional to the square of the r.m.s. electric field strength normal to the inside cavity, the cavity resonant mode and the cavity stored energy can be determined using this technique [l]. Probe experiments were done by measuring the electric field strength at various z positions along the cylindrical cavity wall of the microwave plasma reactor during diamond thin film deposition. The microwave excitation frequency was 2.45 GHz, and the cavity height was adjusted to a position near 20.4 cm for a minimum reflected power condition. The result is shown in Fig. 4, where the resonant condition for the microwave cavity during a deposition process matches closely to the ideal TM0r3 mode. Additionally, the electric field was probed circumferential around the cavity. The field strength 2.5,
* Experiment
Theory
I
Q=E
where o is the excitation frequency. The measured Q value at a pressure of 50 Torr and an input microwave power of approximately 1500 W was 100. A Q value of 100 gives an indication as to why the microwave energy pattern with the plasma load present is similar to the pattern of an empty resonant cavity structure.
5. Simulation results The microwave cavity plasma reactor loaded with a H2 discharge for diamond film deposition was simulated and investigated by the self-consistent numerical model which has been described in the previous sections. The input parameters for this model include the pressure and input microwave power. A set of experimentally determined empirical equations developed by King [16] are used to establish the neutral temperature for this simulation. These empirical equations were obtained by doing discharge diagnostics across a parameter space including pressure and input microwave power variation. The empirical equations used for the translational temperature of H, gas and the discharge volume are [ 161: Translational
*
*
(19)
temperature (K)
= 228.6 + 374.3 x Incident power (kW) + 16.5 x Pressure (Torr) k 94.2
(20)
Plasma volume (cm3) = 449.7 + 116.2 x Incident power (kW) - 18.1 x Pressure (Torr) + 57.1 x [Incident power (kW)]* I
+ 0.25 x [Pressure (Torr)12 - 5.4 0.
0
a
2
* 4
6
8
10
12
14
16
18
20
Cavity Height (cm) Fig. 4. Microwave electric field strength (radial directed) vs. cavity height for microwave cavity plasma reactor during the diamond deposition process.
x Pressure (Torr) x Incident power (kW) + 15.4 (21) The translational temperature was determined using Doppler broadening measurements of the H atom emission and the plasma discharge volume was determined
W. Tan, T.A. GrotjohnjDiamond and Related Materials 4 (1995) 1145-1124
using photographs of the plasma. Eqs. (20) and (21) are valid for the pressure range from 35 to 6.5Torr and incident microwave power range from 1.4 to 2.6 kW. The plasma discharge volume serves as a boundary for the portion of the plasma simulation region where the neutral temperature (T,) is high. Inside the discharge volume region, the translational temperature given by Eq. (20) is the neutral temperature. Outside the volume, the neutral temperature is equal to the temperature assigned on the boundaries (T, = 1000 K). Between these two regions, a linear temperature change profile is used to prevent an abrupt change of neutral temperature. This is done by assigning the boundary temperature (T, = 1000 K) to the region several grids, namely three grids, away from the plasma volume in both the r- and z-directions. The temperature assigned on the grids between these two regions are linearly interpolated from the neutral temperature given by Eq. (20) and boundary temperature (T, = 1000 K) as shown in Fig. 5. For the FDTD electromagnetic field model, the technique used to excite the electromagnetic field in this study is to select grid points on a cross-section plane of the coaxial probe as source points and assign the timevarying electric field component at these points based at theoretical transverse electromagnetic (TEM) wave solutions in a coaxial structure, as shown in Fig. 4. The electromagnetic wave then propagates down into the cavity region where power is absorbed by the discharge. Any reflected electromagnetic wave will propagate to the end of the open boundary of the input power probe and be terminated by the non-reflecting boundary. The simulation results for the microwave electric field distribution are shown in Figs. 6 and 7 for the E, and E, components respectively. The input condition for this simulation is 50 Torr in pressure and 1500 W input microwave power. The electric field patterns basically follow the TMo13 mode electromagnetic field distribution. The effect of the electric field caused by the presence of discharge can be observed. It reduces the amplitude
Neutral temperature 18CUl
1600
14ca 0
1200
1090 0
Fig. 5. Neutral temperature region.
distribution
in the plasma simulation
Xi04
Sliding Short
1151
Quartz Dome
0 0
Plasma Discharge
\
Substrate
Fig. 6. E, field distribution at a pressure of 50Torr and an input microwave power of 1500 W.
XiOP
Sliding Short
Plasma Discharge
Quartz Dome
Substrate
Fig. 7. E, field distribution.
of the radial directed electric field component, E,, in the plasma discharge region since electromagnetic power is absorbed in the discharge region. The discontinuity of the electric field resulting from the difference of dielectric constant between the quartz and air can also be seen. It can also be observed that in the plasma simulation region the E, field has the highest magnitude just above the substrate. The conductivity of the plasma can be determined using the solved density and effective collision frequency. Fig. 8 shows that the maximum conductivity of the discharge is around 0.4 mho m. Fig. 9 shows the power absorption pattern versus I and z and indicates that the microwave power is highly absorbed in the central portion of the plasma region close to the substrate. The simulated electric field strength in Figs. 6 and 7, as well as the field measured in Section 4, indicate that the microwave electric field strengths are of the order of lo4 V m-‘. If this electric field strength is assumed to
W. Tun, T.A. Grotjohn/Diumond and Reluted Muterials 4 (1995) 1145-1154
1152
2
4
51
> 0 0
2 0 0
(a)
r(m)
(4
em)
10 6
(b)
r(m)
10 6
10 6
Substrate
Plasma Discharge
Fig. 10. Simulation results of a microwave cavity plasma reactor for pressure = 50 Torr and input microwave power = 1500 W. (a) Electron temperature; (b) plasma potential; (c) plasma density. Note that z =0 is located at the substrate position for the plasma simulations. (This same z location is at 4.5 cm for the electromagnetic field simulations.)
Fig. 8. Conductivity vs. r and z.
2 81 0 30
10
x lop' 10
/
0 0
Su&rate
Plasma Discharge
?
E 7
J
Fig. 9. Absorbed microwave power vs. r and z.
contribute to the heating of the substrate through direct electromagnetic energy absorption, then the amount of this direct electromagnetic heating can be estimated. Specifically, if a conductivity is assumed for the substrate holder similar to that of graphite, the absorbed power in the substrate and substrate holder is calculated to be less than 1% of the power entering the plasma source. Hence in the resonant cavity microwave system described in this paper the microwave electric fields are not strong enough to contribute substantially to the substrate heating; rather, the substrate is heated by the plasma. The simulation result of the plasma characteristics of a H2 discharge at 50 Torr and 1500 W input power is shown in Fig. 10. The plasma density is highest near the centre (r=O) of the discharge. The electron temperature is also highest at the centre of the discharge. The spatial variations of ionization and recombination rate are shown in Fig. 11. The ionization rate is higher near the substrate due to the electron temperature being higher
5
x
0 0 0
Fig. 11. Spatial variation of ionization rate (top) and recombination rate (bottom) for pressure = 50 Torr and input microwave power = 1500 w.
at that region. The recombination rate profile basically follows the electron density variation. The characteristics of Hz discharges and cavity quality factor Q were analysed for various pressures and absorbed powers by this numerical model. The neutral temperature and plasma volume determined by Eqs. (20) and (21) are the input data for simulations and their variations for different incident power and pressure are shown Fig. 12. The simulation results of average electron temperature vs. absorbed microwave power are shown in Fig. 13. These values are calculated by averaging the electron temperature across the plasma volume region. The average electron temperature is about 1.5 eV. The electron temperature increases with higher absorbed microwave power. This occurs because the measured neutral gas temperature increases with microwave
W. Tan, T.A. Grotjohn/Diamond and Related Materials 4 (1995) 1145-1154
1153
1.50
*: 40 Tom 0: 50 Torr x: 60 TOIT
x
I
x
40 Torr 50 Torr 0: 60 Torr *:
140-
0
x:
Y
130120-
2
G
o: 50
Torr x: 60 Tom
I_.*
,i 9
x
loo
B a
50
.z
0 1600
1500
x 1700
--_+__L 1800
ioco 1900
21W
2cwJ
Fig. 12. Experimentally measured neutral temperature volume vs. incident power for various pressures.
1500
2ow
J
2500
Power absorbed (watt)
and plasma
Fig. 14. Simulated Q value vs. absorbed power for three different pressures. The experimentally measured Q at a pressure of 50 Torr and an input microwave power of 1500 W is also indicated.
nique described in Section 4. The Q value determined experimentally is about 100, which is close to the simulation result (Q= 107). This close agreement is an indication of the validity of the model.
*: 40 Torr o: 50 Torr x: 60 Torr II
*
6. Summary
0
0.9
0
x
1600
2000
x
1
I
1200
1400
1603
2200
2400
Absorbed Power (Watt) Fig. 13. Simulated electron temperature vs. absorbed power for three different pressures.
power, as indicated in Eq. (20). The higher neutral gas temperature leads to a reduced neutral density and, hence, to a reduced electron-neutral collision rate at a given electron temperature. Also, the electron temperature drops with increases of pressure. This can be understood since higher pressure causes a higher collision frequency for electrons with neutrals. Therefore electrons more easily transfer their energy to neutral particles and reduce their kinetic energy. The results for Q vs. absorbed microwave power are shown in Fig. 14. It shows that when the input microwave power is increased the Q factor decreases. The pressure has only a small effect on the Q factor. At the 50 Torr pressure and 1500 microwave power condition, the Q factor of the microwave cavity plasma reactor with an Hz discharge was experimental determined by using the tech-
A self-consistent numerical model has been developed to simulate the electromagnetic excitation of discharges inside microwave cavity plasma reactors. This software includes an electromagnetic field model and a fluid plasma model. The microwave cavity plasma reactors simulated by this numerical software were loaded with a H, discharge and were used for diamond thin film deposition. The complex reactor geometry and input power coupling probe structure were included in the simulation. The input parameters for the microwave cavity plasma reactor simulation included the pressure (40-60 Torr) and input microwave power. The electromagnetic behaviour of the discharge loaded resonant cavity, such as the electric field distributions, power absorption patterns, and cavity quality factor Q, were studied and analysed for various input conditions. The characteristics of H, discharges, including the plasma density, electron temperature and plasma potential, were also investigated and studied. The electromagnetic mode and cavity Q value of an Hz discharge loaded reactor were simulated and shown to be in agreement with measured experimental values. The numerical model presented in this study allows the calculation of the power absorption profile in microwave discharges, which has often been a difficulty in plasma simulations. The understanding of the interactions between the microwave fields and the plasma
1154
W. Tan, T. A. GrotjohnlDiamond and Related Materials 4 ( 1995) 1145-I 154
discharges allows the reactor design to be analysed for improvement and optimization of such quantities as uniformity of microwave power absorption, uniformity of the plasma species and kinetic energy of the plasma species. Additionally, the simulation and understanding of the electromagnetic excitation of the discharge is expected to be important for microwave discharge diamond deposition machine control.
c31 J.P. Boeuf, Phys. Rev. A, 36 (1987) 2782. c41 F. Young and C-H. Wu, IEEE Trans. Plasma Sci., 21 (1993) 312. c51 C.R. Koemtzopoulos, D.J. Economou and R. Pollard, Diamond Relat. Mater., 2 (1993) 25. IL61K.S. Yee, IEEE Trans. Antennas Propagat., I4 (1966) 302. c71 T.A. Grotjohn et al., Rev. Sci. Instrum., 65 (1994) 1761. C81 W. Tan and T.A. Grotjohn, .I. Vat. Sci. Technol., 12 (1994) 1216. c91 B.E. Cherrington, Gaseous Electronic and Gas Lasers, Pergamon,
New York, 1979. Cl01 J.D.P. Passchier and W.J. Goeheer, .I. Appl. Phys., 74 (1993) 3744. Cl11 M. Surendra et al., J. Appl. Phys., 71 (1992) 5189. Cl21 MS. Barnes, T.J. Cotler and M.E. Elta, J. Camp. Phys., 77
(1988) 53.
References
Cl31 G. Mur, IEEE Trans. Electromagn. Compat., 23 (1981) 377. Cl41 D.M. Manos and D.L. Flamm, Plasma Etching, Academic Press,
San Diego, CA, 1989. J. Zhang, B. Huang, D.K. Reinhard and J. Asmussen, J. Vat. Sci. Technol. A, 8 (1990) 2124. [2] J. Asmussen, J. Vat. Sci. Technol. A, 7 (1989) 883. [l]
Cl51 R.K. Janev et al., Elementary Processes in Hydrogen-Helium Plasmas, Springer, Berlin, 1987. Cl61 G. King, Ph.D. Dissertation, Michigan State University, 1994.
Diamond
and Related
Materials
4 (1995) 1145-1154
Modelling the electromagnetic field and plasma discharge in a microwave plasma diamond deposition reactor W. Tan, T.A. Grotjohn Department of Electrical
Engineering, Michigan State University, East Lansing, MI 48824, USA
Received 7 September 1994; accepted in final form 29 March 1995
Abstract A numerical model which includes an electromagnetic field model and a fluid plasma model has been developed for a microwave cavity plasma reactor used for diamond thin film deposition. The microwave plasma reactor modelled is a cylindrical, single mode excited cavity with an input power probe for coupling the microwave energy into the cavity. The time-varying electromagnetic fields inside the resonant cavity are obtained by applying the finite difference time domain (FDTD) method to solve Maxwell’s equations. The microwave electric field interactions with the plasma discharge are described using a finite difference solution of the electron momentum transport equation. The characteristics of the discharge are simulated using a fluid plasma model which solves the electron and ion continuity equations, electron energy balance equation, and the Poisson equation. The spatial electric field patterns, power absorption patterns, and quality factor of the cavity loaded with a hydrogen discharge are investigated in the moderate pressure range (40-60 Torr). The physical behaviour of the diamond deposition discharge, such as plasma density, electron temperature, and plasma potential, are also simulated and analysed for various input conditions. The simulated results are compared with experimental data. Keywords:
Microwave discharge; Plasma modelling; Diamond deposition reactor; Electromagnetic
1. Introduction Microwave excited plasmas have demonstrated good potential in moderate pressure (l-150 Torr) plasma processing applications including diamond film deposition Cl]. One way to create a microwave excited plasma is utilizing a single mode microwave resonance cavity which couples microwave energy into the discharge via joule heating [2]. The resonant electromagnetic field energy is imparted to the electron gas, where the energy is then transferred to ions and neutrals by elastic and inelastic collisions. In order to understand the behaviour of the plasma inside a cavity reactor source, the microwave fields, gas transport and energy transport, as well as their relations, need to be investigated carefully. The coupling among these is complex and analytic models are hard to apply to handle this problem. Thus, a systematic numerical model is required to solve the plasma behaviour and the exciting microwave fields selfconsistently. One type of numerical model for plasma discharge simulation is the fluid (continuum) plasma model, which solves the moments of the Boltzmann equation, including Elsevier Science S.A. SSDI 0925-9635(95)00291-X
modelling
the continuity equation, momentum transport equation and energy transport equation. This technique has been widely used to investigate the behaviour of the discharges inside parallel plate r.f. plasma reactors [ 3,4]. The electric fields used for plasma excitation inside the reactor are solved by Poisson’s equation with the boundary conditions assigned on the electrodes. In terms of microwave discharge simulations for diamond deposition, the microwave electric field distribution has been assumed uniform [ 51. In this paper, a numerical model is developed which self-consistently solves the electromagnetic fields as well as the excitation mechanism and characteristics of the diamond deposition discharges inside a cylindrical microwave cavity plasma reactor. The electromagnetic fields are solved by using a finite difference time domain (FDTD) [6] method. Previous work successfully adapted the FDTD method to solve the electromagnetic field distribution inside a compact ECR source and also compared the FDTD method solution directly to analytical solutions on simple cavity geometries [7]. Moreover, the FDTD method has been applied to study the electromagnetic fields inside a microwave cavity
1146
W. Tan, T.A. Grotjohn/Dinmond
und Related Mat&u/s
plasma reactor used for diamond thin film deposition by assuming a uniform plasma load inside the cavity [S]. In this paper this uniformity assumption is removed, and the characteristics of the discharge are simulated by a fluid plasma model which solves the electron and ion continuity equations, the electron energy balance equation, and the Poisson equation. A final electromagnetic field and plasma discharge solution is obtained by iteratively solving the electromagnetic field model and the fluid plasma model. The details of coupling the FDTD electromagnetic field model and the fluid plasma model are provide in the following sections. In particular, a TMO13 mode, 17-78 cm i.d. microwave cavity ,plasma reactor loaded with a hydrogen discharge is simulated using this numerical model. The TM,,, mode is an azimuthally symmetric mode. The hydrogen plasma example was chosen because diamond film deposition processes often consist of high percentages of hydrogen in the discharge. And, as shown by Koemtzopoulos et al. [S], adding small percentages (e.g. 1%) of methane to a hydrogen discharge has only a minimal effect on the electron energy distribution function. The method used to excite electromagnetic waves inside the cavity and the boundary condition assignments for both the electromagnetic field model and fluid plasma model will be discussed. The methods to evaluate the cavity quality factor Q will be presented. The spatial electric field patterns, power absorption patterns and quality factor of the cavity loaded with a hydrogen discharge are investigated in the moderate pressure range of 40-60 Torr. The physical behaviour of the hydrogen discharge, such as plasma density, electron temperature and plasma potential, are also simulated and analysed for various input conditions. The simulated results are compared with experimental data.
2. Model description The numerical model developed in this study includes an FDTD electromagnetic field model and a fluid plasma model. These numerical models and the coupling of these two models are described in this section. 2.1. FDTD electrtzmagnetic$eId model
4 (1995) 1145-1154
density. The FDTD method is formulated by discretizing Eq. (1) with a centred difference approximation in both the time and space domains. In this study, Eq. (1) is discretized in two-dimensional (P and z) cylindrical coordinates [S]. One of the major motivations to solve Maxwell’s equations is to investigate the power dissipation inside the microwave plasma source due to the presence of plasma discharges. The power dissipation density, Pa&t), with a power absorbing load (such as a discharge) present is Pa&t)
= J(r,WW)
(2)
The current density, J, which is induced by the microwave fields, can be determined by solving the momentum transport equation ofelectrons. For example, if ion motion is ignored, the momentum transport equation of electrons can be written as m, $ v(r,t) = - eE(r,t) - mev&r,t)Y(r,t)
(3)
which is also called the Langevin equation [9], and the current density as J (r,t) = - @&Wr,r)
(4)
where v is the average electron velocity, e is the electron charge, m, is the electron mass, veff is the effective collision frequency [!I] and n, is the electron density. Eqs. (3) and (4) are solved by the finite difference method [S]. From the above equations, the spatially dependent current density J is determined by the spatially dependent discharge characteristics, such as electron density n, and effective collision frequency veff. Moreover, the effective collision frequency which describes the momentum transfer due to electron-neutral elastic collisions is a function of the electron energy distribution (electron temperature). For a complete microwave plasma excitation solution, the electron density and electron temperature must be solved self-consistently by using a fluid plasma model. 2.2. Fluid plasma mode!
The electromagnetic fields inside the cavity can be described by the time-dependent Maxwell’s equations, which are written as
0) where E and H are the electric and magnetic fields, p is the permeability, E is the permittivity and J is the current
The behaviour of the charged particles in a weakly ionized gas can be described by the particle, momentum and energy balance equations for electrons and ions, which are obtained from the moments of the Boltzmann equation. Moreover, these equations are often combined with Poisson’s equation to provide a self-consistent space charge field. Fur the moderate pressures (I150 Torr) of this study, the momentum balance equations of electrons and ions are given by Eqs. (8) and (9) according to Ref. [lo). In the steady state, these equa-
W. Tan, T.A. Grotjohn/Diamond and Related Materials 4 (1995j 1145-1154
tions are written as [ 111: (5) V.J, = nennkion - apine
(6)
V.Ji = nennkion - a,nin,
(7)
- D,Vn,
(8)
J,= -n,p,E
Ji = nil”iE - Di V ni V.q, = - 4 J,.E -2
-
(9)
nenn(Ei&on + Eextkxt + Edd
2
i ksl&in,
qe= ; kgzJe-
c
; k,D,n,
(10)
>
UT,
(11)
Eq. (5) is the Poisson equation, where Y is the potential (E = - Vlu). Eqs. (6) and (7) are the electron and ion continuity equations, respectively. Eqs. (8) and (9) are in essence the momentum balances for electrons and ions, respectively. Eq. (10) is the electron energy balance equation, with the total electron energy flux given by Eq. (1 1), where the thermal energy is assumed to be greater than the kinetic energy. In the above equations, n, and ni are the electron and ion densities, respectively; J, and Ji are the electron and ion fluxes, respectively; T, is the electron temperature; qe is the electron energy flux; kionr kext and kdis are the inelastic rate constants, with threshold energies &ion,E,,, and Edis for ionization, excitation and dissociation, respectively; and CI, is the recombination rate constant. Note that elastic losses have been included in Eq. (lo), where \I,,, is the electron-neutral momentum transfer frequency and m, is the species mass. Moreover, D,,i and ~~,i are the electron and ion diffusivities and mobilities, respectively. The fluid plasma model describes the characteristics of discharges by solving Eqs. (5)-( 11). Finite difference methods and staggered mesh techniques [12] are used to discretize these equations in cylindrical coordinates. Since the resonant electromagnetic modes of the microwave plasma sources used for diamond thin film deposition in this study are normally TMoln mode, the electric field spatial distributions are C$symmetric in nature and the E+ component is zero. Therefore the discharge behaviour can be assumed to be 4 symmetric and the discretization of the equations reduced to a two-dimensional problem. Thus, the simulation region remains in the r-z plane only. Note that the heating term qJ*E in Eq. (10) primarily represents the microwave power absorbed by the plasma, which is determined by the microwave
1147
electric field and the microwave induced current density using the FDTD electromagnetic field model described earlier. Therefore, for simplicity, this heating term given in Eq. (2) as Pabs can be viewed as an input parameter for the electron energy balance equation. The flow chart of the microwave cavity plasma simulation is shown in Fig. 1. The solution of the discretized fluid plasma equations, Eqs. (5)-( 1l), on a grid of N nodes involves a total of 4N unknowns for the fluid model since at each point the unknowns are Y, n,, ni and T,. To solve the system of 4N non-linear type discretized equations, both direct method (Newton’s method) and iterative method techniques are used. In particular, the Poisson equation, electron continuity equation and ion continuity equation are tightly coupled and solved by Newton’s method. The unknowns solved at each grid point are the electron density, ion density and plasma potential, and the electron and ion flux (steady state d.c. flux) can be determined by Eqs. (8) and (9) using the solution of these unknowns. Then, the electron density and flux are coupled into the discretized electron energy balance equation, Eq. (lo), to solve the electron temperature. The calculated electron temper-
Solve electron energy balance equation to determine electron temperature(steady
state)
Solve Poisson equation and continuity equations to determine particle densities, flux and plasma potential (steady state) by Newton’s method
I
Solve Maxwell’s equation by FDTD model to determine the power absorption
‘ig. 1. Flow chart for microwave
cavity plasma
reactor
simulation.
1148
W. Tun, T.A. GrotjohnJDiamond and Related Materials 4 (1995) 1145-1154
ature is then coupled back to the continuity equations and Poisson equation to modify the reaction rates and parameters and update the electron density, ion density and plasma potential. The final stable solution of density and temperature is achieved by iteratively solving the continuity/Poisson equations and the energy balance equation. As described above, an input term in the discretized electron energy balance equation is the absorbed microwave power density Pabs, which comes from the FDTD electromagnetic field model. The FDTD model which solves the electromagnetic fields inside the discharges provides output information, such as the microwave power absorption of the discharges, by simultaneously solving the microwave electron momentum transport equation, Eq. (3). The discharge characteristics, such as the plasma density and electron temperature, are determined in the steady state under this power absorption condition. The discharge characteristics information is then coupled to the FDTD model to modify the plasma conductivity and calculate a new discharge power absorption. Therefore, in this iterative manner, the power absorbed by the discharges will converge to a stable value, and the electromagnetic fields inside the reactor source and the plasma characteristics are solved selfconsistently.
3. Implementation for a microwave cavity plasma reactor
which the microwave energy is coupled into the cavity through a coaxial input probe. The specific resonant mode is determined by the geometrical size of the cylindrical cavity, and it can be adjusted by the movable sliding short. In this study, the height of the cavity is adjusted for TM0i3 mode excitation with a 2.45 GHz input microwave frequency. The coaxial input probe and movable short position are adjusted for maximum power transfer to the plasma load, as indicated by minimum reflected power from the cavity. A hydrogendominated discharge is excited and sustained inside the quartz dome by the input microwave power. The volume and the temperature of the discharge are dependent on the input power and pressure. The substrate and substrate holder are placed on the top of a quartz tube which is used to adjust the height of the substrate. A metal screen is placed at the bottom of the cavity to prevent the microwave field from propagating into the downward region. The simulation structure for the FDTD model and the plasma fluid model is shown in Fig. 3. The sliding short, the base plate, the substrate holder, the coaxial probe (both the inner and outer conductor) and the cavity sidewalls form the boundary of the FDTD simulation region. The boundary conditions for the electric fields on these surfaces are that only the normal components of the electric fields exist and the tangential electric fields on these surfaces are zero. This is based on assuming that these boundaries are formed by perfect conductors. At the input end of the coaxial probe is an
In this study, the numerical model is applied to simulate a microwave cavity plasma reactor used for diamond thin film deposition, as shown in Fig. 2. The cavity is a cylindrical, 17.78 cm i.d. microwave cavity in
Wave Excitation /Plane (TEM Wave)
Perfect Conductor
kFig. 2. Configuration of a microwave cavity plasma reactor diamond film deposition. The cavity diameter is 17.78 cm.
used for
Fig. 3. Simulation structure reactor with quartz tube.
cross-section
Quartz Tube (q~3.75) for microwave
cavity plasma
W. Tan, T.A. Grotjohn/Diamond
and Related Materials 4 (1995) 1145-1154
open boundary where the region in which the field has to be computed is unbounded. At this open boundary a truncation method is used to prevent any artificial reflection of outgoing waves. This non-reflecting boundary condition has been investigated extensively in solving electromagnetic wave scattering problems [ 131. Basically, this boundary condition allows electromagnetic waves arriving at the open boundary of the coaxial probe to leave the simulation region without any artificial reflections at the termination of the grid structure. The non-reflecting boundary condition which is applied in this model can be written as: (E:+‘(i,k,)
1149
+ ~C~%41’21- WlNW(w~1’21 + Ad)
AzrAr2rr) are approximately equal. For the FDTD electromagnetic field model, the number of grids used in the simulation is 30 x 83 in the r- and z-directions respectively, and the time step is 0.5 ps. The grid spacing for the FDTD model is uniform in both r- and zdirections. Since the grid structure spacing for the FDTD model is different from that of the fluid plasma model, a linear interpolation technique is used to interpolate the absorbed power to each grid point of the plasma model, and the plasma density and electron temperature to each grid point of the FDTD model. For the partially ionized and partially dissociated H2 plasma discharge simulation, the major particle interaction processes are the electron-H, molecule inelastic collision, electron-H, molecule elastic collision and electron-hydrogen ion recombination. The electron-H, inelastic collisions include the H2 molecule ionization, excitation and dissociation processes. The rate coefficients in Eqs. (6), (7) and ( 10) for these collision processes can be expressed using the Arrhenius relationship [ 111 as
x [Ez+‘(i,(k, - 1)) - Ez(i,k,)]
kion= &n exp ( - aion/KBT,)
= E:(i,(k,-
1))
+ KWW’21
- WlK~~/W”“1
x [E:“(i,(k,-
1)) - E:(i,k,)]
+ AZ)) (12)
and E$+‘(i,k,,) = E$(i,(ko - 1))
(13)
where i denotes the grid location in the r-direction and k, is the grid terminating point in the z-direction. Additionally, the quartz dome and the quartz tube, where the dielectric constant is equal to 3.75&,, are included in this model. The boundaries of the fluid plasma simulation are the quartz dome and the substrate holder. The diameter of the substrate holder is assumed to be the same as that of the quartz tube. The boundary conditions of the fluid plasma model are, at the substrate: n,=ni=O T,= T,
(14)
Y=O
kext= A,,, exp ( - GJGT,) his = &is exp ( - EsdKBTe)
where &ion, eext and ldis are the threshold energies for H,
molecule ionization, excitation and dissociation and Ai,,, A,,, and Adis are the pre-exponential factors, which are obtained by approximating the rate constant data at low electron temperatures to these relationships [11,15]. For simplicity, only the reactions with higher rate coefficients are considered in this study for H, discharges. The types of inelastic collisions and their corresponding rate parameters are summarized in Table 1. The collision frequency for electron-H, molecule momentum transfer is relatively independent of the electron temperature, so it can be written as [9]
and at the quartz wall: n,=n,=O
v,(H2)=
1.44 x 1012 x
Pressure (Torr) (17)
T,(K) Table 1 H, reaction
The boundary conditions at the quartz wall and the substrate wall assume that recombination occurs when electrons and ions reach the walls and that thermal equilibrium exists at the walls. The substrate and the substrate holder can be either grounded or floating. In the simulations of this paper the substrate potential is assumed grounded. The potential is assumed to be at the floating potential Yr on the quartz disk [ 141. The grid structure used in the plasma simulation is in the Y-and z-directions. In the z-direction a uniform grid spacing is used, and in the r-direction the grids are constructed such that all the unit cell areas (or volumes,
(16)
rates
Reaction
Expression
Ionization
e+H*-+e+H:+e H; + H,+H:
+ H
Threshold energy
Pre-exponential factor or rate coefficient (m3 s-j)
15.4 eV OeV
1.0 X 10-l” H+ is the diminant ion species
Excitation
e+H,-+Hf+e
12.0 eV
6.5 x lO_”
Dissociation
e+H,+e+H+H
10.0 eV
1.0 x lo-‘4
e + ion-neutral
OeV
1.0 x lo-‘4
Recombination
W. Tan, T.A. GrotjohnjDiamond and Related Materials 4 (1995) 1145-1154
1150
where T, is the neutral temperature, which can be represented by the translational temperature of H2 gas. Since v, is relatively independent of the electron temperature, the effective collision frequency veff is set equal to V It should be noted that the only neutral species ctnsidered in the ion and electron simulations was the Hz species and that the dominant ionic species in the plasma is Hz [5,11-J. Hence these simulations assume a low hydrogen dissociation percentage.
showed no significant variation along the #-direction as expected for the TM,,s mode. Once the electromagnetic field strength and mode is determined inside the cavity, the cavity quality factor Q of the plasma discharge loaded cavity reactor can be found by calculating the electromagnetic energy U stored in the cavity and the power absorbed in the cavity Pa,,*. The energy stored in the cavity U is expressed as U
1E12dV
=E
s 4. Experimental results
(18)
and the Q factor is given by
In order to verify the cavity resonant mode, the radial component of electric field strength at the cavity outer wall was measured by calibrated micro-coax electrical probes inserted through the cavity wall versus z and 4 [ 11. The centre conductor of the probe is inserted 2 mm beyond the inside surface of the cavity side walls, and the other end is connected to a microwave power meter. Since the probe power reading is proportional to the square of the r.m.s. electric field strength normal to the inside cavity, the cavity resonant mode and the cavity stored energy can be determined using this technique [l]. Probe experiments were done by measuring the electric field strength at various z positions along the cylindrical cavity wall of the microwave plasma reactor during diamond thin film deposition. The microwave excitation frequency was 2.45 GHz, and the cavity height was adjusted to a position near 20.4 cm for a minimum reflected power condition. The result is shown in Fig. 4, where the resonant condition for the microwave cavity during a deposition process matches closely to the ideal TM0r3 mode. Additionally, the electric field was probed circumferential around the cavity. The field strength 2.5,
* Experiment
Theory
I
Q=E
where o is the excitation frequency. The measured Q value at a pressure of 50 Torr and an input microwave power of approximately 1500 W was 100. A Q value of 100 gives an indication as to why the microwave energy pattern with the plasma load present is similar to the pattern of an empty resonant cavity structure.
5. Simulation results The microwave cavity plasma reactor loaded with a H2 discharge for diamond film deposition was simulated and investigated by the self-consistent numerical model which has been described in the previous sections. The input parameters for this model include the pressure and input microwave power. A set of experimentally determined empirical equations developed by King [16] are used to establish the neutral temperature for this simulation. These empirical equations were obtained by doing discharge diagnostics across a parameter space including pressure and input microwave power variation. The empirical equations used for the translational temperature of H, gas and the discharge volume are [ 161: Translational
*
*
(19)
temperature (K)
= 228.6 + 374.3 x Incident power (kW) + 16.5 x Pressure (Torr) k 94.2
(20)
Plasma volume (cm3) = 449.7 + 116.2 x Incident power (kW) - 18.1 x Pressure (Torr) + 57.1 x [Incident power (kW)]* I
+ 0.25 x [Pressure (Torr)12 - 5.4 0.
0
a
2
* 4
6
8
10
12
14
16
18
20
Cavity Height (cm) Fig. 4. Microwave electric field strength (radial directed) vs. cavity height for microwave cavity plasma reactor during the diamond deposition process.
x Pressure (Torr) x Incident power (kW) + 15.4 (21) The translational temperature was determined using Doppler broadening measurements of the H atom emission and the plasma discharge volume was determined
W. Tan, T.A. GrotjohnjDiamond and Related Materials 4 (1995) 1145-1124
using photographs of the plasma. Eqs. (20) and (21) are valid for the pressure range from 35 to 6.5Torr and incident microwave power range from 1.4 to 2.6 kW. The plasma discharge volume serves as a boundary for the portion of the plasma simulation region where the neutral temperature (T,) is high. Inside the discharge volume region, the translational temperature given by Eq. (20) is the neutral temperature. Outside the volume, the neutral temperature is equal to the temperature assigned on the boundaries (T, = 1000 K). Between these two regions, a linear temperature change profile is used to prevent an abrupt change of neutral temperature. This is done by assigning the boundary temperature (T, = 1000 K) to the region several grids, namely three grids, away from the plasma volume in both the r- and z-directions. The temperature assigned on the grids between these two regions are linearly interpolated from the neutral temperature given by Eq. (20) and boundary temperature (T, = 1000 K) as shown in Fig. 5. For the FDTD electromagnetic field model, the technique used to excite the electromagnetic field in this study is to select grid points on a cross-section plane of the coaxial probe as source points and assign the timevarying electric field component at these points based at theoretical transverse electromagnetic (TEM) wave solutions in a coaxial structure, as shown in Fig. 4. The electromagnetic wave then propagates down into the cavity region where power is absorbed by the discharge. Any reflected electromagnetic wave will propagate to the end of the open boundary of the input power probe and be terminated by the non-reflecting boundary. The simulation results for the microwave electric field distribution are shown in Figs. 6 and 7 for the E, and E, components respectively. The input condition for this simulation is 50 Torr in pressure and 1500 W input microwave power. The electric field patterns basically follow the TMo13 mode electromagnetic field distribution. The effect of the electric field caused by the presence of discharge can be observed. It reduces the amplitude
Neutral temperature 18CUl
1600
14ca 0
1200
1090 0
Fig. 5. Neutral temperature region.
distribution
in the plasma simulation
Xi04
Sliding Short
1151
Quartz Dome
0 0
Plasma Discharge
\
Substrate
Fig. 6. E, field distribution at a pressure of 50Torr and an input microwave power of 1500 W.
XiOP
Sliding Short
Plasma Discharge
Quartz Dome
Substrate
Fig. 7. E, field distribution.
of the radial directed electric field component, E,, in the plasma discharge region since electromagnetic power is absorbed in the discharge region. The discontinuity of the electric field resulting from the difference of dielectric constant between the quartz and air can also be seen. It can also be observed that in the plasma simulation region the E, field has the highest magnitude just above the substrate. The conductivity of the plasma can be determined using the solved density and effective collision frequency. Fig. 8 shows that the maximum conductivity of the discharge is around 0.4 mho m. Fig. 9 shows the power absorption pattern versus I and z and indicates that the microwave power is highly absorbed in the central portion of the plasma region close to the substrate. The simulated electric field strength in Figs. 6 and 7, as well as the field measured in Section 4, indicate that the microwave electric field strengths are of the order of lo4 V m-‘. If this electric field strength is assumed to
W. Tun, T.A. Grotjohn/Diumond and Reluted Muterials 4 (1995) 1145-1154
1152
2
4
51
> 0 0
2 0 0
(a)
r(m)
(4
em)
10 6
(b)
r(m)
10 6
10 6
Substrate
Plasma Discharge
Fig. 10. Simulation results of a microwave cavity plasma reactor for pressure = 50 Torr and input microwave power = 1500 W. (a) Electron temperature; (b) plasma potential; (c) plasma density. Note that z =0 is located at the substrate position for the plasma simulations. (This same z location is at 4.5 cm for the electromagnetic field simulations.)
Fig. 8. Conductivity vs. r and z.
2 81 0 30
10
x lop' 10
/
0 0
Su&rate
Plasma Discharge
?
E 7
J
Fig. 9. Absorbed microwave power vs. r and z.
contribute to the heating of the substrate through direct electromagnetic energy absorption, then the amount of this direct electromagnetic heating can be estimated. Specifically, if a conductivity is assumed for the substrate holder similar to that of graphite, the absorbed power in the substrate and substrate holder is calculated to be less than 1% of the power entering the plasma source. Hence in the resonant cavity microwave system described in this paper the microwave electric fields are not strong enough to contribute substantially to the substrate heating; rather, the substrate is heated by the plasma. The simulation result of the plasma characteristics of a H2 discharge at 50 Torr and 1500 W input power is shown in Fig. 10. The plasma density is highest near the centre (r=O) of the discharge. The electron temperature is also highest at the centre of the discharge. The spatial variations of ionization and recombination rate are shown in Fig. 11. The ionization rate is higher near the substrate due to the electron temperature being higher
5
x
0 0 0
Fig. 11. Spatial variation of ionization rate (top) and recombination rate (bottom) for pressure = 50 Torr and input microwave power = 1500 w.
at that region. The recombination rate profile basically follows the electron density variation. The characteristics of Hz discharges and cavity quality factor Q were analysed for various pressures and absorbed powers by this numerical model. The neutral temperature and plasma volume determined by Eqs. (20) and (21) are the input data for simulations and their variations for different incident power and pressure are shown Fig. 12. The simulation results of average electron temperature vs. absorbed microwave power are shown in Fig. 13. These values are calculated by averaging the electron temperature across the plasma volume region. The average electron temperature is about 1.5 eV. The electron temperature increases with higher absorbed microwave power. This occurs because the measured neutral gas temperature increases with microwave
W. Tan, T.A. Grotjohn/Diamond and Related Materials 4 (1995) 1145-1154
1153
1.50
*: 40 Tom 0: 50 Torr x: 60 TOIT
x
I
x
40 Torr 50 Torr 0: 60 Torr *:
140-
0
x:
Y
130120-
2
G
o: 50
Torr x: 60 Tom
I_.*
,i 9
x
loo
B a
50
.z
0 1600
1500
x 1700
--_+__L 1800
ioco 1900
21W
2cwJ
Fig. 12. Experimentally measured neutral temperature volume vs. incident power for various pressures.
1500
2ow
J
2500
Power absorbed (watt)
and plasma
Fig. 14. Simulated Q value vs. absorbed power for three different pressures. The experimentally measured Q at a pressure of 50 Torr and an input microwave power of 1500 W is also indicated.
nique described in Section 4. The Q value determined experimentally is about 100, which is close to the simulation result (Q= 107). This close agreement is an indication of the validity of the model.
*: 40 Torr o: 50 Torr x: 60 Torr II
*
6. Summary
0
0.9
0
x
1600
2000
x
1
I
1200
1400
1603
2200
2400
Absorbed Power (Watt) Fig. 13. Simulated electron temperature vs. absorbed power for three different pressures.
power, as indicated in Eq. (20). The higher neutral gas temperature leads to a reduced neutral density and, hence, to a reduced electron-neutral collision rate at a given electron temperature. Also, the electron temperature drops with increases of pressure. This can be understood since higher pressure causes a higher collision frequency for electrons with neutrals. Therefore electrons more easily transfer their energy to neutral particles and reduce their kinetic energy. The results for Q vs. absorbed microwave power are shown in Fig. 14. It shows that when the input microwave power is increased the Q factor decreases. The pressure has only a small effect on the Q factor. At the 50 Torr pressure and 1500 microwave power condition, the Q factor of the microwave cavity plasma reactor with an Hz discharge was experimental determined by using the tech-
A self-consistent numerical model has been developed to simulate the electromagnetic excitation of discharges inside microwave cavity plasma reactors. This software includes an electromagnetic field model and a fluid plasma model. The microwave cavity plasma reactors simulated by this numerical software were loaded with a H, discharge and were used for diamond thin film deposition. The complex reactor geometry and input power coupling probe structure were included in the simulation. The input parameters for the microwave cavity plasma reactor simulation included the pressure (40-60 Torr) and input microwave power. The electromagnetic behaviour of the discharge loaded resonant cavity, such as the electric field distributions, power absorption patterns, and cavity quality factor Q, were studied and analysed for various input conditions. The characteristics of H, discharges, including the plasma density, electron temperature and plasma potential, were also investigated and studied. The electromagnetic mode and cavity Q value of an Hz discharge loaded reactor were simulated and shown to be in agreement with measured experimental values. The numerical model presented in this study allows the calculation of the power absorption profile in microwave discharges, which has often been a difficulty in plasma simulations. The understanding of the interactions between the microwave fields and the plasma
1154
W. Tan, T. A. GrotjohnlDiamond and Related Materials 4 ( 1995) 1145-I 154
discharges allows the reactor design to be analysed for improvement and optimization of such quantities as uniformity of microwave power absorption, uniformity of the plasma species and kinetic energy of the plasma species. Additionally, the simulation and understanding of the electromagnetic excitation of the discharge is expected to be important for microwave discharge diamond deposition machine control.
c31 J.P. Boeuf, Phys. Rev. A, 36 (1987) 2782. c41 F. Young and C-H. Wu, IEEE Trans. Plasma Sci., 21 (1993) 312. c51 C.R. Koemtzopoulos, D.J. Economou and R. Pollard, Diamond Relat. Mater., 2 (1993) 25. IL61K.S. Yee, IEEE Trans. Antennas Propagat., I4 (1966) 302. c71 T.A. Grotjohn et al., Rev. Sci. Instrum., 65 (1994) 1761. C81 W. Tan and T.A. Grotjohn, .I. Vat. Sci. Technol., 12 (1994) 1216. c91 B.E. Cherrington, Gaseous Electronic and Gas Lasers, Pergamon,
New York, 1979. Cl01 J.D.P. Passchier and W.J. Goeheer, .I. Appl. Phys., 74 (1993) 3744. Cl11 M. Surendra et al., J. Appl. Phys., 71 (1992) 5189. Cl21 MS. Barnes, T.J. Cotler and M.E. Elta, J. Camp. Phys., 77
(1988) 53.
References
Cl31 G. Mur, IEEE Trans. Electromagn. Compat., 23 (1981) 377. Cl41 D.M. Manos and D.L. Flamm, Plasma Etching, Academic Press,
San Diego, CA, 1989. J. Zhang, B. Huang, D.K. Reinhard and J. Asmussen, J. Vat. Sci. Technol. A, 8 (1990) 2124. [2] J. Asmussen, J. Vat. Sci. Technol. A, 7 (1989) 883. [l]
Cl51 R.K. Janev et al., Elementary Processes in Hydrogen-Helium Plasmas, Springer, Berlin, 1987. Cl61 G. King, Ph.D. Dissertation, Michigan State University, 1994.