Plasma Dynamics and Electron Transport in a Hall-Thruster-Representative Configuration with Various Propellants: II—Effects of the Magnetic Field Topology

Abstract

We investigate the effects of the magnetostatic ($B$) field topology on the plasma behavior in a 2D collisionless simulation setup that represents an axial–azimuthal cross-section of a Hall thruster. The influence of the $B$-field topology is assessed in terms of two principal design properties of the field in a typical Hall thruster, i.e., the field’s peak intensity along the axial direction, and the field’s axial distribution. The effects of the field’s intensity are investigated for three propellants—xenon, krypton, and argon. Whereas, the effects of the axial profile of the magnetic field are studied only for the xenon propellant as an example. We primarily aim to understand how the changes in the $B$-field topology affect the spectra of the resolved instabilities as well as the electrons’ transport characteristics and the contributions of various momentum terms to transport. The numerical observations on the instabilities’ characteristics are compared against the relevant existing theories to determine the extent to which the simulated and the theoretically predicted characteristics are consistent across the studied parameter space. It was, most notably, found that modes related to ion acoustic instability are dominantly present across the simulation cases. The ion transit time instability additionally develops at the highest $B$-field intensities as a long-wavelength structure. The main influence of the axial profile of the $B$ field on the plasma discharge was observed to be in terms of the electrons’ transport characteristics. Where possible, the insights from the simulations are discussed with respect to the relevant experimental observations available in the literature.

1. Introduction

Hall thrusters comprise a class of plasma propulsion technologies for spacecraft. Hall thrusters ionize a neutral, often inert, propellant gas injected into their discharge channel to create plasma and accelerate the plasma ions to generate thrust [1]. The ionization process is typically via DC electron bombardment and the acceleration process is electrostatic, as the ions are mostly influenced by an electric field [1,2]. Compared to the other “electrostatic” electric propulsion (EP) technologies on the market, such as Gridded Ion Engines, Hall thrusters feature a relative ease of manufacturing and offer a unique combination of favorable operational characteristics, including high thrust-to-power ratio, and high operational versatility. These features have resulted in a great deal of interest in the Hall thruster technology over the past two decades. Today, Hall thrusters are prime candidates for the next generation of space applications, from on-orbit servicing near the Earth to exploration missions to the Moon and Mars [3].

Since the introduction of Hall thrusters to the Western world in the early 1990s, this plasma technology has been the subject of rigorous and continued academic research in order to understand the physics of plasma associated with a Hall thruster’s operation. The scientific interest in Hall thrusters is because of several reasons: first, the plasma configuration in these devices, which consists of a perpendicular orientation of the applied electric ($E$) and magnetic ($B$) fields and is often referred to as a cross-field or $E$ × $B$ configuration, features rich underlying physics that is of great interest across various plasma science domains, owing, in part, to the presence of instabilities and turbulence that extend across a broad spectrum of spatial and temporal scales [2,4,5]. Second, the Hall thruster’s $E$ × $B$ configuration closely matches many other industrially relevant plasma sources, like magnetrons used in the manufacturing industry (e.g., for materials’ etching and coating). Third, a Hall thruster features relatively convenient accessibility and affordability for numerical and experimental physics studies, the insights from which can inform the research across the different plasma regimes and applications.

As the applied magnetic field is one of the two main forces acting on the plasma in a Hall thruster, the topology and the magnitude of the $B$-field plays a significant role in determining the plasma behavior and the operation and performance of the thrusters. For instance, the intensity and the distribution of the self-consistent electric field, which directly influence the acceleration process and, hence, the thrust efficiency, in a Hall thruster, is in fact tied to the efficacy with which the applied $B$ field can hinder the cross-magnetic-field motion of the electrons toward the anode element of the thruster [1]. Additionally, the overlap between the ionization and the acceleration zones within the Hall thruster domain is in part driven by the topology of the magnetic field. The extent of this overlap can affect the energy distribution of the ion species, hence influencing the efficiency of the acceleration process [1,2,6].

An important aspect of the physics of plasma in Hall thrusters that is not yet fully understood regards the influence of the $B$ field on the plasma instabilities and on the contribution that the instabilities in turn have to the phenomenon of “anomalous” electron cross-field transport [2,4,5,7,8,9,10]—a phenomenon that considerably limits the devices’ performance [1,2]. The complex, highly coupled nature of the plasma processes in a Hall thruster poses a great challenge in isolating and distinguishing the mechanisms/pathways through which the variations in the $B$ field may alter the underlying phenomena and the global behavior of the plasma. This is where parametric studies in controlled simulation environments become highly beneficial. Through reasonable simplifications to the overall complex picture of the problem, and by performing high-fidelity parametric simulations in the simplified setups, we will be able to set the foundational knowledge on the effects of the parameters of interest, such as the $B$ field in this case. Building upon the insights derived from the simplified parametric simulations, we can then methodically progress toward setups that are increasingly more representative of the real world, ultimately discovering the parametric dependencies of plasma processes in realistic settings.

There are two crucial enabling elements here: (a) suitable simulation setups for parametric studies of simplified, controlled problem definitions, and (b) suitable high-fidelity simulation tools that allow cost-effective extensive parametric investigations. Regarding (a), the available benchmark cases in the literature [11,12] serve as viable solutions. This is because of two reasons: (1) these benchmarks represent well-accepted cases within the community, and (2) they enable assessing the parametric variations of the underlying processes with respect to the rigorously studied baseline settings of the benchmark problem.

