Simulation study on electron heating characteristics in magnetic enhancement capacitively coupled plasmas with a longitudinal magnetic field

1.   Introduction

The parameters of capacitively coupled plasmas (CCPs), such as plasma density and electron energy distribution function (EEDF), etc, play important roles in material etching or thin film deposition, and they usually depend on how electrons gain energy from driving electric fields, which is traditionally known as electron heating [1].

Ohmic heating is an important mechanism in electron heating, which origins from the collision between electrons and neutral particles [2]. In this case, the synchronous relationship between the motion of the colliding electrons and the driving electric fields is destroyed, and therefore the electrons obtain net energy. In some cases, however, the strong power deposition in CCPs cannot be explained well by Ohmic heating, and there must be some additional 'collisionless' mechanisms to heat electrons, as pointed out by Godyak et al [3, 4]. A classical model described collisionless heating is known as hard-wall model, which proposed by Lieberman [5]. In this model, the electric filed in plasma zone is ignored, and the collisionless heating is considered as a result of the random interactions between the moving electrons and the oscillating sheath electric filed. Some further researches show that the strong electric filed reversals in the region between the sheath and plasma zone can also introduce additional collisionless heating [6, 7], which is known as ambipolar electron heating [8]. In addition, the electrostatic waves originated from the energetic electron beam propagating from the sheath edge toward plasma zone can also contribute to the collisionless heating in some cases [9]. Another classical model proposed to understand the collisionless heating is known as the pressure heating model. In this model, the collisionless heating is considered to be related to the pressure action in electronic fluid [10, 11].

To analyze the dependence of collisionless heating on discharge parameters, some analytical models have been proposed to describe the level of collisionless heating [6, 12]. These models, however, often show some errors when compare their results with that from particle-in-cell (PIC) simulations [13], and cannot describe all the physical mechanisms of collisionless heating completely. To determine the level of the heating accurately, Surendra and Dalvie proposed a quantitative method, in which a dynamic fluid model is used to analyze the heating, and parameters of this model are taken from simulations [14]. Through this method, heating components from different physical mechanisms can be distinguished, and thereby which components are important can be determined. They further revealed that the pressure tensor of electrons is the physical origin of collisionless heating, and time-varying electron temperature is very important for this heating. Based on such a method, people have a clearer understanding of electronic heating in CCPs. For example, Lafleur et al found that the total heating in their most CCPs tested cases consists mainly of Ohmic heating and pressure heating, with the additional collisionless heating associated with the electron inertia being negligible [15]. Liu et al found that the thick dielectric sidewall in CCPs can enhance the level of the axial Ohmic heating and result in the increase of axial electron density [16]. Schulze et al found that the strong bipolar electric fields are very important for pressure heating in electropositive CCPs [17], and Proto et al further expanded this conclusion to electronegative CCPs recently [18, 19].

If one introduces a transverse magnetic field which oriented parallel to the electrodes into CCPs, then the electron heating as well as plasma characteristics of the CCPs may be modulated [20]. For example, transverse magnetic field can enhance the level of Ohmic heating significantly, and induce the Ohmic heating to become dominant in total heating [21, 22]. In addition, the electrons near the magnetized plasma sheath edge can undergo periodic cyclotron motion, and in some cases the couple between the oscillating sheath and the cyclotron electron can enhance the electron heating, which is known as electron bounce-cyclotron resonance heating [23–25]. On the other hand, some quantitative studies about the influences of transverse magnetic field on electron heating in CCPs have also been reported. For example, Wang et al found that with the increase of transverse magnetic field, the dominant heating mechanism in magnetized oxygen CCP will translate from electronegative drift-ambipolar heating to electropositive stochastic electron heating [26]. Zheng et al revealed that the transverse magnetic field induced Hall currents play important roles in enhancing the level of Ohmic heating in CCPs [27].

In summary, a series of accurate quantitative researches have discussed the effects of transverse magnetic fields on the electron heating in CCPs, few studies, however, report the effects of longitudinal magnetic fields on that in CCPs. In this work, the electron heating characteristics in a magnetic enhancement CCP with a longitudinal uniform magnetic field are explored through both theoretical and numerical calculations. This paper is structured as follows. In section 2, a dynamic fluid model used for analyzing the electron heating characteristics is reviewed. The simulation model used to perform numerical calculations and its parameter setting is described in section 3. The effects of longitudinal magnetic fields on electron heating characteristics are stated in section 4. The last part is a summary of the full text.

