Behavior of electrons in a dual-magnetron sputter deposition system: a Monte Carlo model

A Monte Carlo model has been developed for investigating the electron behavior in a dual-magnetron sputter deposition system. To describe the three-dimensional (3D) geometry, different reference frames, i.e. a local and a global coordinate system, were used. In this study, the influence of both closed and mirror magnetic field configurations on the plasma properties is investigated. In the case of a closed magnetic field configuration, the calculated electron trajectories show that if an electron is emitted in (or near) the center of the cathode, where the influence of the magnetic field is low, it is able to travel from one magnetron to the other. On the other hand, when an electron is created at the race track area, it is more or less trapped in the strong magnetic field and cannot easily escape to the second magnetron region. In the case of a mirror magnetic field configuration, irrespective of where the electron is emitted from the cathode, it cannot travel from one magnetron to the other because the magnetic field lines guide the electron to the substrate. Moreover, the electron density and electron impact ionization rate have been calculated and studied in detail for both configurations.


Introduction
Magnetrons have found wide applications in the sputter deposition of different types of industrially important coatings, such as metals, oxides and nitrides [1]. In a magnetron discharge, in addition to an electric field, a magnetic field is also present, causing the electrons to be trapped close to the cathode. This leads to a more efficient use of the electrons for ionizing the gas atoms, resulting in enhanced sputtering of the target in comparison with a non-magnetized discharge. Moreover, the pressure can be lower. Hence, by adjusting the pressure, one can control the level of thermalization of the sputtered atoms to tune the thin film properties.
In the last three decades, many types of magnetrons have been developed and the characteristics of magnetron discharges have been studied intensively (see e.g. [2]- [4]). Each of these magnetrons has its own specific advantages and applications. Two single magnetrons can also be combined into a dual configuration. Different reasons can be found for using this configuration. One application of the dual-magnetron configuration is to enhance the process stability in industrial coaters, where the power is switched between two cathodes consisting of the same material [5]. Another application of dual magnetrons is the so-called facing target sputtering. In this setup, two cathodes, mostly consisting of the same material, are facing each other. In this way, a well-defined plasma is generated between the two sources. The substrate is placed outside this region to eliminate substrate bombardment with high energetic species and substrate heating [6]. For this reason, this configuration is widely used for depositing sensitive materials, such as hightemperature superconductors [7,8]. Recently, this configuration was studied for sideways deposition of thin films using dual-HiPIMS [9]. A final, but perhaps most important, reason for using a dual-magnetron configuration is related to materials research of complex materials. Indeed, to modify the stoichiometry in a flexible way, a dual-magnetron configuration is often chosen. By changing the target-substrate distance and/or the target power for both cathodes independently, a large range of compositions can be investigated [10]- [13]. This flexibility comes at a cost, especially in reactive magnetron sputtering. Indeed, both cathodes will react independently on addition of the reactive gas, resulting in possible target poisoning for one target. Nevertheless, using this configuration can speed up materials research by the combinatorial deposition of thin films. This is not limited to two sources but sometimes 3 extended to a multisource approach [14][15][16]. Some companies are offering this kind of setup 4 .
An important issue in a dual-magnetron configuration is the magnetic field design. Two designs are possible, i.e. closed field and mirror field. In the mirror-field configuration, the magnets of both magnetrons are mounted in the same way (i.e. two magnetrons with the same magnetic polarity). On the other hand, in the closed-field configuration, the magnets of both magnetrons are mounted in the opposite way, so that the field lines close in across the chamber. Electrons following these field lines can give rise to ionizing collisions and, thus, maintain a high plasma density in the vicinity of the substrate. In the mirror field, the electrons are directed towards the chamber walls, resulting in a lower plasma density. This is nicely demonstrated by Musil and Baroch [17,18]. As such, the choice of the magnetic field design can have a strong influence on the thin film growth.
A certain amount of experimental work on the plasma characterization of these configurations has already been performed, but there is still a lack of a detailed description. As dual-magnetron configurations are becoming more popular in materials research of complex materials, there is a need to deeply investigate this configuration. Modeling can assist in understanding this configuration. However, to the authors' knowledge, no modeling has been performed so far on dual-magnetron configurations.
There are different modeling approaches in the literature for properly describing the processes in a magnetron discharge. In an analytical model [19], mostly semi-empirical formulae describe the relation between certain plasma quantities and macroscopic parameters (e.g. voltage, current and pressure). This method can rapidly predict the plasma behavior, but it is only an approximation, valid for a limited range of conditions. A fluid model [20] is based on the continuity equations of particle density, momentum and energy, usually coupled to Poisson's equation to calculate the electric field distribution selfconsistently. It is also rather simple and fast, but it assumes that the plasma species are more or less in equilibrium with the electric field, i.e. the energy gain from the electric field is more or less balanced by the energy loss due to collisions. This is not always true, e.g. at low gas pressure (where the mean free path of charged particles exceeds the characteristic length of the discharge) and for the fast electrons in regions characterized by a strong electric field.
A Boltzmann model takes into account the non-equilibrium behavior of the plasma species by solving the Boltzmann transport equation. However, this method can become mathematically complicated, certainly when solving in more than one dimension. In the works of Porokhova et al [21,22], this model has been used for cylindrical postmagnetrons, where the magnetic field is one dimensional (1D) and constant in the largest part of the discharge region, which simplifies the Boltzmann equation.
Monte Carlo (MC) simulations [23], on the other hand, are mathematically simple and account for the non-equilibrium behavior of the plasma species correctly. Indeed, they describe the trajectories of individual particles by Newton's laws, and treat the collisions by random numbers. However, in order to reach statistically valid results, a large number of individual particles need to be simulated. Hence, MC simulations can be quite time consuming, especially for slow-moving particles. Moreover, an MC model on its own is not self-consistent because it requires a certain electric field distribution as input.