Concerning (b), following the discussions in part I [13], the reduced-order particle-in-cell (PIC) scheme [14,15] represents a breakthrough. This is because the remarkable cost-effectiveness of this PIC scheme compared to the traditional PIC methods, which has been extensively verified in various cases [14,15,16,17,18,19], allows exploring wide parameter spaces within practical timeframes.

Noting the above remarks, in this part II paper, we have once again adopted a simulation setup largely similar to that of the axial–azimuthal Hall thruster benchmark reported in Ref. [11] and have performed a series of parametric studies using Imperial Plasma Propulsion Laboratory’s (IPPL’s) reduced-order quasi-2D kinetic PIC code, IPPL-Q2D. IPPL-Q2D has been used to verify the reduced-order PIC scheme in various 2D setups [14,18,19] and has already been deployed for a broad range of physics studies in our prior publications [20,21,22,23]. Our quasi-2D axial–azimuthal parametric simulations here are performed using an approximation order of the 2D problem that corresponds to 40 “regions” ($M$) along the axial direction and 20 horizontal regions ($N$) along the azimuthal direction. As was discussed in detail in part I of this article [13], this specific approximation order (domain decomposition) offers a reasonable balance between the computational gain and the accuracy of the simulations.

The main objective of this part II paper is to evaluate the variations in plasma behavior due to the changes in the properties of the $B$-field topology, namely, the field’s peak intensity along the axial direction and the field’s axial gradients. Concerning the latter, as will be seen in Section 2, we have parametrically changed the $B$-field’s axial gradients (profile) as a way to mimic varying degrees of overlap between the ionization and the acceleration zones to isolate the corresponding potential effects on the discharge behavior. The variations in the plasma behavior are assessed in terms of the macroscopic plasma properties, the characteristics of the plasma fluctuations, the electrons’ cross-field transport, and the properties of the ion species. With regard to the variations in the characteristics of the plasma fluctuations, we additionally aim to assess the degree to which the theories of the pseudo-saturated electron drift instability (EDI)/ion acoustic instability (IAI) are consistent with the simulated variations in the characteristics of the resolved dominant instability modes. Finally, similar to part I, we have also been interested in evaluating the impacts that the ion mass has on a collisionless plasma discharge. This aspect is assessed for the simulations with various peak intensities of the $B$-field using three propellant gases—xenon, krypton, and argon.

In terms of the positioning of this work in the context of the previous research, our methodologies and results in this paper build upon a series of publications that studied the physics of the plasma instabilities and the electrons’ transport along the axial–azimuthal coordinates of a Hall thruster [7,8,11,24,25,26,27,28,29,30,31,32]. In particular, we would mention the work by Charoy et al. in Ref. [26], who adopted an axial–azimuthal simulation setup with an imposed ionization source to assess how the changes in the peak intensity of the magnetic field in a xenon plasma discharge affect the characteristics of the electrons’ transport, including the contributions of various force terms in the electrons’ azimuthal momentum equation. Furthermore, Petronio [31] and Reza et al. [32] carried out axial–azimuthal simulations of different Hall thruster geometries with a self-consistent description of the ionization process through the Monte Carlo Collisions (MCC) scheme [33], and by treating the neutral species as a fluid described by the 1D Euler system of equations [34]. Although pursuing different approaches and numerical tools, both studies were aimed primarily at evaluating the impacts of the operating condition, including the magnetic field peak intensity, on the global operation of the Hall thruster and the breathing mode characteristics [31,32].

The present paper goes beyond the prior relevant publications in the literature in three main respects: (i) the broad extent of the investigated parameter space, which is primarily characterized by several values of the quantities that define the magnetic field topology (Section 2), (ii) the extensive comparisons of the simulated characteristics of azimuthal instabilities against the theoretically predicted characteristics for different magnetic field topologies and ion masses (with various $B$-field intensities), and (iii) the adoption of the novel Dynamic Mode Decomposition (DMD) technique [35,36] to reliably identify the dominant instability modes for each simulation case and to simultaneously derive the spatial–temporal characteristics of the dominant modes.

2. Overview of the Simulations’ Setup and Conditions

The setup of the simulations was described at length in part I of the article [13]. The simulation domain is a 2D Cartesian ($x-z$) plane representative of an axial–azimuthal cross-section of a Hall thruster, with $x$ along the axial direction and $z$ along the azimuth. The computational and physical parameters, including the domain’s lengths along the axial and azimuthal directions, the cell size, the timestep, the total simulation time, and the initial plasma conditions are exactly as mentioned in part I [13]. The discharge voltage is $200 V$ and the imposed ion current density is $400 A{m}^{-2}$ corresponding to the baseline conditions.

An important point to emphasize concerning the simulations’ setup here is that the electrons’ boundary condition at the cathode side of the domain follows the quasi-neutrality approach [15,25]. This is the main difference of our setup compared to that of the axial–azimuthal benchmark in Ref. [11] and the work of Charoy et al. in Ref. [26], for both of which the electrons’ cathode boundary condition was based on the current-continuity approach [15,25]. As was detailed in part I [13], we adopted the quasi-neutrality cathode condition in this effort because this boundary condition is less sensitive to the chosen temperature of the electrons re-injected into the domain from the cathode boundary to sustain the discharge [25].

