M. J. Goeckner, J. Goree, and T. E. Sheridan

Department of Physics and Astronomy

The University of Iowa

Iowa City, Iowa 52242

A simulation of ion dynamics in a planar magnetron discharge is performed using separate three-dimensional Monte Carlo codes for the electrons and ions. First, to predict the ionization sites, the orbits of energetic electrons are simulated for prescribed dc electric and magnetic fields, subject to collision with neutrals at random intervals. In the second code, the predicted sites are used as the starting positions of ion trajectories. The ion trajectories are followed taking into account collisions with neutrals, turbulent electric fields, and the dc fields. We report results for ion impact on the cathode and substrate anode surfaces (energy, angle, and spatial distribution) and ion parameters in the plasma (density, drift velocity, random energy, and transit time). To test these results, we compare them to several previously reported experiments, and in most cases find good agreement. These simulation methods not only are useful for gaining an understanding of magnetron plasma operation, but may also aid in designing magnetrons.

I. INTRODUCTION

Sputtering magnetrons are used for thin film deposition and sputter etching. [1] In these devices an electric sheath and an external magnetic field are configured to trap electrons. [2,3] Typically the magnetic field strength is weak enough that the electrons are magnetized, but the ions, due to their larger mass, are unmagnetized. [3] Most of the ions hit the target cathode, having gained several hundred eV of energy in the sheath, sputtering material from the surface [2,4] and causing secondary electron emission. [5] A small fraction of the ions travel toward the substrate anode, where they may become embedded in the deposited thin film of sputtered atoms. [2] It is worthwhile to develop models of magnetrons not only to gain an understanding of their operation, but also to aid in their design.

The greatest recent progress has been in the development of electron models. Wendt [6] used a Hamiltonian approach to calculate the location of the electron trap. We [3] then presented a model of electron transport which assumes that energetic electrons are scattered from the trap by collisions with neutrals. This electron model has been implemented with a Monte Carlo code [3,7] to predict the spatial distribution of ionization in the plasma; several tests of this code have shown good agreement with experimental measurements [3,8,9] from two different magnetron designs. Guimaraes, Almeida, and Bretagne [10] have reported a Boltzmann equation model for electrons in a planar magnetron that differs from the Monte Carlo simulation by taking into account lower energy electrons and their collisions but with a schematic treatment of electron orbits and their losses.

In this paper we report a Monte Carlo simulation of ion dynamics in our magnetron. The simulation includes an electron code to predict the ionization locations and then a separate code to follow ion trajectories. For electrons, we use a particle simulation rather than a continuum model because their complicated trajectories dominate magnetron physics. For ions, the trajectories are calculated assuming collisions with neutrals, dc electric and magnetic fields, plus a turbulent electric field. In the next section, we describe the magnetron that we modeled and detail the assumptions made in the simulation. In Sec. III we present simulation predictions of ion impact location, energy, and impact angle on the electrodes. Additionally, we compare and find good agreement between simulation results and earlier experimental measurements of ion density, drift velocity, random energy, [11] and transit time. [12] We also find that the dc electric field and ion-neutral collisions strongly influence the ion trajectories while the influence of turbulence is weak.

II. SIMULATION

A. Magnetron model parameters

For the work presented here, we have performed simulations of the cylindrically-symmetric planar magnetron described in Ref. 11. In Fig. 1 we sketch this device, indicating the location of the electron trap, and defining the cylindrical coordinate system (r,[[theta]],z). These and other variables are listed in Table I. A set of permanent magnets located behind the copper cathode provide a magnetic field with a purely radial component of 245 G on the surface at r = 1.7 cm, which is the same location as the deepest point in the etch track erosion profile. [3]