4
This problem can be overcome by the particle-in-cell/Monte Carlo collisions (PIC/MCC) model that couples MC simulations for the behavior of ions and electrons, to the Poisson's equation for a self-consistent electric field calculation. The PIC/MCC model does not need many a priori physical assumptions except for the classical nature of particle trajectories, being very useful for non-equilibrium processes. Hence, it is the most accurate technique applied to 2D [24,25] or 3D [26] simulations. A drawback of the PIC/MCC model, however, is that it requires a much longer computational time, certainly for larger reactors, such as a dualmagnetron configuration.
Finally, hybrid models combine different modeling approaches, achieving more reasonable computation times. For instance, in [27,28], fast electrons (such as electrons, emitted at the cathode), which are not in equilibrium with the electric field, are treated by an MC model, whereas for slow electrons and ions, which can be considered in thermal equilibrium, the fluid equations are used.
Based on the above-mentioned considerations, in this work, we develop an MC model to investigate the behavior of the fast electrons in an Ar/O 2 direct current dual-magnetron discharge. Indeed, the calculated electron properties will already give an idea of the plasma behavior near the target. Moreover, in the future, the electron MC model will be combined with other models to be developed (e.g. MC models for other species and surface models for the plasma-surface interactions) into a hybrid modeling approach, in order to describe all of the plasma species in a dual-magnetron configuration and to investigate the sputter deposition processes of complex oxide layers. The results from the electron MC model (e.g. electron impact ionization rate) can then be used as input values in the other models. At this stage, we assume that both targets are made of titanium. Later on, when the models for sputtering are included, this will be changed to different targets in order to describe the deposition of complex oxide films.