2.   Theoretical background

For a magnetized plasma, the first momentum Boltzmann equation for electrons is the momentum balance equation, read as [28]

where $n_{\mathrm{e}},$ ${\mathit{u}}_{\mathrm{e}}$ and ${\mathbf{\Gamma }}_{\mathrm{e}}$ represent the density, mean velocity and flux of electrons, respectively. E and B represent driving electric fields and applied magnetic fields, respectively. $\stackrel{⃡}{{\mathit{P}}_{\mathrm{e}}}$ represents the pressure tensor of electrons and ${\mathrm{\Pi }}_{\mathrm{c}}$ represents the momentum change resulting from electron–neutron collisions. Since elastic collision is dominant in the electron–neutral collisions, the ${\mathrm{\Pi }}_{\mathrm{c}}$ can be obtained by subtracting the electron momentum before elastic collision from that after elastic collision. In this work, we take electron positions at integral time of $t$ and electron velocities at half integral time of $t+\mathrm{\Delta }t/2$ to calculate ${\mathrm{\Pi }}_{\mathrm{c}},$ a more detailed discussion of ${\mathrm{\Pi }}_{\mathrm{c}}$ can be found in Lafleur et al [15]. Multiplying equation (1) with the mean velocity ${\mathit{u}}_{\mathrm{e}},$ the new equation can be expanded as follows,

where the subscript $i=x,y,z$ and ${\mathit{E}}_{\mathrm{H}}={\mathit{u}}_{\mathrm{e}}\times \mathit{B}$ represents Hall fields. $p_{\mathrm{e},xi}=m_{\mathrm{e}}{n}_{\mathrm{e}}\left(〈v_x{v}_i〉-u_x{u}_i\right)$ represents the elements of the pressure tensor, where $v$ represents velocities of the electrons. The term $-e{n}_{\mathrm{e}}{u}_{\mathrm{e},i}{E}_i=P_{\mathrm{abs}}$ represents the level of the electron heating. In the case of a one-dimensional (1D) plasma with a driving electric field of $\mathit{E}=\left(E_x,0,0\right),$ $P_{\mathrm{abs}}$ obtained from equation (2) can be expressed in two forms, i.e. the heating in x-direction including magnetic field term

or the heating in global space excluding magnetic field term, $P_{\mathrm{abs}}=P_{\mathrm{in}}+P_{\mathrm{press}}+P_{\mathrm{ohmic}},$ where

represent the level of inertia heating, pressure heating and Ohmic heating, respectively. In the following studies, $P_{\mathrm{in}}$ is ignored since the time averaged $P_{\mathrm{in}}$ is usually unimportant [15].

3.   Simulation descriptions

In this work, an implicit 1d3v PIC/MCC code is used to perform all the simulations, and the code has been benchmarked with the results from Zheng et al [24] and the results from Yang et al [29]. We consider a magnetic enhancement CCP discharge system operated in argon, as shown in figure 1. In this model, we assume the working pressure and gas temperature are 100 mTorr and 300 K, respectively. The transverse magnetic field $B_y$ is a fixed uniform magnetic field, while the longitudinal magnetic field $B_x$ is a tunable uniform magnetic field. The driving electric field $E_x$ is generated by a voltage source having a waveform of $\tilde{V}=V_0\cos (2\pi ft),$ where $V_0=300$ V and $f=13.56$ MHz. The gap between the electrodes is 2.5 cm, and we have 80 equal cells in this gap.

Figure  1.  The model of the magnetic enhancement CCP operated in argon.

DownLoad: Full-Size Img PowerPoint