The plasma in this device has been characterized thoroughly in a series of previously reported experiments. Using a Langmuir probe we have measured the plasma potential, shown in Fig. 2, in the presheath from z = 0.5 to 4.0 cm. [11] These measurements were made for a -400 Vdc, 1 Pa argon discharge, which had an electron density of 2.5 x 10¹⁰ cm^-3 and temperature of 3 eV, measured at r =1.7 cm, z =1.0 cm. To reduce sputter contamination of the probe and vacuum windows, this discharge was operated at a discharge current of 0.1 A. This corresponds to an average current density at the cathode of 9 mA/cm², which is 10 to 50 times less than might be used for efficient film deposition. Sheridan and Goree measured the low-frequency turbulence in the electron trap and found that the rms density fluctuation was 2.7 % with a corresponding rms electric field of 0.46 V/cm.[12]

We have also reported time-averaged laser-induced fluorescence (LIF) measurements of ion velocity distributions at z = 1 cm. [11] We found that there was no measurable average drift velocity parallel to the cathode, that the ion density is peaked near r = 1.7 cm, and that the random ion energy parallel to the cathode ranged from room temperature at the edge of the discharge to 0.32 eV at the center.

B. Model

To model ion trajectories in this magnetron, we make a number of assumptions based upon experimental measurements. We assume that the argon ions are born predominately in the trap region at sites which can be accurately predicted by an electron Monte Carlo simulation. [3] The ions are accelerated from these sites by the dc and turbulent electric fields, while possibly undergoing collisions with neutrals. Below, we evaluate a number of physical processes and determine which to include in our simulation.

The dc electric field, a major influence on the ion's trajectory, is included in the simulation. In the sheath, where ions gain most of their energy, the dc electric field is in excess of 800 V/cm. In the presheath, our previously reported measurements of electric potential contours [11] , shown in Fig. 2, indicate that the dc electric field is typically 2 V/cm along the r and z directions. Coupled with the spatial distribution of the ionization events, this process could bring about the random energy in the radial direction that was measured in the LIF experiment. However, the LIF data also showed a large random energy in the [[theta]] direction where, due to symmetry, there is no potential drop. This azimuthal random energy indicates that an additional process must be occurring. Two possibilities are turbulence and collisions.

The low-frequency turbulent electric field is included in the ion simulation. In the electron trap region of the presheath, this field is weaker than the dc field by a factor of four. [12] However, the large random energy in the azimuthal direction might be due to the turbulent field.

Collisions with neutral gas atoms are also included because they can scatter an ion's trajectory. Both elastic and charge exchange collisions, which have comparable cross sections, may be important. The significance of elastic scattering can be seen as follows. Based on Langmuir probe measurements of the density in the discharge described above, Sheridan and Goree reported that the ion confinement time is 0.8 us. [12] This time can be compared to the 2.3 us time between elastic collisions, [13] indicating that an ion has a 35% chance of having an elastic collision between its creation and its collection on the cathode.

The magnetic field does not significantly influence ion trajectories, [3] however we retain it in the ion simulation for completeness.

There are a number of processes that we do not take into account in the simulation. We assume that the neutral gas is at room temperature throughout the discharge and ignore any hot neutrals. [14] We also ignore ion-ion collisions because the energy transfer rate is much slower than the rate of ion loss. [15] Finally, we neglect collisions with the sputtered atoms because their number density is less than the neutral gas density by 3 to 4 orders of magnitude.

C. Numerical implementation

To implement our physical assumptions, we first use an electron Monte Carlo simulation to establish the ionization sites and then an ion simulation to follow the ion trajectories. In these simulations, we assume that the electrons and ions move in the electric and magnetic fields, subject to collisions at random intervals. To do this, we must provide expressions for the fields, follow the orbits based on the equation of motion, and provide accurate scattering cross sections. In this section, we describe the simulations and list the approximations that they entail.

1. Initial Conditions

In the electron simulation, which is documented in Ref. [3], only the electrons starting at rest on the cathode are followed. This electron code was used previously to predict the spatial distribution of ionization in the plasma, and it showed good agreement with experimental data. [3,7]

In the ion simulation, the ions are assumed to start from the ionization sites predicted by the electron simulation. The initial ion velocities are randomly chosen from an isotropic room-temperature Maxwellian distribution.

2. Integrating the Orbits

Particle orbits x(t) are computed by integrating the equation of motion

(1)