Considering the objectives of this part II paper as described in Section 1, we have shown in Figure 1 and Figure 2 the axial profiles of the radial magnetic field used for various simulation cases. Figure 1 illustrates how the variation in the peak intensity of the field from 5 to 30 mT changes the $B$-field profile imposed in the simulations.

Figure 2a,b depict the axial distributions of the $B$-field for the simulations with varying axial gradients of the magnetic field in terms of the location of the field’s peak value ($x_p$) [plot (a)] and the standard deviations ($\sigma$) in the Gaussian field profile (Equation (1)) [plot (b)]. The axial profile of the imposed ionization source [11,24,25] in the simulations is shown on both plots. This illustrates the extent of overlap across various simulation cases between the ionization zone, as defined by the ionization source’s distribution, and the acceleration zone, which typically corresponds to the region of the highest magnetic field intensities.

$$B\left(x\right)=A_k\exp \left(-{\displaystyle \frac{{\left(x-x_p\right)}^2 }{2{\sigma }^2}}\right)+C_k , k=\mathrm{1,2}$$

(1)

In Equation (1), $k=1$ for $x\le x_p$, $k=2$ for $x>x_p$, and $\sigma =0.625 cm$. $A_k$ and $C_k$ are constant values, which are obtained given $B\left(0\right)=6 mT$, $B\left(L_x\right)=1 mT$, and $B\left(x_p\right)=10 mT$. Referring to Figure 2b, we observe that, by increasing the value of $\sigma$, the $B$-field profiles gradually approach the shape given by the case with $\sigma =inf$, which is associated with an arbitrarily large standard deviation value (set here to $100$). The baseline value of $\sigma$ (${\sigma }_0$) is $0.625 cm$, which is the same value as that of the benchmarking simulation case in Ref. [11].

It is noted that $A_k$ and $C_k$ coefficients are determined in such a way as to constrain the values of the $B\left(x\right)$ profile at the two boundaries of the domain, which are different from $B\left(x_p\right)$. As a result, even at very large values of $\sigma$ (resembling $\sigma$ approaching infinity), the $B(x$) profile is not flattened. Also, while we refer to σ in this case as approaching infinity, in practice, its value is not infinite but rather set to a relatively large number, such as 100 here. In fact, the shape of the profile does not change appreciably for values of σ beyond a sufficiently large value, e.g., 10.

Before moving on, we would emphasize that the simulations in this work do not account for the radial direction and the wall losses. As a result, the location of the channel’s exit plane, which is referred to throughout the text and in the figure captions, merely represents the location of the peak of the magnetic field profiles, except for the simulations with different $x_p$ values. The terms inside the “channel” and in the “plume”, hence, refer to the axial locations before and after the B-field’s peak, respectively, for the simulations with various field intensities and various $\sigma$ values.

3. Results

In this Section, we look at the variations that changing the properties of the magnetic field topology induces in the time-averaged macroscopic plasma properties and the various current terms, namely, the discharge current, the electron current, and the ion current. The results related to the effects of the $B$-field’s peak intensity are presented in Section 3.1, whereas Section 3.2 shows how the axial profile of the $B$-field affects the time-averaged properties of a plasma discharge with xenon propellant.

3.1. Effects of Varying the Magnetic Field’s Peak Intensity on Macroscopic Plasma Properties

Figure 3 illustrates the axial distributions of several plasma properties from the simulations with different magnetic field peak intensities ($B$). The profiles of the plasma properties are averaged over the time interval of $20–30 \mu s$, which is after the simulations arrive at the quasi-steady state.

We observe from Figure 3 that, for the magnetic field peak intensities in the range of $10$ to $20 mT$, the plasma profiles exhibit minor variations. However, the cases with the lowest and the highest magnetic field peak values show different behaviors. At $B=5 mT$, the density peak is much higher than in the other cases. Also, the acceleration zone, as represented by the region of the highest axial electric field values, is widened and shifted further downstream, resulting in a relatively lower overlap with the ionization region. This, in turn, causes the ions introduced into the domain in the ionization zone to experience a larger portion of the potential drop, leading to their acceleration to higher velocities on average and, consequently, lower ion densities in the plume. The observed behavior in this case is most likely due to weaker azimuthal instabilities and, hence, lower instability-induced electrons’ axial transport upstream of the channel’s exit. The downstream shift in the acceleration region at low magnetic field intensities is also observed experimentally [37,38,39]. In Ref. [39], it is reported that, at low $B$-field intensities, the high $E_x$ region occurs entirely outside the channel, and it shifts upstream as $B$ increases. Chaplin et al. [37] and Hargus et al. [38] also observed that most of the potential drop happens outside the channel. These experimental observations are in line with our numerical results.

At $B=30 mT$ (and, to a lower extent, at $B=20 mT$), the discharge becomes rather turbulent (as will be demonstrated later in Section 4.2 and Section 4.3). At this $B$ level, the magnetic field barrier is very strong and, thus, it severely limits the electrons’ axial transport. As was also observed by Boeuf and Garrigues [24], in the conditions where the electron conductivity is insufficient to balance the ion beam current, axial oscillations corresponding to ion transit time instability (ITTI, refer to Section 4.3) are developed. In the 30 mT case, the turbulent nature of the discharge affects the acceleration of the ions, which is reflected in the density profiles in the plume. From the $E_x$ and the $n_i$ profiles in Figure 3, it appears that these effects are stronger for argon than for xenon and krypton. These effects were observed to significantly impact the electric potential, which included the appearance of a small bump in the potential profile downstream of the exit plane. This small bump manifests itself as a dip in the $E_x$ profiles. The turbulence together with a long-wavelength azimuthal instability, which will be discussed in Section 4.2 and Appendix A, seems to be causing substantial electron heating, leading to higher electron temperatures in the 30 mT case compared to the cases with lower $B$ values.