In all simulations, we consider electron–neutral collisions and argon ion–neutral collisions, where the collision types and the collision cross sections are shown in figure 2. The cross sections of electron–neutral collisions are taken from [30], in which for elastic electron–neutral collision the momentum transfer cross-section is used, and all scattering events are treated as isotropic, as referred as [17]. The cross-sections of argon ion–neutral collisions are taken from [31], and elastic ion–neutral collision is isotropic, while the ion after the charge exchange collision (backscattering collision) has a velocity equal to the velocity of the neutral atom that collides with the ion. The velocities of particles after collision events are calculated by the method in [32].

Figure  2.  The collision cross sections used in simulations, which are taken from [30, 31].

DownLoad: Full-Size Img PowerPoint

The initial electrons and ions have Maxwell distributions with temperatures of 2 eV and 0.026 eV, respectively, and have a uniform profile in the gap. When the simulations reach a steady state, there are about 60 000–80 000 super-electrons (the number of super-ions is similar). We use a double-time-scale algorithm in all simulations, that is, the time step of electrons and that of ions are $1\times {10}^{-10}$ s and $4\times {10}^{-10}$ s, respectively. The secondary electron emission process is not considered here to simplify the analysis and it will not have a significant impact on our main conclusions. All simulations run for 5000 radio frequency (RF) periods and take about 3000 periods to reach a steady state. The macroscopic quantities, such as $n_{\mathrm{e}},$ $p_{\mathrm{e}}$ and electron temperature $T_{\mathrm{e}},$ etc, are calculated by a linear weighting scheme [33], and their time averaged value in 1500 RF periods are taken as the input parameters for the fluid model described in section 2.

4.   Results and discussions

Evolutions of $n_{\mathrm{e}}$ at different uniform $B_x$ are shown in figure 3(a), where the fixed uniform $B_y$ is 50 G. It can be seen that the global $n_{\mathrm{e}}$ decreases gradually with the increased $B_x,$ which implies a weakened electron heating. Figure 3(b) further demonstrates evolutions of the peak $n_{\mathrm{e}}$ and the averaged $n_{\mathrm{e}}$ with $B_x.$ When $B_x$ increases from 0 to 40 G, the peak $n_{\mathrm{e}}$ decreases from $2.65\times {10}^{16}$ to $1.05\times {10}^{16}$ m⁻³, while the averaged $n_{\mathrm{e}}$ decreases from $1.4\times {10}^{16}$ to $0.47\times {10}^{16}$ m⁻³, both electron densities decrease. As $B_x$ continues to increase, however, the rates of decrease of both electron densities slow down gradually, which indicates that there is a limit to the effect of $B_x$ on the heating.

Figure  3.  Evolutions of $n_{\mathrm{e}}$ versus different $B_x,$ where $B_y$ is 50 G. (a) Global electron density, (b) peak electron density as well as averaged electron density.

DownLoad: Full-Size Img PowerPoint

Figure 4 shows evolutions of the total heating profile $〈E_x·J_x〉$ with different $B_x,$ where $J_x=-e{n}_{\mathrm{e}}{u}_{\mathrm{e},x}$ represents current density in x-direction and $e$ is unit charge. The figure also shows $P_{\mathrm{abs}}$ as well as the different heating components calculated from equation (4). It can be seen that the results from the theoretical calculations agree very well with those directly from numerical calculations in all cases. Meanwhile, $P_{\mathrm{press}}$ is negative in the plasma zone, and oscillates at the sheath boundary, while $P_{\mathrm{ohmic}}$ makes an important contribution to the heating in plasma zone, and it reaches a maximum near the sheath region. As $B_x$ increases continuously, the level of all these heating components decreases gradually, and the heating in plasma zone tends to be negligible.

Figure  4.  Evolutions of both the total heating profile and each heating component profile versus different longitudinal magnetic fields. (a) $B_x=0$ G, (b) $B_x=10$ G, (c) $B_x=20$ G.

DownLoad: Full-Size Img PowerPoint