using a fourth-order Runge-Kutta integrator [16] and the prescribed electric and magnetic fields. For most of the trajectory a fixed time step [[Delta]]t, which is smaller than the inverse collision frequency, is chosen so that energy is conserved. In the ion simulation, when a particle approaches the cathode, [[Delta]]t is reduced to provide accurate predictions of the impact energy.

Particle orbits are stopped when they hit the cathode (z = 0), the substrate anode (z = 4 cm), or escape out the side (r = 4 cm). In the electron simulation, the negligible fraction of orbits that are still confined after 5 us are also stopped. We then begin a new particle in the Monte Carlo simulation.

We prescribed dc electric and magnetic fields on a 40 x 40 grid in the r-z plane, and used bilinear interpolation between the grid points. The dc electric field is composed of separate presheath and sheath portions joined at z = 0.5 cm, where we made the electric potential continuous. The presheath was based on the experimental measurements shown in Fig. 2, while the cathode sheath was calculated from a power-law potential. [17] While this two-dimensional electric field does not include an anode sheath, it is more realistic then the one-dimensional model used in previous magnetron simulations. [3,7]

The two-dimensional magnetic field for our cylindrically symmetric planar magnetron was computed from Poisson's equation [3] based on the magnet configuration. [12]

A turbulent electric field included in the ion simulation is based on Langmuir probe measurements made by Sheridan and Goree[12]. These experimenters reported that the fluctuations are probably ion-acoustic waves that grew from a current-driven instability in the plasma. Using a transient recorder, sampling at 200 ns intervals, they recorded a digitized waveform of the electron density fluctuations at two locations in the plasma. We chose to use the data from the location where the fluctuations were larger, based on their rms value. This density fluctuation waveform was then scaled to produce an electric field waveform, adjusted in amplitude to have the same rms fluctuation value reported by Sheridan and Goree. In the ion code, each ion is given a different turbulent electric field chosen from data points in this waveform, and this is added to the dc field. The experimenters found that time scale of the fluctuations was much less than the transit time of an ion across the discharge, and so we assume in the code that the turbulent electric field is actually dc for a given ion. Because the experimental method was unable to determine the direction of the electric field, in our simulation we assume that it points in a random direction that is different for each ion.

3. Collisions

At each time step during an orbit, we ascertain whether a collision has taken place. This is done by generating a random number that is compared to the probability per unit time of a collision. When a collision takes place, the scattering process (elastic, etc.) is determined according to the cross sections. Then, in a manner consistent with the differential scattering cross section, the particle's velocity direction is randomly scattered and its energy is reduced. This method, as described in Ref. [3], takes into account elastic, excitation, and ionizing collisions in the electron simulation.

In the ion simulation, elastic and charge exchange collisions with the neutral gas atoms are included. We model the elastic scattering cross section [[sigma]]el in SI units using

[[sigma]]el = [c + d ln (v)]^f . (2)

We computed the coefficients c = - 2.27 x 10³⁶, d = 3.08 x 10³⁵, and f = - 0.52 by fitting (2) to Cramer's experimental measurements of the cross section, which were made for v >= 3800 m/s [18]. Below that velocity, we assume that the cross section is constant. To evaluate the charge exchange cross section, we use a similar expression reported by Iovitsu and Ionescu-Pallas [19] ,

[[sigma]]cx = ^{[a + b ln(v)]2}, (3)

who found that a = 1.51 x 10^-9 and b = -9.53 x 10^-11 in SI units for v >= 100 m/s. We assume a constant cross section for slower ions.

A hard sphere model is used when a collision takes place to decrement the ion's energy and to scatter its velocity direction. The neutral gas atom is assumed to be stationary. It is traditional to use the hard-sphere approximation for elastic collisions, while for charge exchange collisions at high energies it is customary to assume purely forward scattering. However, most of our collisions are at low energy, 0 to 5 eV, where we do not know of any experimental measurements of differential cross sections, so we use the hard-sphere model for all collisions. A differential cross section with a stronger forward peak would result in less randomization of ion velocities, while isotropic scattering would randomize them more.