Description of the model
The behavior of fast electrons in the dual-magnetron discharge is studied by the MC method. The term 'fast electrons' is used here for those electrons with kinetic energies above the threshold for inelastic collisions (3 eV), which are present in the discharge at positions where the absolute value of the electric field is above 1000 V m −1 ; the latter is done to prevent electrons with a kinetic energy below 3 eV from being removed from the simulation, as they can be accelerated again. A detailed overview of the MC collision algorithm can be found in [29]- [31]. In this paper, only a brief description is given, but special attention will be paid to the use of different reference frames. The presented model is 3D in both the coordinate and the velocity space. The real particles in the discharge are represented by a limited ensemble of superparticles (SPs), with a certain weight corresponding to the number of real particles per SP. In our simulations, the initial number of electron SPs, emitted from the cathodes is taken as 40 000.
During successive time steps, the SPs' trajectory under the influence of both the electric and magnetic fields is calculated by Newton's equations of motion, based on the Lorentz force, where q is the elementary charge, v is the velocity, m is the electron mass and E and B are the electric and magnetic fields, respectively. The equation is discretized using the central finitedifference method (the so-called leap-frog method) [32]. The v×B rotation term is treated according to the algorithm suggested by Boris [33].
In the middle of each time step, the collisions are treated by random numbers; that is, the total collision probability, summed over all possible collision processes, is calculated based on the collision cross sections and the target species densities, and compared with a random number uniformly distributed between 0 and 1. If the probability is lower than this random number, no collision occurs. If it is higher, a collision takes place, and the new energy and direction after the collision are determined.
The following collisions are considered in this work: collisions with argon atoms (elastic collision, excitation to metastable and radiative state and ionization), ionization of titanium atoms, collisions with argon metastable atoms (transfer to a radiative state and ionization), collisions with oxygen molecules (ionization, excitation to a radiative state, dissociative attachment, ion pair formation, dissociative ionization and dissociation) and collisions with oxygen ions (recombination, neutralization and dissociative recombination). All of the reactions and their corresponding cross sections have been adopted from [31]. The effect of Bohm (or anomalous) diffusion is also included in the model [34,35].
The dual-magnetron setup is schematically represented in figure 1. Both magnetrons are positioned at angles of 45 • with respect to the substrate. For convenience, we use two different frames of reference: the collisions are treated in the local coordinate systems (LCS), the position (i.e. the trajectory) and the velocity of the electrons are calculated in the LCS as well, but they are also transferred to the global coordinate system (GCS), so that the electrons are able to move from one magnetron to the other (in other words, from one LCS to the other). The transformation 6 from the LCS to the GCS is performed in the following way: Left magnetron: Coordinates Right magnetron: Coordinates where x global , y global , z global (v x,global , v y , global , v z,global ) and x local , y local , z local (v x,local , v y,local , v z,local ) are the coordinates (velocity components) of the electron in the GCS and LCS, respectively; x center1 , y center1 , z center1 and x center2 , y center2 , z center2 are the global coordinates of the origins of the LCS (or in other words the global coordinates of the target centers) of the left and right magnetrons, respectively; and finally, ω 1 and ω 2 are the angles (ω 1 = ω 2 = 45 • ) between the LCS and the GCS of the left and right magnetrons, respectively (see figure 1). The electric and magnetic fields are needed as input in the model. The electric field has been adopted from PIC-MCC simulations [31] in one magnetron; hence, it is defined in the LCS. Accurate values of the magnetic field can be obtained from experiments, but can also be calculated analytically (see e.g. [36]- [38]) or by using finite-element models (see e.g. [39]- [41]). The analytical models are fast and give a smooth variation in the magnetic field, but they can only be applied in relatively simple magnetic configurations. On the other hand, finite-element models can deal with very complex magnetic configurations. There exist several packages in the literature (e.g. [40], [42]-[44]). In the present work, we use the open-source finite-element solver GetDP [45] to calculate the magnetic field distribution in a full 3D configuration, together with the finite-element mesh generator Gmsh [46]. Hence, the magnetic field is calculated in the GCS and transferred to the LCS.
In this paper, we study both a closed (i.e. the magnetic field lines go from one magnetron to the other) and a mirror magnetic field configuration (i.e. the magnetic field lines are reflected 7 from each other) (see [47] and see below). The polarity of the magnetron configuration becomes important when several magnetrons are combined in one vacuum chamber. In a closed-field configuration, where the magnets within the magnetrons are arranged such that alternating poles are next to each other, resulting in the linkage of the field lines, the electrons are prevented from escaping to the chamber walls, resulting in much higher ion current density. This can result in much denser, harder and well-adhered thin films [48]. The mirror-field configuration is not so widely used, but in larger systems it is sometimes combined with two sources in closedfield configuration. In this way, electrons can be guided towards the center of the deposition chamber.

Operating conditions
The direct current dual-magnetron configuration under study is a sputter deposition device with Ti targets with a radius of 25 mm. It operates in an Ar/O 2 gas mixture at 300 K with partial pressures of 1 Pa Ar and 0.24 Pa O 2 . Under these conditions, the targets are completely poisoned. The electric field is calculated self-consistently with a PIC/MCC model [31] and taken as input in the MC simulation. This electric field gives rise to a cathode current of 0.2 A and a cathode voltage of −312 V. The distributions of the components of the electric fields (E x and E z ) are presented in figure 2. In the sheath in front of the cathode, E z is approximately ten times higher than the corresponding E x (and E y ) values, which is logical because of the strong potential drop from the cathode to the bulk plasma. Since the cathode is powered by a negative voltage, the cathode sheath is much more pronounced than the sheath in front of the substrate (= anode). Indeed, the maximum value of E z (figure 2(a)) in front of the cathode is −400 kV m −1 , whereas in the substrate sheath it is at maximum 10 kV m −1 , hence two orders of magnitude lower. The same is true for the distribution of E x ( figure 2(b)). Moreover, the values in the bulk area are still much lower both for E r (where E 2 r = E 2 x + E 2 y ) and E z (∼ 100 V m −1 ). Note that the values of the electric field were only calculated in the PIC-MCC model till 24 mm in the z-direction and 28 mm in the r-direction of the LCS. Since the reactor under study here is 70 mm in the z-direction and 55 mm in the r-direction (or in the xand y-directions, r 2 = x 2 + y 2 ) in the LCS, the electric field values in the bulk plasma are extrapolated based on the PIC-MCC values. This approach is justified as the electric field values in the bulk do not influence the processes in the magnetron a lot, since most of the plasma processes occur in the regions near the cathodes. However, in the sheath in front of the substrate, the exact values as calculated from the PIC-MCC simulations were used as input values.
As mentioned above, the magnetic field distributions are calculated with the open-source finite-element solver GetDP [45] and presented in figure 3. Figure 3(a) presents a closed and figure 3(b) a mirror magnetic field configuration. It is clearly observed in figure 3(a) that the magnetic field lines go from the left magnetron to the right one (characteristic for a closed magnetic field configuration), whereas in figure 3(b) the magnetic field lines are reflected from each other in the center of the dual-magnetron configuration (characteristic for a mirror magnetic field configuration). Furthermore, in both magnetic field configurations, the maximum radial magnetic field strength in front of the cathode is around 0.1 T, as is found at a radial position of 14 mm from the center of the target. In the bulk area of the