Concerning other notable insights from Figure 3, the peak $E_x$ is generally lower at higher $B$-field intensities ($B\ge 15 mT$) compared to the peak $E_x$ at $B=10 mT$ (except for argon at $B=20 mT$). Additionally, the width of the high $E_x$ region (acceleration zone) increases from $B=10 mT$ to $B=30 mT$. These behaviors are consistent with the experimental observations. In this regard, Hargus et al. [38] observed that the peak of the electric field is higher for the reduced $B$-field condition. Gawron et al. [40] noticed that the acceleration zone expands with increasing $B$-field intensity. They explained that the expansion of the acceleration zone is due to the enlargement of the region over which the $B$-field is strong enough to effectively confine electrons [40]. Contrary to this explanation, in Section 3.2, we will show that increasing the extent of the magnetized region does not appreciably alter the width of the acceleration zone.

Note that the electron temperatures predicted by the model represent an overestimation of the real values in a typical Hall thruster. This is due to not capturing the energy losses to the wall as well as the ionization energy losses.

The plots in Figure 4 present the variations in the mean values of the ion number density and the ions’ axial drift velocity ($V_{ix}$) in the plume region ($1.5 cm<x<2.5 cm$) with the magnetic field peak intensity for the three studied propellants. On the $V_{ix}$ plot, the superimposed horizontal dashed lines show the theoretical (ideal) values of the ions’ axial drift velocity (exhaust velocity) from Equation (2) for the different propellants and the $B$ intensities.

$$v_{exhaust}=\sqrt{{\displaystyle \frac{2e{V}_d}{m_i}}}.$$

(2)

In Equation (2), $e$ is the elementary charge, and $m_i$ is the ion mass. It is recalled from the discussions in part I [13] that the kinetic energy losses in the ion population due to, for instance, an overlap between the acceleration and the ionization zones or the interactions between the ions and the instabilities’ waves can cause a deviation in the simulated ions’ axial drift velocity from the ideal exhaust velocity given by Equation (2). The sheath established near the anode (appearing in Figure 3 as a bump followed by a dip in the $E_x$ profiles upstream of the channel’s exit) acts an opposing factor, which can partly compensate for the acceleration losses by increasing the potential drop experienced by the ions during acceleration [13].

With the above in mind, we observe from the $V_{ix}$ plot that, for xenon and krypton and for the $B$ intensities of 5 and 10 mT, the simulated $V_{ix}$ is quite consistent with the theoretical $v_{exhaust}$. This is because, in these conditions, the acceleration inefficiencies, mainly due to the overlap of the acceleration and the ionization zones, are almost balanced against the increased potential due to the formation of a rather large anode sheath (as seen in Figure 3).

As $B$ increases to 15 and 20 mT, the deviation between the $V_{ix}$ and the $v_{exhaust}$ grows. A major discrepancy between the simulated $V_{ix}$ and the ideal $v_{exhaust}$ occurs at $30 mT$. These deviations are increasing substantially from xenon to argon. The disparity at $30 mT$ is believed to be due to the turbulent nature of the discharge under this extreme $B$ condition, which distorts the electric potential profile and leads to an incomplete acceleration of the ions. Nevertheless, in practice, operating at higher magnetic field intensities often corresponds to higher operating voltages. The results here illustrate the consequences of not increasing the voltage when operating at higher magnetic fields.

Figure 5 shows the variations with the $B$-field intensity of the several current ratios, namely, the ratio of the electron current toward the anode to the nominal total ion current density ($I_0$) (left-hand-side plot), the ratio of the ion current toward the cathode to $I_0$ (middle plot), and the ratio of the electron current to the ion current (right-hand-side plot). More detailed definitions of the electron and the ion current terms, as well as how these were obtained from the simulations, can be found in part I [13].

Starting with the variations trend of the normalized $I_e$ vs. the $B$-field intensity, we expected to see that the electron current decreases with $B$ for all the propellants due to stronger electron confinement. It is also noticed that the most notable changes in the $I_e/I_0$ occur between $5$ and $10 mT$ and between $15$ and $20 mT$. $I_e/I_0$ seems to have reached an asymptotic value between 2$0$ and $30 mT$ for xenon and krypton. However, it shows a slight increase from $20$ to $30 mT$ for argon. The largest disparity in $I_e/I_0$ among the propellants occurs at $B=5 mT$. Next, the $I_i/I_0$ has an increasing trend with the $B$ intensity until it almost reaches a plateau around $15 mT$ for all the propellants, indicating a threshold level for the fraction of the ions that can be extracted at the cathode boundary. The variation trend of $I_e/I_i$ vs. $B$ is very similar to that of $I_e/I_0$ across the three propellants.

3.2. Effects of Varying the Axial Profile of the Radial Magnetic Field on Macroscopic Properties of Xenon Plasma