We term the ions created by ionization collisions "primary ions" to distinguish them from the ions created by charge exchange collisions. An asterisk will denote results based solely on the primary ions.

D. Ion simulation diagnostics

1. General diagnostics

Ion impact density, angle, and energy are recorded when an ion hits the cathode or the substrate anode.

The mean ion transit time [[tau]]transit is calculated from the confinement time tk for each ion,

, (4)

where tk is the time from an ionization event until that ion, or a subsequent ion born by charge exchange, leaves the simulated region.

The spatially-resolved ion density in the plasma ni(r,z) is calculated by accumulating after each full time step the bilinearly-weighted ion location in 40 x 40 bins on the r-z plane.

2. Diagnostics related to LIF data

To compare the simulation results with our earlier experimental LIF data, [11] we count only the primary ions in the simulation. We used LIF to measure velocity distributions of metastable excited-state ions in the presheath. Both the excited-state ions and the ground-state ions are predominately produced by electron impact ionizations in this discharge. However, we suspect that charge exchange collisions involving metastable excited-state ions will most likely produce ground state ions and excited-state neutrals. Thus, we compare the metastable excited-state ions in the experiment only to the primary ions in the simulation. We assume that the cross sections are the same for all ions.

To compare to the reported LIF data, velocity distributions of the primary ions f^*(v) are collected from the simulation in 50-m/s wide linearly weighted bins. In the experiment, 1-cm diameter, 0.84-cm long cylindrical regions of the presheath were probed by the laser beam. [11] In the simulation, we record the velocities of primary ions when they are in toroidal regions approximating those cylindrical regions. From f^*(v) we calculate the primary ion density

, (5)

the average drift velocity

, (6)

and the random energy

, (7)

where Mi is the ion mass. The radial and azimuthal components of the velocities are collected separately so that two components of the drift velocity ([v^*r , [v^*[[theta]]) and the random energy (dE^*r , dE^*[[theta]]) can be calculated.

III. RESULTS

First, we performed an electron Monte Carlo simulation to establish the starting positions of the ions. We followed 2000 electrons in 50 ps time steps and recorded 24833 ionization events. Next we performed the ion simulation, starting an ion at the location of each ionization event. Ion orbits were integrated using 10 ns time steps until they escaped from the 4 x 4 cm boundary. If an ion approached the cathode, the time step was reduced to 1 fs just before impact. Both the electron and ion simulations made use of the two-dimensional magnetic and dc electric fields described in Sec. II.

To gauge the influence of collisions and turbulence, we repeated the ion simulation for three different cases:

(1) no collisions and no turbulent electric field,

(2) collisions but no turbulent electric field, and

(3) both collisions and a turbulent electric field.

We will emphasize case 3 because it includes the most physics and should be the most accurate of the simulations.

A. Surface impact results

Here we report the ion impact energy, angle and spatial distribution predicted by the ion simulation. These parameters determine the probability of sputtering and other ion-surface interactions. [2] They have been used elsewhere to model the transport of sputtered cathode particles. [20]

We find that most of the ions strike the cathode surface. In case 3, 99.5% of the ions (24699 out of 24833) hit the cathode while less than 0.5% (119) hit the substrate anode. The remainder were lost radially toward the vacuum vessel wall.

The cathode impact angle distribution, shown in Fig. 3, indicates that ions strike almost perpendicular to the cathode surface. This is a direct consequence of the large electric field in the sheath and the long mean free path between collisions. The mean impact angle, defined with respect to the surface normal, is 7.95 +/- 0.14 degrees. Ions striking the cathode with angles large than this have suffered at least one collision.

The cathode impact energy distribution, Fig. 4, shows that most of the ions have less than the maximum energy. The average impact energy is 293.02 +/- 1.34 eV at the cathode. Ions slower than the average are born in the sheath (84%) or have lost energy in a collision (43%). Thirty-three percent of the ions born in the sheath also suffer at least one collision. Ions entering the substrate anode sheath typically had less than 0.2 eV of energy. When considering these results, recall that the cathode bias is -400 V.

One may wish to know the angle and energy distributions for ions at impact on the cathode for various discharge currents and neutral pressures. Ions born in the presheath that do not undergo a collision in the sheath will strike the cathode with the full energy at nearly normal incidence. Ion collisions in the sheath account for the tail of the angular distribution (Fig. 3), and in combination with ionizations in the sheath, they also account for the tail of the energy distribution (Fig. 4). The number of ions born in the sheath will vary linearly with the sheath thickness d, which scales approximately as d [[proportional]] . The fraction of the ions that undergo collisions in the sheath will scale approximately as the product of the neutral density and the sheath thickness, Pd [[proportional]] P , provided that a large fraction of the ions do not undergo multiple collisions in the sheath. By using these scaling laws to adjust the proportions of the tails, one can construct the angle and energy distributions for any choice of neutral density and discharge current.

At the cathode, the radial profile of the impact density is strongly peaked at r = 1.7 cm, while at the substrate anode it is uniformly distributed. In Fig. 5 a), the radial profile of the impact density on the cathode is compared to the measured etch track depth, which is formed after many hours of use at various pressures and biases. Figure 5 b) shows the radial profile that would result from ions falling directly to the cathode without radial deflection. [3,7] Both profiles are peaked where the etch track is deepest. The use of the electron simulation unexpectedly offers a more accurate prediction of the radial profile then does the more complicated ion simulation.