Calculated trajectories for one electron
3.2.1. Closed magnetic field configuration. As mentioned before, the electrons are emitted from the cathode and will be accelerated away from the cathode by the electric field. However, they are also affected by the magnetic field: the electrons are trapped along the lines of this magnetic field and can hence travel in the transverse direction. Moreover, the electrons have a low mass, which causes them to circulate around the magnetic field lines with a small Larmor radius. As a consequence, the electrons stay longer in the discharge and therefore they have a higher probability for ionization of the background gas than in a non-magnetized plasma.
The calculated trajectory of an electron, which is emitted in the middle of the strong magnetic trap, i.e. at the race track area of the target, is shown in figure 4. Nearly all of the electrons are created here. It is clear that the electron circulates around the center of the cathode with a spiral trajectory and can also move further away from the cathode target due to elastic collisions with the background gas atoms, as well as Bohm diffusion. Elastic collisions 'knock' the electron from the magnetic trap to another magnetic field line. Bohm diffusion also results in enhanced electron mobility, since the electric field fluctuations push the electron away from its magnetic trap. Figure 4 also illustrates that in a closed magnetic field configuration, the electron can escape from one magnetron to the other. However, it is worth noting that only 12% of the electrons do cross from one magnetron to the other, as predicted by our simulations (see figure 5; 'crossing electrons'). This is caused by the fact that the energy of the electrons, once they reach the bulk plasma, is typically not high enough to reach the second magnetron. From our calculations, it was observed that electrons can cross the center (in between both magnetrons) up to 18 times (see figure 5, inset).
In order to show some statistics about the electron trajectories, a histogram of the fraction of electrons with corresponding end positions in the z-direction of the LCS is plotted in figure 5. It is clear that most of the electrons (almost 70%) will be absorbed again at the target, after one or several circulations around the magnetic field lines. It should be mentioned that electron impact ionization also takes place in the cathode sheath, resulting in a multiplication of the electrons, and the multiplication factor can even reach two orders of magnitude. The peak around 10 mm  corresponds to the electrons, which reach higher positions above the target but cannot escape from one magnetron region to the other. As mentioned above, about 12% of the electrons can escape from one magnetron to the other. In figure 6, the trajectory is illustrated for an electron that is released at the center of the cathode, where the radial magnetic field is very weak (cf figure 3(a)) and strong trapping does not occur. Moreover, at the cathode, the electric field is strong, so the electron has a high initial potential energy. Figure 6 clearly shows that the electron can now easily escape from one target to the other and can also return to its original side. Our calculations predict that this occurs for all electrons released at the cathode center. However, it is important to note that only 0.1% of all electrons are emitted near the target center. Although this percentage is low, it is interesting to investigate this process to gain a better understanding of the behavior of all electrons in the dual-magnetron system.

Mirror magnetic field configuration.
In the case of a mirror magnetic field configuration, the calculated trajectory of an electron, which is emitted in the middle of the strong magnetic trap, is virtually the same as in the case of a closed magnetic field configuration, but it should be noted that in this case the electron cannot escape from one magnetron region to the other due to the mirror magnetic field configuration.
A histogram of the fraction of the electrons with corresponding end positions in the z-direction of the LCS, for the mirror magnetic field configuration, is shown in figure 7. Similar to the closed magnetic field configuration, most of the electrons (more than 70%) will be absorbed again at the target, after one or several circulations around the magnetic field lines. The peak around 10 mm again corresponds to the electrons reaching higher positions above the target. However, in contrast to the closed magnetic field configuration, there is another peak between 50 and 60 mm, which represents the electrons being drawn towards the substrate in the mirror-field configuration.
If an electron is released at the center of the cathode in the mirror-field configuration (see figure 8), it can also not travel from one magnetron to the other because the mirror magnetic field configuration guides the electron towards the substrate. In summary, we can conclude that in the case of the mirror magnetic field configuration, the electrons, independent of where they are emitted from the cathode, cannot escape from one target to the other but are drawn towards the substrate.