Figure 6 shows the changes in the time-averaged axial profiles of the $n_i$, the $E_x$, and the $T_e$ for different axial distributions of the radial magnetic field. Additionally, Figure 7 presents the corresponding variations in the spatiotemporal mean of the $n_i$ and the $V_{ix}$ in the plume region.

A couple of interesting observations can be made from Figure 6. First, changing the standard deviation of the $B$-field’s profile over a rather large range is seen to minimally affect the axial profiles of the plasma properties. The case with $\sigma$ = 0.25${\sigma }_0$ presents the most difference among all the cases, particularly in terms of the ion number density and the axial electric field. The plasma profiles are almost invariant with $\sigma$ from $0.5{\sigma }_0$ to $inf$. These observations are significant because they highlight that, when the ionization process is modeled through a temporally invariant source, the variations in the width of the region of the highest magnetic field intensities do not majorly affect the distributions of the plasma properties. This contradicts the explanations given in Ref. [40] to justify the experimental observation that the acceleration zone expands with the increasing $B$-field intensity.

As a result, the often-reported significant alterations of the plasma profiles and the performance of a Hall thruster when varying the width of the B-field profile in a real thruster setting [41,42,43] may primarily be attributed to the variations in the ionization process, including the efficiency of the ionization and the relative location of the ionization zone with respect to the acceleration zone. Alternatively, variations in the axial gradient of the magnetic field can change the curvature and the radial gradients of the field topology. As was demonstrated in Ref. [21], the $B$-field’s curvature can significantly influence the plasma processes and conditions. Hence, the changes in the curvature of the magnetic field could also be another contributing factor to the experimentally observed changes in the performance when varying the width of the high magnetic field region [41].

Second, the changes in the axial location of the peak $B$-field intensity are evident to have a more notable impact on the profiles of the plasma properties. In this respect, we see that, when $x_p$ is varied from 0.5 cm to 1.5 cm, the extent of the overlap between the ionization zone and the acceleration zone is continuously reduced as the acceleration zone shifts progressively downstream toward the plume (as inferred from the $E_x$ plot in Figure 6). This downstream movement of the acceleration zone leads to consistently larger ion velocities in the plume (see Figure 7, bottom row), as well as an increase in the peak $n_i$ value and a gradual reduction in the peak $T_e$ values together with a downstream shift in the axial location of the peak $T_e$ (Figure 6).

Nevertheless, the extent of the downstream movement of the acceleration zone does not match the extent of the shift in the location of the magnetic field peak. This is, in part, because the $B$-field value is fixed on the left boundary, thereby the $B$-field intensity is maintained at relatively high values upstream of the acceleration zone. As a result, the region of high $B$-field intensity does not significantly move with the change in $x_p$.

Despite the minor changes observed in the plasma profiles, the electron current ($I_e$) and, thus, the discharge current ($I_d$) vary substantially with changes in the $B$-field profile, as is shown in Figure 8. While the ion current ($I_i$) remains almost invariant with $\sigma$ and $x_p$, the electron current decreases significantly as either $\sigma$ or $x_p$ increases. Indeed, as will also be seen in Section 4.1, the consistent increase in the width of the high magnetic field region with increasing $\sigma$ or $x_p$ leads to a continuous reduction in the overall electron transport.

4. Discussions

In the following Section, we present and discuss the outcomes of the more detailed analyses we carried out on the results of the various parametric simulations. Section 4.1 delves into how the variations in the peak magnitude of the magnetic field as well as in its axial profile change the electrons’ transport and the contributions of the different force terms in the electrons’ azimuthal momentum equation to this phenomenon.

In Section 4.2, using 2D fast Fourier transform (FFT) analyses and, informed by the 1D FFT plots and the DMD spectral analyses [35,36] in Appendix A, we look at the changes in the characteristics of the resolved fluctuations with the $B$-field’s peak intensity. In this Section, where relevant, we compare the characteristics of the main numerically resolved fluctuations with the corresponding theoretical predictions from the available instability theories (Appendix B).

Finally, in Section 4.3, we discuss how the changes in the intensity of the $B$-field affects the ion population in the plasma discharges of the three studied propellants. We focus on the effects of the $B$-field intensity only since the variations in the axial profile of the magnetic field were not found to majorly influence the ion population. The $B$-field intensity’s impacts on the ion population for the different propellants are analyzed by assessing the ions’ distribution functions, the axial profiles of the ion temperature, and the time evolution signals of the ion current.

4.1. Variations in the Electrons’ Transport and the Contributing Momentum Force Terms

Here, we first assess, for the three studied propellants, how the peak intensity of the $B$-field changes the electrons’ transport characteristics. Then, for the xenon propellant as an example, we evaluate the effects that the axial profile of the $B$-field in terms of the $\sigma$ and the $x_p$ parameters have on the transport.

As was discussed in part I [13], the electrons’ azimuthal momentum equation can be presented in terms of a balance equation—shown by Equation (3)—between the magnetic force ($F_B$) on the left-hand side, and the sum of four momentum terms on the right-hand side, namely, the temporal inertia force ($F_t$), the convective inertia force ($F_I$), the viscous force ($F_{\Pi }$), and the electric force ($F_E$). These terms represent various contributing factors to the electrons’ axial transport. The definition of these force terms is given by Equations (4)–(8).

$$F_B=F_t+F_I+F_{\Pi }+F_E.$$

(3)

$$F_B=-q{n}_e{v}_{e,x}B,$$

(4)