B. Results in the plasma

In this section, we compare simulation results to experimental measurements of the magnetron plasma.

First, we compute [[tau]]transit using (4) and compare it to the experimental measurement made by Sheridan and Goree. [12] The distribution of ion lifetimes, shown in Fig. 6, yields [[tau]]transit = 0.259 us, with an uncertainty of +/- 0.004 us, which should be compared to the experimental value of 0.8 us. [12] This distribution has two distinct peaks, due to the birth locations of the ions. Those born in the approximately 5 to 6-mm wide sheath hit the cathode sooner. The minimum in the transit time at 0.9 us might be considered as the time required to cross the sheath. On average, ions born above 6 mm require more than 1 us to reach the cathode while ions born below 5.9 mm require less than 0.8 us.

Second, we compare several ion parameters predicted by the simulation to our LIF measurements. [11] In Fig. 7, we find good agreement between the radial profiles of ni^* and the measured ion density. Experimentally, the density is highest near the center of the electron trap, r = 1.7 cm, and lowest near the symmetry axis, r = 0. All three cases of the simulation agree with the experimental data. Including collisions (case 2) noticeably improves the agreement, while the addition of the turbulent electric field (case 3) has little effect.

We also find good agreement between [v^*[[theta]] and the measured azimuthal drift velocities, as shown in Fig. 8. This agreement is improved with the addition of collisions (case 2) while the turbulent electric field (case 3) has little influence. In contrast to the azimuthal direction, none of the simulations exhibited good agreement with the measured radial drift velocity, which was measured to be less than 93 m/s. [11] The simulations predict [v^*r > 850 m/s in the trap.

Comparing dE^* from the simulation and the random ion energy from the LIF experiment, we find good agreement for the radial but not the azimuthal velocity direction, as shown in Fig. 9. Experimental measurements indicate that near the center of the electron trap the random energy in the azimuthal direction is 0.32 eV while in the radial direction it is 0.13 eV. Outside of the trap, near the magnetron axis, the ions were measured to be room temperature. Including collisions in the simulation (case 2) yields good agreement between dE^*r and the experimental data. This is a major improvement over case 1 where the random energy is twice experimental values. The addition of turbulence (case 3) again has little effect.

A more accurate model of the differential cross section in the simulation might lead to better agreement with the radial drift velocities and the random energy in the azimuthal direction. The significant slowing of [v^*r in case 2 (Fig. 8) suggests that collisions are very important. A less forward-peaked scattering cross section might further slow down the radial drift and increase the azimuthal random energy in the region of the electron trap.

We can conclude that turbulence does not significantly influence ion transport. We found only small differences between the results of case 2 and case 3, despite the fact that the role of low-frequency turbulence was exaggerated by modeling it with the strongest field measured. [12]

Contours of the ion density ni(r,z), shown in Fig. 10, are markedly different then the ionization event location density from the electron simulation, shown in Fig. 11. Both the ionization events and the ion density are localized to elliptical regions in the r-z plane, with the ionization events concentrated near z = 0.25 cm and the ion density peaked near z = 0.7 cm. The difference is due to the high ion velocity in the cathode sheath, which reduces the density there for a given ion flux.

IV. SUMMARY

We have reported a Monte Carlo simulation of ion trajectories in our magnetron. We calculated the trajectories including ion-neutral collisions, a dc magnetic field, a dc electric field, and a turbulent electric field. To gauge the influence of collisions and turbulence, we repeated the ion simulation for three different cases: (1) no collisions and no turbulent electric field, (2) collisions but no turbulent field, and (3) both collisions and a turbulent electric field. We found that the dc electric field and collisions play major roles in ion transport, but the turbulent electric field does not.

We predicted the energy, angle, and spatial distribution of ion impact on the electrodes. Additionally, we tested the accuracy of the simulation by comparing to previously reported experimental measurements in the plasma. We found that the simulation successfully predicts: the etch track shape, ion transit time, radial profile of ion density, azimuthal drift velocity, and random ion energy for the radial velocity component. The simulation provides inaccurate predictions of the radial drift velocity and the azimuthal ion energy.

The ion Monte Carlo simulation reported here is a logical extension of our earlier electron transport model. These models enhance our understanding of magnetrons and may be useful tools for designing improved devices.

ACKNOWLEDGEMENTS

This work was funded by a grant from the Iowa Department of Economic Development and an IBM Faculty Development Award.

TABLE I

LIST OF VARIABLES

^{_________________________________________________________________________________________________________________________________________________}

Symbol Meaning

^{_________________________________________________________________________}

B magnetic field

E electric field

f^*(v) ion velocity distribution

Mi ion mass

ni(r,z) ion density

ni^* primary ion density

r radial coordinate

tk confinement time of a single ion

v velocity

[v^*r average radial drift velocity

[v^*[[theta]] average azimuthal drift velocity

z axial coordinate

dE^*r random energy in the radial direction

dE^*[[theta]] random energy in the azimuthal direction

[[Delta]]t integration time step

[[theta]] azimuthal coordinate

[[sigma]]el ion-neutral elastic collision cross section

[[sigma]]cx charge exchange collision cross section

[[tau]]transit ion transit time

* ^{denotes primary ions}

________________________________________________________________________

REFERENCES

[1] John A. Thornton and Alan S. Penford, in Thin Film Processes, ed. J. L. Vossen and W. Kern, (Academic, New York, 1978), p. 75.

[2] Brian Chapman, Glow Discharge Processes, (Wiley-Interscience, New York, 1980), p. 265 and pp. 177-184.

[3] T. E. Sheridan, M. J. Goeckner, and J. Goree, "Model of energetic electron transport in magnetron discharges," J. Vac. Sci. Technol., vol. A8, pp. 30-37 Jan./Feb. 1990.

[4] John A. Thornton, "Magnetron sputtering: basic physics and application to cylindrical magnetrons," J. Vac. Sci. Technol., vol. 15, pp. 171-177, Jan./Feb. 1978.

[5] Mark A. Lewis, David Glocker, and Jacob Jorne, "Measurements of secondary electron emission in reactive aluminum and titanium nitride," J. Vac. Sci. Technol., vol. A7, pp. 1019-1024, May/June 1989.

[6] A. E. Wendt, Ph.D. thesis, University of California, Berkeley, 1988 and A. E. Wendt, M. A. Lieberman, and H. Meuth, "Radial current distribution at a planar magnetron cathode," J. Vac. Sci. Technol., vol. A6, pp. 1827-1831, May/June 1988.

[7] J. E. Miranda, M. J. Goeckner, J. Goree, and T. E. Sheridan, "Monte Carlo simulation of ionization in a magnetron plasma," J. Vac. Sci. Technol., vol. A8, pp. 1627- 1631, May/June 1990.

[8] A. E. Wendt, M. A. Lieberman, and H. Meuth, "Radial current distribution at a planar magnetron cathode," J. Vac. Sci. Technol., vol. A6, pp. 1827-1831, May/June 1988.

[9] Lan Gu and M. A. Lieberman, "Axial distribution of optical emission in a planar magnetron discharge," J. Vac. Sci. Technol., vol. A6, pp. 2960-2964, Sept./Oct. 1988.

[10] F. Guimaraes, J. Almeida, and J. Bretagne, "Modeling of the energy deposition mechanisms in a argon planar discharge," submitted to J. Vac. Sci. Technol., June 6, 1990.

[11] M. J. Goeckner, J. Goree, and T. E. Sheridan, "Laser-induced fluorescence characterization of ions in a magnetron plasma," J. Vac. Sci. Technol., vol. A8, (in press 1990).

[12] T. E. Sheridan and J. Goree, "Low-frequency turbulent transport in magnetron plasmas," J. Vac. Sci. Technol., vol. A7, pp. 1014-1018, May/June 1989.

[13] Earl W. McDaniel, Collision Phenomena in ionized Gases, (Wiley, New York, 1964), p. 164.

[14] S. M. Rossnagel, "Gas density reduction effects in magnetrons," J. Vac. Sci. Technol., vol. A6, pp. 19-24, Jan. 1988.

[15] David L. Book, Physics Vade Mecum 2nd Ed., ed. Herbert L. Anderson, (American Institute of Physics, New York, 1989), pp. 280-281.

[16] C. E. Roberts, Ordinary Differential Equations: A Computational Approach, (Prentice-Hall, Englewood Cliffs, NJ, 1979).

[17] T. E. Sheridan and J. Goree, "Analytic expression for the electric potential in the plasma sheath," IEEE Trans. Plasma Sci., vol. 17, 884-888 , Nov./Dec. 1989.

[18] W. H. Cramer, "Elastic and inelastic scattering of low-velocity ions: Ne⁺ in A, A⁺ in Ne and, A⁺ in A," J. Chem Phys., vol. 30, pp. 641- 642, March 1959.

[19] I. Popescu Iovitsu and N. Ionescu-Pallas, "Resonant charge-exchange and the kinetics of ions," Sov. Phys.-Tech. Phys., vol. 4, pp. 781-791, Jan. 1960.

[20] A. M. Myers, J. R. Doyle, J. R. Abelson, and D. N. Ruzic, "Monte Carlo simulation of magnetron sputtering particle transport," submitted to J. Vac. Sci. Technol. (1990).

Biographies

Matthew J. Goeckner -- In 1982, he received B.S. degrees in mathematics and physics from Southern Illinois University. In addition to the magnetron simulations described here, he has reported laser-induced fluorescence and Langmuir probe measurements of magnetron and multidipole filament discharges. For this work he will receive a Ph.D. in plasma physics from The University of Iowa.

John A. Goree -- received the B.S. in applied physics from the California Institute of Technology in 1980, and the Ph.D. in plasma physics from Princeton University in 1985 for work on rf waves and far-infrared laser scattering diagnostics of plasmas. He has been an assistant professor of physics at the University of Iowa since 1985, where he has been working on transport processes in glow discharge and magnetron plasmas, laser-induced fluorescence diagnostics, and plasma wave theory. John received an IBM Faculty Development Award in 1986.

Terrence E. Sheridan, Jr. -- received the B.A. in physics from Hiram College in 1983 and the Ph.D. in plasma physics from Dartmouth College in 1987 for measurements of the electron distribution function in a magnetized laboratory plasma. He is currently doing post-doctoral research at the University of Iowa. Research topics include experimental measurement and Monte Carlo simulation of transport in magnetron plasmas and the theory of plasma sheaths.

FIGURE CAPTIONS

Fig. 1. Sketch of our cylindrically symmetric magnetron. [12] The electrons are mainly confined in the trap region near the cathode. While most of the ions strike the copper cathode, a small fraction of them hit the substrate anode where they may be buried in the deposited thin film. The coordinate system (r,[[theta]],z) used throughout this paper is indicated. In the simulations described in this paper, we assumed a 1 Pa argon discharge operated with a cathode bias of -400 V.

Fig. 2. Experimental plasma potential in the presheath. Ions fall down this potential surface from their ionization sites, mainly towards the cathode. These Langmuir probe measurements indicate electric potential drops of several volts in both the r and z directions. Because of symmetry, there is no potential drop in the [[theta]] direction. From these measurements, we calculated the dc electric field in the plasma presheath used in our electron and ion simulations. The cathode sheath, between 0 < z < 0.5 cm, has a 400 V potential drop not shown here. Reprinted from Ref. [11]

Fig. 3. Angular distribution of ion impact at the cathode. Most of the ions strike almost perpendicular to the cathode surface. The mean impact angle, defined with respect to the surface normal, is 7.95 +/- 0.14 degrees. We repeated the ion simulation for three different cases: (1) no collisions and no turbulent electric field, (2) collisions but no turbulent field, and (3) both collisions and a turbulent electric field. Here, and in Figs. 4-6, we show results only from case 3.

Fig. 4. Energy distribution of ion impact at the cathode. A large fraction of the ions are born in the sheath or lose energy in a collision; these strike the cathode with less than the maximum energy. The mean impact energy is 293.02 +/- 1.34 eV for a cathode bias of -400 V.

Fig. 5. Radial profile of the ion impact density at the cathode. In (a), the ion impact density from the ion simulation is compared to the experimental [3] etch track depth. The ion simulation takes into account the radial deflection of the ion trajectories. In (b), the ion impact density profile from the electron simulation is compared to the etch track. This profile assumes that ions fall directly to the cathode without radial deflection. The electron model shows better agreement with experiment than the more complicated ion model.

Fig. 6. Histogram of ion lifetime. This semi-log plot shows that ions born in the sheath (z < 0.5 cm) leave the plasma faster than those born in the presheath. From this histogram we calculate that the transit time is [[tau]]transit = 0.259 +/- 0.004 us. This time can be compared to the 0.8 us transit time measured by Sheridan and Goree. [12] The minimum at 0.9 us corresponds to ions born 6 mm above the cathode.

Fig. 7. Radial profile of ion density in the plasma presheath. Here we compare the ion density profiles from the ion simulation to earlier laser-induced fluorescence (LIF) measurements. [11] To gauge the influence of collisions and turbulence, we show results from all three cases simulated. Case 2 and 3 include collisions with neutrals, and show the best agreement. Including turbulence in the simulation (case 3) has little effect on this or any other result.

Fig. 8. Radial profile of ion drift velocity in the plasma presheath. We compare the drift velocities from the LIF experiment [11] and the ion simulation. In the azimuthal direction (a), all three cases agree with experiment, but in the radial direction (b), none agree. This disagreement is made less severe by including collisions with neutrals (case 2 and 3) in the simulation.

Fig. 9. Radial profile of random ion energy in the plasma presheath. We compare the random ion energies from the LIF experiment [11] and the simulation. These random energies are computed from the velocity components in the (a) azimuthal and (b) radial directions. As with the drift velocity in Fig. 8, we find agreement with only one component. Here, however, it is the radial direction (b) where the random ion energy shows the best agreement, especially when collisions are included in the ion simulation (case 2 and 3).

Fig. 10. Ion density contours on the r-z plane in the discharge. This data from case 3 of the ion simulation shows that the ion density is peaked near r = 1.5 cm, z = 0.7 cm.

Fig. 11. Ionization event density contours. This data from the electron simulation shows that the ionization events are concentrated near r = 1.7 cm, z = 0.25 cm, which is closer to the cathode than the ion density peak in Fig. 10. This difference is attributable to the transit time of the ions, which is longer in the presheath (z > 0.5 cm) than in the cathode sheath.