The phenomenon of $B_x$ -induced reduction in the heating is related to the electron confinement in the z-direction. It can be seen from equation (3) that the total heating can be expressed as a sum of a series of heating components in x-direction, i.e. $P_{\mathrm{abs}}=P_{\mathrm{in},x}+P_{\mathrm{press},x}+P_{\mathrm{ohmic},x}+e{n}_{\mathrm{e}}{u}_{\mathrm{e},x}{E}_{\mathrm{H},x},$ where the term $e{n}_{\mathrm{e}}{u}_{\mathrm{e},x}{E}_{\mathrm{H},x}$ can be considered as a magnetically induced heating [26]. With the increased $B_x,$ the Hall field $E_{\mathrm{H},x}$ decreases since the increased $B_x$ constrains electrons more so that $u_{\mathrm{e},z}$ decreases [2], and this results in a reduction in the level of the heating. On the other hand, as $B_x$ continues to increase, the term $e{n}_{\mathrm{e}}{u}_{\mathrm{e},x}{E}_{\mathrm{H},x}=-e{n}_{\mathrm{e}}{u}_{\mathrm{e},x}{u}_{\mathrm{e},z}{B}_y$ is expected to be negligible since $u_{\mathrm{e},z}$ tends to zero, which can explain why there is a limit to the effect of $B_x$ on the heating.

Figure 5 shows the spatio-temporal evolutions of the total heating as well as that of different heating components at different $B_x$ in one RF periodic. It can be seen from figures 5(a1)–(c1) that the electrons near the sheath boundary lose and gain energy periodically with the oscillating sheath, and the loss of energy is correlated highly with the energy loss caused by pressure heating mechanism during sheath collapse, as shown in figures 5(a2)–(c2). Electrons gain energy in two main ways, that is, energy gain caused by pressure heating mechanism during sheath expansion and that caused by Ohmic heating mechanism during entire RF cycle. In addition, the Ohmic heating reaches maximum near the sheath boundary when the sheath begins to expand or the sheath collapse is about to end. Finally, with the increase of $B_x,$ the level of all heating decreases. Interestingly, when a small $B_x$ of 10 G is introduced, the Ohmic heating in front of the sheath boundary is enhanced when the sheath begins to expand, and this phenomenon disappears with the further enhancement of $B_x,$ as shown in figures 5(a3)–(c3).

Figure  5.  The spatio-temporal evolutions of both the total heating and different heating components at different $B_x$ in one RF periodic. (a) $B_x=0$ G, (b) $B_x=10$ G, (c) $B_x=20$ G.

DownLoad: Full-Size Img PowerPoint

To further analyze the influences of $B_x$ on the heating, figure 6 shows the evolutions of both the spatio-temporally averaged total electron heating (includes $〈\mathit{E}·\mathit{J}〉$ and $P_{\mathrm{abs}}$ ) and each heating component. It can be seen that the total heating decreases with the increased $B_x,$ and tends to be stable, which is consistent with the change in $n_{\mathrm{e}}$ shown in figure 3. Meanwhile, the decreased total heating is correlated highly with the weakened Ohmic heating, as shown in figure 6(b). In addition, as $B_x$ keeps increasing, the level of the Ohmic heating decreases gradually while that of the pressure heating shows an opposite trend, which implies a transition in the dominant heating mechanisms from Ohmic heating to pressure heating when $B_x$ increases gradually.

Figure  6.  The evolutions of the spatio-temporally averaged heating components with different $B_x.$ (a) Pressure heating, (b) Ohmic heating.

DownLoad: Full-Size Img PowerPoint

On the other hand, figure 6 also shows how $B_x$ affects the heating components. It can be seen that when $B_x=0,$ $P_{\mathrm{press}}$ and $P_{\mathrm{ohmic}}$ are zero in the y-direction. When $B_x>0,$ however, net heating in the y-direction appears, which arises from $B_x$ -induced Hall electric field $E_{\mathrm{H},y}.$ In addition, as $B_x$ increases, $P_{\mathrm{press},x}$ and $P_{\mathrm{press},y}$ show an upward trend, while $P_{\mathrm{press},z}$ shows a downward trend. It is clear that the varying pressure heating is correlated highly with $P_{\mathrm{press},x}.$ Meanwhile, $P_{\mathrm{ohmic},x}$ changes little with $B_x,$ while $P_{\mathrm{ohmic},y}$ first increases and then decreases, and $P_{\mathrm{ohmic},z}$ shows a downward trend. The varying Ohmic heating is correlated highly with $P_{\mathrm{ohmic},z},$ and the weakened $P_{\mathrm{ohmic},z}$ is attributed to the reduced Hall currents along 'E × B' direction [27]. Here the 'E × B' direction is referred as a direction which is perpendicular to both $E_x$ and $B_y.$ Now, by analyzing the level of different heating components, it is confirmed that $B_x$ affects the electron heating mainly by changing pressure heating along the longitudinal direction and by changing Ohmic heating along the 'E × B' direction.