$$F_t={\displaystyle \frac{\partial }{\partial t}}\left(m{n}_e{v}_{e,z}\right),$$

(5)

$$F_I={\displaystyle \frac{\partial }{\partial x}}\left(m{n}_e{v}_{e,x} v_{e,z}\right),$$

(6)

$$F_{\Pi }={\displaystyle \frac{\partial }{\partial x}}\left({\Pi }_{e,xz}\right),$$

(7)

$$F_E=-q{\tilde{n}}_e{\tilde{E}}_z.$$

(8)

In Equations (4)–(8), $q$ is the unit charge, m is the electron mass, $n_e$ is the electron number density, $v_{e,x}$ and $v_{e,z}$ are the electrons’ axial and azimuthal drift velocities, and ${\tilde{n}}_e$ and ${\tilde{E}}_z$ are the number density and the azimuthal electric field fluctuations, respectively.

The time-averaged axial profiles of the force terms in Equation (3) are presented in Figure 9 for the various $B$-field intensities and the three studied propellants. The $F_t$ profile is not shown because it was found to be negligible across all the simulations. Figure 10 shows the variations in the peak values of the $F_B$ and the $F_E$ terms with the $B$ intensity as well as the variations in the axial locations of these force terms, respectively denoted as $x_{F_{B,max}}$ and $x_{F_{E,max}}$.

The axial profiles of the force terms in Figure 9 and the variation trends of the $F_B$ peak and the $F_E$ peak in Figure 10 show a large degree of variability across $B$ and ion mass, not showing distinct, generally applicable trends.

In particular, the cases with $B=5 mT$ and $B=20–30 mT$ represent two extremes in terms of the characteristics of the momentum force terms. At 5 mT, all the force terms—$F_I$, $F_{\Pi }$, and $F_E$—have comparable magnitudes. The $F_{\Pi }$ term is dominant near the exit plane and exhibits a large peak, whereas the peak of the $F_E$ occurs further downstream in the plume. The profiles of force terms in this case are nearly similar to those of the $V=600 V$ case presented in part I [13].

The $F_{\Pi }$ term is less significant at the high B-fields ($B=20–30 mT$). The force terms at $B=30 mT$ condition, moreover, exhibit the largest disparity among themselves across the different propellants. Notably, in the case of argon, the $F_E$ term has disproportionately large values compared to the other gases at $B=30 mT$.

From Figure 10, as $B$ increases, the $x_{F_{B,max}}$ and the $x_{F_{E,max}}$ gradually shift upstream for the xenon propellant, with the $x_{F_{B,max}}$ moving continually further away from the channel’s exit plane. For krypton and argon, the upstream movement of the $x_{F_{B,max}}$ and the $x_{F_{E,max}}$ reverses downstream after $B=15–20 mT$.

Now, referring to Figure 11, we discuss how the variations in the axial profile of the $B$-field affect the electrons’ transport characteristics. Looking at the left-column plots, it is noticed that as the parameter $\sigma$ increases, the axial distribution of the $F_B$ (or equivalently the RHS term, which is the sum of the $F_E$, $F_{\Pi }$ and $F_I$ terms), tends to become more uniform. This behavior is related to the consistent reduction in the axial gradient of the $B$-field’s profile. Additionally, the peak magnitude of the $F_B$, along with the peak of the individual contributions ($F_E$, $F_{\Pi }$ and $F_I$), all consistently decrease with increasing $\sigma$.

Considering the plots in the right column of Figure 11, when the axial location of the peak $B$-field intensity moves from $0.5 cm$ to $1.5 cm$, the profile of the $F_B$ is seen to become increasingly more flattened, with its upstream values dropping notably. At $x_p\ge 1 cm$, the $F_B$ profile begins to develop a peak at the $x_p$ position, which is nearly as prominent as the upstream peak. For the large $x_p$ values ($x_p\ge 1.25 cm$), the $F_E$ term as well as the $F_{\Pi }$ and the $F_I$ terms show oscillatory axial distributions across the domain.

4.2. Numerical and Theoretical Analyses of the Instabilities’ Characteristics

Here, we analyze how the characteristics of the resolved fluctuations in our simulations vary with the peak intensity of the magnetic field and compare the numerical trends with the theoretical ones. These analyses are performed for the three propellants. We discuss the effects of the $B$-field peak intensity only because this is the relevant physical parameter that yields meaningful insights regarding the variations in the instabilities’ characteristics.

To begin the discussions, we compare in Figure 12 the numerical dispersions of the azimuthal electric field fluctuations from the various $B$ simulations against the theoretical nonlinear dispersion relation of the IAI (given by Equation (A1) of Appendix B) at two axial locations, one inside the channel and one in the plume.

From Figure 12, for $B\le 15 mT$, the numerical dispersion of the instabilities inside the channel (at $x=0.375 cm$) matches quite well the IAI’s nonlinear dispersion relation. At high $B$-field intensities ($B\ge 20 mT$), however, the dispersion map is localized around a single mode with a small wavenumber and a relatively low frequency.

In the plume (at $x=1.625 cm$) and for $B\ge 10 mT$, the linear dispersion relation of the ion sound waves is overall well fitted to the instabilities’ dispersion maps. Nonetheless, there is an increase in the instabilities’ dispersion across a wider range of frequencies and azimuthal wavenumbers as the $B$ increases, implying the co-existence of several modes. At $B=5 mT$, the instabilities’ dispersion in the plume is close to the nonlinear dispersion relation of the IAI rather than the linear one.