Fast electron density and electron impact ionization rate
3.3.1. Closed magnetic field configuration. Important results of this model include the electron impact ionization rate and the fast electron density, because they give information about the regions of strongest plasma activity. Furthermore, in a hybrid modeling approach, these output data can be used as input for the treatment of other plasma species.
The calculated 2D electron density profile, for the case of the closed magnetic field configuration, is shown in figure 9(a). Note that this profile reflects only the so-called fast electrons, which are treated with this MC model. The thermal electrons are not taken into account in this model in order to save computation time. Most of the (fast) electrons are located near the target due to the strong magnetic trap, and consequently the electron density reaches its maximum there. The electron density agrees well with PIC/MCC results in comparable conditions [31], showing that the effect of Bohm collisions is taken into account in our MC model in a realistic way.
The maxima of the profiles are somewhat tilted towards the centers of both cathodes, which is caused by the similar tilted shape of the magnetic field (see figure 3). Moreover, the profile spreads out more in the z-direction than in the r-direction because the z-component of the electric field is much larger than its r-component, causing most of the fast electrons to move in the z-direction and remain fast.
As mentioned in section 3.2, electrons might be able to move from one target to the other because of the closed magnetic field configuration. It is also clear from figure 9(a) that the electron density profiles of both magnetrons overlap at x = 0. Figure 9(b) presents the electron impact ionization rate of Ar gas atoms. Most of the ionizations occur where the fast electron density is at maximum (see figure 9(a)). Similar to the fast electron density profile, the ionization rate profile is somewhat tilted towards the centers of both cathodes, as well as more spread out in the z-direction than in the r-direction and the profiles of both magnetrons overlap at x = 0. Between the targets, less ionization takes place compared to the target sheath, because of the lower electron density.

Mirror magnetic field configuration.
In the case of a mirror magnetic field configuration, similar profiles are found for the fast electron density (figure 10(a)) and electron impact ionization rate of Ar gas atoms (figure 10(b)) in the regions close to the cathode, where there is a strong magnetic trap. However, in both figures, the profiles are aligned towards the substrate rather than towards the center of the dual-magnetron reactor, in contrast to the case of the closed magnetic field configuration (see figures 9(a) and (b)). As mentioned above (section 3.2.2), due to the mirror magnetic field configuration, the electrons cannot escape from one magnetron to the other; they mostly move towards the substrate. Therefore, the fast electron density in between both targets is very low and consequently no ionization takes place here. This result is important because it will also affect the behavior of the other plasma species to some extent, and therefore it can explain the different plasma behavior and application field of the closed-and mirror-field configurations.

Conclusion
To simulate the behavior of the fast electrons in a dual-magnetron sputter deposition device, an MC model was developed. A total of 18 types of electron collisions are included in this model, such as collisions with background gas atoms (elastic scattering, excitation to metastable and radiative states and ionization), with target atoms (ionization) and with reactive gas molecules and ions.
The electron trajectory was thoroughly studied. It was found that most of the electrons are created and absorbed again at the target, after one or several circulations around the magnetic field lines in both magnetic field configurations. In the case of the closed magnetic field configuration, 12% of all created electrons can escape from one magnetron region to the other, and some of them can even cross the center between the two magnetrons a number of times. The other 88% of all electrons cannot escape from one magnetron region to the other because their energy has become too low. In the mirror magnetic field configuration, no electrons can escape from one magnetron region to the other because the magnetic field lines guide them towards the substrate.
The calculated electron density and the electron impact ionization rate were also presented. Due to the strong magnetic trap, most of the electrons reside near the cathode and most of the collisions occur in this region. The closed magnetic field configuration and the resulting electron trajectories cause the fast electron density profile to overlap in between the two-magnetron regions. In the mirror configuration, this profile is aligned towards the substrate and almost does not overlap in between the two magnetron regions. These results can give insights into the behavior of complex deposition setups where several magnetrons are combined.
In the future, this electron MC model will be coupled to other models, to treat all of the plasma species and to study the relevant plasma processes in the dual-magnetron sputter deposition system. Moreover, surface models will be included to describe the processes occurring at the target and the substrate, to investigate the sputter deposition process of complex oxide thin films.