The influences of $B_x$ on plasma characteristics are also explored, as shown in figure 7, which shows the evolutions of $T_{\mathrm{e}}$ as well as that of global EEDF at different $B_x.$ As can be seen from figure 7(a), when a small (10 G in this work) $B_x$ is introduced, $T_{\mathrm{e}}$ in plasma zone increased sharply, and then decreased gradually with the continuously increased $B_x,$ which results in the spatial profile of $T_{\mathrm{e}}$ evolving from bimodal-like to inverted U-like and then back to bimodal-like.

Figure  7.  Evolutions of $T_{\mathrm{e}}$ as well as that of the global EEDF versus different $B_x.$ (a) Electron temperature distribution, (b) global EEDF.

DownLoad: Full-Size Img PowerPoint

Meanwhile, the global EEDF shows a similar change. When a small $B_x$ of 10 G is introduced, the proportion of a considerable part of low energy (about 2.5–15 eV in this figure) electrons in the EEDF increases, which means that electrons become 'hotter'. As $B_x$ further increases to 40 G, however, the proportion of most low energy (about 2.5–9 eV) electrons decreases significantly, while the proportion of some lower energy electrons (about 0–2.5 eV) increased slightly, as shown in the blue dotted line in figure 7(b). As a result, the electron temperature decreases.

To understand the dramatic change in both $T_{\mathrm{e}}$ and EEDF caused by the small $B_x,$ figure 8 shows evolutions of loss rates of both the spatio-temporally averaged electron density and the spatio-temporally averaged electron absorbed power at different $B_x,$ where the benchmark data are the plasma parameters at $B_x=0.$ It can be seen that when $B_x>10$ G, the two loss rates rise at almost the same speed and tend to be stable. However, when $B_x\le 10$ G, the change in the loss rate of the averaged absorbed power is not as intense as that of the averaged electron density. For example, when $B_x=10$ G, the loss rate of the averaged absorbed power rises from 0% to 9.113% (9.113%/10 G), while that of the averaged electron density rises from 0% to 22.14% (22.14%/10 G). In this case, electrons are lost quickly, while the total electron heating remains at a relatively high level, and results in most of the electrons gaining more energy. As a result, the electron temperature and proportion of some low energy electrons in EEDF increase.

Figure  8.  Evolutions of loss rates of both the spatio-temporally averaged electron density and the spatio-temporally averaged electron absorbed power versus different $B_x.$

DownLoad: Full-Size Img PowerPoint

5.   Conclusions

In summary, the electron heating characteristics of a magnetic enhancement CCP in presence of both longitudinal and transverse uniform magnetic field have been explored by using a combination of theoretical and numerical calculations. It is found that the electron density as well as level of the total heating decreases gradually and tends to be stable with a continuously increased longitudinal magnetic field, which results from the confinement of electrons in the direction which is perpendicular to both the driving electric field and the applied transverse magnetic field (a direction referred to as 'E × B' direction in this work). Further analysis shows that longitudinal magnetic field affects the heating mainly by changing pressure heating along the longitudinal direction and by changing Ohmic heating along the 'E × B' direction. As longitudinal magnetic field increases continuously, the level of the Ohmic heating decreases gradually, and the pressure heating tends to be the dominant heating mechanism. Finally, it is worth mentioning that when a small longitudinal magnetic field is introduced, the electron temperature as well as the proportion of some low energy electrons in the EEDF increases. This is because the small longitudinal magnetic field causes a significant decrease in electron density but remains a relatively high level of the electron heating, which results in more energy per electron. Although the analyses in this paper are based on 1D simulations, these research conclusions, however, are still expected to be effective in 2D cases [34], and can provide guidance for optimizing the magnetic field configuration of some discharge systems with both transverse and longitudinal magnetic fields introduced.