The 1D spatial and the 1D temporal FFT plots of the azimuthal electric field from the simulations with various magnetic field intensities and ion masses are provided in Appendix A. We have used these spatial and temporal FFT spectra to derive the azimuthal wavenumber and the frequency of the dominant instability modes at the two investigated axial locations. The FFT spectra in Appendix A also exhibited features that demonstrate the turbulent nature of the plasma at $30 mT$, in particular, but also at $20 mT$.

In Figure 13, the numerically obtained characteristics of the dominant modes as described above are compared against the theoretical $k_z$ and $f$ (frequency in Hz) values associated with the fastest-growing ion acoustic mode (Equation (A3) of Appendix B) at two axial positions, one corresponding to the mid-location within the channel and the other midway in the plume. The theoretical data points in Figure 13 are calculated from Equation (A3) using the values of the time-averaged plasma properties from the simulations at the probed axial locations.

From Figure 13, we primarily observe that, inside the channel, the $k_z$ and the $f$ of the dominant mode from the simulations show stronger variations with the $B$ intensity than what is predicted by the theory. The closeness of the simulated characteristics of the instabilities to the theory varies significantly as well across various $B$-field values and the different propellants.

Moving from inside the channel into the plume, we can notice a variation in the characteristics of the instabilities. In this regard, Equation (A3) suggests that the azimuthal wavenumber of the IAI’s fastest-growing mode is proportionally correlated with the $n_e/T_e$ ratio. Accordingly, the most remarkable change in the $k_z$ is observed at $B=5 mT$, which corresponds to the largest variation in the $n_e/T_e$ ratio between the channel and the plume regions according to Figure 3. This is in line with the relevant observations reported in Ref. [25].

It is also interesting to observe that, at $B=30 mT$, the wavenumber of the instabilities for xenon increases significantly from within the channel to the plume, whereas, for the cases of krypton and argon, the wavenumber remains almost constant at this $B$ value.

Looking at the ions’ azimuthal phase-space plots in Figure 14, we overall notice the trapping of the ions in the potential field of the azimuthal fluctuations. These plots suggest that ion-wave trapping is a saturation mechanism for the instabilities across almost all the simulation cases. However, at the highest field intensity ($B=30 mT$), the waves seem to not be strong enough to significantly affect the ions’ azimuthal velocity distribution.

To conclude the discussions here, Figure 15 compares the numerical and the theoretical predictions of the variations with the $B$ intensity in the RMS amplitudes of the azimuthal electric field fluctuations ($\delta E_{RMS}$) and the electron number density fluctuations ($\delta n_{e,RMS}$) within the plume (spatiotemporally averaged across $x=1.25$ to $2.25 cm$ and over $20–30 \mu s$). The theoretical data points are obtained from Equation (A4) in Appendix B. As was pointed out in part I [13] as well, the theoretical RMS amplitudes correspond to the saturated amplitudes of the IAI calculated by assuming the ion-wave trapping as the saturation mechanism. Thus, the comparisons presented in Figure 15 show the extent to which these theoretical oscillation amplitudes are consistent with the numerical observations.

The salient point to elaborate on from the plots in Figure 15 is that the theoretical and the simulated $\delta E_{RMS}$ values (for all the propellants) exhibit different trends with the $B$ intensity. While the theoretical values for each propellant show a slight decrease from 5 to 10 mT and then an increasing trend for the higher $B$ intensities, the simulated $\delta E_{RMS}$ values exhibit an overall decreasing behavior with $B$. At high $B$ intensities ($B=20, 30 mT$) in particular, the substantial disparity could partly be due to the presence of a long-wavelength instability, especially upstream of the channel exit plane, and the impacts of this instability on the downstream waves.

The presence of the instabilities of long wavelength (with their wavelengths almost equal to the domain’s azimuthal size) is detected by the DMD analysis for the cases with $B=20$ and $30 mT$ upstream of the channel exit (see the DMD plots in Appendix A). The frequency of these instability modes is about $1–1.5 MHz.$ The spatial structure and the frequency of the modes are similar to those of the long-wavelength mode observed for the low ion-current-density conditions in part I [13] and are reminiscent of the ITTI [44,45,46,47].

4.3. Variations in the Distribution Function and the Properties of the Ion Species

For the last part of our discussions in this paper, we examine the variations in the ions’ velocity distribution function (IVDF) with the peak intensity of the magnetic field. To this end, we have presented in Figure 16 and Figure 17, respectively, the IVDFs in the plume integrated over the axial extent of $1.5 cm<x<2.5 cm$, and the 1D1V ions’ distribution functions along the axial direction and the axial velocity component.

The behaviors at $B=30 mT$ are particularly notable. Indeed, as was anticipated from the analyses of the results in the preceding Sections, the rather turbulent state of plasma at this $B$ intensity is seen from Figure 16 and Figure 17 to have majorly affected the distribution functions for all the propellants, especially those of the lighter mass gases—krypton and argon. The peaks of the IVDFs at $30 mT$ have shifted toward lower energies and the distributions feature strong broadening (Figure 16). Additionally, at $30 mT$, the ion beam is significantly perturbed as is noticed from Figure 17. The notable spread in the ions’ axial velocities as well as the strongly perturbed ion beam can be associated with the presence of a strong axial ITTI and the consequent “wave-riding” mechanism [24,44,46].

The substantial scattering of the ions at $B=30 mT$ leads to a strong heating of the ion population, which is evident from the time-averaged axial distributions of the ion temperature ($T_i$) in Figure 18. The heating effects become more pronounced from xenon to argon.

The observations made above concerning the influence of the instability spectra on the ion population are also evident from the time evolution plots of the ion current shown in Figure 19. Referring to this figure, for the $B$ intensities of $30$ and $20 mT$, which were consecutively associated with the highest levels of turbulence in the plasma behavior, and for which the presence of the ITTI waves was identified, the $I_i$ signals feature strong, high-amplitude periodic oscillations, in line with the observations reported in Ref. [24]. The amplitude of the oscillations in the $I_i$ for all the propellants gradually decreases with decreasing $B$ intensity from $15$ to $10 mT$ and, finally, the $I_i$ is almost steady for the $5 mT$ case.

5. Conclusions

We investigated in this paper how the variations in the peak intensity and the axial profile of the externally applied magnetostatic ($B$) field affect the plasma discharge in a collisionless axial–azimuthal Hall-thruster-representative simulation case. The effects of the $B$-field topology were assessed in terms of the changes in the macroscopic plasma properties, the instabilities spectra, the characteristics of the electrons’ cross-field transport, and the properties of the ion population. The influence of the $B$-field intensity was studied for three propellants, namely, xenon, krypton, and argon. The changes in the plasma behavior due to the axial profile of the magnetic field were evaluated only for the xenon propellant. The variations in the simulated characteristics of the dominant ion acoustic instability (IAI) related modes with the $B$-field topology were compared against the theoretical predictions of the characteristics of the fastest-growing IAI.

Increasing the peak intensity of the $B$-field implies that the electrons would “feel” an increased level of resistance against their cross-field motion toward the anode side of the domain. Accordingly, we observed that increasing $B$ intensity results in a continuous decrease in the $I_e$ until it reaches an asymptotic value between 20 and 30 mT. $I_i$; however, increases with $B$ until it nearly plateaus around 15 mT.

At the low $B$-field extreme ($B=5 mT$), the acceleration zone is widened and shifts downstream, reducing its overlap with the ionization zone. Furthermore, at this field intensity, the electrons’ transport within the channel is primarily driven by the $F_{\Pi }$ term in the electrons’ azimuthal momentum equation. Whereas, in the plume, the $F_{\Pi }$, $F_I$, and $F_E$ terms are of comparable magnitudes. In any case, the balance between of the $F_{\Pi }$ and the $F_I$ in the plume leaves the $F_E$ term as the dominant transport mechanism in this region.

At the high $B$-field extreme ($B=30 mT$, and to a lower extent at $B=20 mT$), the discharge becomes turbulent, leading to ions’ drift velocities much lower than the ideal $v_{exhaust}$. The non-idealities in acceleration substantially increase from xenon to argon. The turbulent fluctuations cause significant heating of the electron and ion populations and induce remarkably high electron transport in the case of argon. The $30 mT$ case exhibited the highest level of disparity in the plasma behaviors among the studied propellants.

Concerning the variations in the waves’ characteristics with the $B$-field intensity, the wavenumber and the frequency of the dominant modes are found to decrease with $B$ inside the channel, while they show no general, distinct trend in the plume. The most notable change in terms of the waves’ characteristics between the channel and the plume regions is seen to occur for the wavenumber and the amplitude of the oscillations at $B=5 mT$.

Finally, at the low $B$-field values (especially, $B=5 mT$), the ions have a highly beam-like velocity distribution. Whereas, at the high $B$ intensities (particularly, $B=30 mT$), the ions’ distribution functions exhibit “swirling” features as a result of the ions’ interactions with turbulent fluctuations.

Regarding the effects of the $B$-field’s axial distribution, the profiles of the plasma properties are minimally affected by varying the standard deviation ($\sigma$) of the defined $B$-field’s Gaussian profile. This limited sensitivity of plasma behavior to $\sigma$ is, in part, attributed to prescribing a fixed ionization source in the simulations. Our observations on this point suggest that the well-acknowledged significant changes in the plasma profiles and the performance of Hall thrusters with variations in the width of $B$-field’s profile [41,42,43] can be primarily due to the variations in the ionization process and/or the curvature of the magnetic field, in contrast to some previously suggested explanations in the literature [40].

Shifting the axial location ($x_p$) of the $B$-field peak intensity toward the plume reduces the overlap between the ionization and the acceleration zones, albeit to a lesser extent than the downstream shift in the $x_p$ itself. The impact of shifting the $x_p$ is less than may have been expected. This is believed to be because of the fixed and relatively high $B$-field values in the near-anode region of our simulations together with the imposed invariant ionization source.

Despite minor changes in the plasma profiles and an almost invariant ion current ($I_i$) with the $\sigma$ and the $x_p$, the electron current ($I_e$) and, thus, the discharge current ($I_d$) are found to decrease substantially with changes in the $B$-field’s profile. This is due to the consistent increase in the axial extent over which the magnetic field is relatively high.

Lastly, in addition to the reduction in the overall electrons’ transport (as represented by the $F_B$ term in the electrons’ azimuthal momentum equation), the profile of the $F_B$ tends to become more uniform with increases in either $\sigma$ or $x_p$.