Effect of auxiliary gas injection on the operation of a Hall current plasma accelerator

Optimization of magnetic and electric fields has been central concern for the design of a Hall current plasma accelerator since its inception decades ago. However, neutral flow dynamics in the discharge channel may have as much impact on the accelerator performance, operation stability and lifetime as the magnetic and electric fields due to its strong coupling with plasma properties. In this article, auxiliary gas injection is numerically studied for a low-power accelerator using a two-dimensional fully kinetic particle-in-cell code. Gas injection through the discharge channel sidewalls increases lifetime of the accelerator, but also degrades thrust performance suggesting that there is an optimum gas injection ratio. Although reduction in the maximum erosion rate is substantially lower than that predicted by a two-dimensional hybrid model for a high-power accelerator [14], extension of lifetime by approximately 20% appears to be possible with little impact (2%) on the thrust. The anode efficiency analysis supported by the simulated plasma properties clarifies that reduction in voltage utilization is the main cause of the observed alterations in the plasma properties and thrust performance deterioration.


Introduction
With the acceptance of reliability and cost benefits, there has been a paradigm shift towards the use of electric propulsion (EP) systems for spacecraft orbit control [1]. Hall current plasma accelerator, also referred to as Hall-effect thruster is a highly competitive form of EP systems due to its advantages such as large thrust-to-power ratio (>50 mN kW −1 ), high efficiency (>50%) at moderate specific impulse (>1500 s), design Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. simplicity (e.g. no grids) and outstanding in-space experience (flown since the 1970s) [2]. However, Hall current accelerator has a distinctive principle of operation involving intricate and non-linear physics leading to challenges in performance improvement and scaling-down for small satellites.
In a Hall current accelerator, the motion of free electrons is forced into spiral orbits by quasi-radial magnetic field (B) lines in an annular discharge channel with the radius of gyration (r g = m|u × B|/qB 2 ) much smaller than characteristic length of the discharge channel (L c ). The hot electrons with large pitch angles (u ⊥ > u/√r mirror ) bouncing forward and backwards between magnetic mirror points (r mirror = B max /B min > 1) ionizes most of the cold neutral atoms by electron-impact ionization collisions as the ionization mean free path (λ i = v n /n e ⟨σ i v e ⟩) is much less than unity (λ i ≪ L c ). Axial electric field applied between cathode and anode crossed with the radial magnetic field forces the electrons to drift in the E × B direction forming the azimuthal Hall current. The Hall current, which is much larger than the axial current (ω e τ e ≫ 1) interacts with the radial magnetic field to produce Lorentz force (F = ρ∇ϕ + J × B) accelerating quasi-neutral plasma mostly in the axial direction. Electric field equipotential lines or socalled 'thermalized potential' lines (ϕ ≈ ϕ 0 − T e ln[n e /n e,0 ]) approximately line up with magnetic field lines as electron mobility in the direction parallel to B is much larger than that in the vertical direction (µ || ≫ µ ⊥ ), and electron temperature (T e ) is constant along the magnetic field lines (∇ || T e ≈ 0) [3]. Magnetic field topography may therefore be expected to have direct influence on the accelerator performance and lifetime. Magnetic field, electric field (e.g. discharge voltage) and cathode coupling physics have been the main topics of research since the inception of Hall current accelerators in 1950s, whereas the neutral gas flow dynamics in the discharge channel have typically been underemphasized. Nevertheless, electron transport across the magnetic field is coupled with electric field and neutral flow dynamics.
Neutral flow dynamics may be described by density and thermal velocity distributions of the neutral atoms. Electron density (n e ) is a function of neutral density (n n ) and electron thermal velocity (v e ) as total collision frequency (v e = n n v th S total ) is proportional to neutral density and electron temperature. Fundamentally, a higher neutral gas density and a lower neutral thermal velocity in the discharge channel will improve ionization or propellant utilization efficiency (η m ) and overall performance of the accelerator. One can increase neutral density in the discharge channel by narrowing its geometry without changing propellant mass flow rate. Arhipov et al and Raitses et al reduced the discharge channel cross sectional area with a ceramic spacer, and found that ionization efficiency enhances due to increased local neutral number density and probability of ionization [4,5]. Neutral axial velocity is related to neutral atom residence time in the discharge channel (τ r = L c /v n ). The residence time should be much larger than the ionization time (τ r ≫ τ i = v n λ i ) so that the neutral atoms can be ionized, and accelerated by the Lorentz force before escaping the discharge channel. An increase in the neutral atom temperature at a given operating condition may lead to a direct reduction in the propellant utilization efficiency. Wilbur et al cooled the discharge channel wall by liquid nitrogen, and showed that reduction in the neutral atom temperature induces improvements in extracted ion beam current [6]. Book and Walker observed that an actively cooled anode achieves higher thrust and anode efficiency compared to a passively cooled anode [7]. The fact that temperature of the anode wall and discharge channel walls affect the neutral atom velocity (v th = √8k B T/mπ) indicates a relationship between the neutral atom velocity and the accelerator performance. However, this relationship may not be linear because the neutral atom temperature might differ from the mean wall temperature due to the perturbing effects induced by atomic collisions with electrons and ions [6] even when neutral atoms are in thermal equilibrium with the walls.
In addition to neutral atom density and temperature, flight paths of the neutral atoms may be altered through injection schemes in order to investigate effect of neutral flow dynamics. Miyasaka et al reduced the neutral atom velocity (v th = √8k B T/mπ) at the discharge channel centreline by increasing the anode orifice diameter, and observed decreased discharge current oscillations [8]. Vial et al demonstrated that changing radial position of injection holes on the anode surface may affect plume behaviour to reduce divergence angle and discharge channel wall erosion rate [9]. Recently, Ding et al studied the effect of a tangential injection scheme similar to one used in vortex reactors [10], and found that propellant utilization and thruster performance increases remarkably at low discharge voltage operations [11]. A more detailed literature review regarding neutral flow dynamics may be found elsewhere [12].
Despite the low attention it attracted, neutral flow dynamics plays a prominent role on Hall current accelerator operation, and remains a topic of interest for researchers. One promising method is auxiliary gas injection. Goebel et al found that auxiliary gas injection may quench the plasma potential oscillations in hollow cathode discharge plume reducing highly energetic ion production responsible for keeper electrode erosion [13]. Garrigues et al showed with a twodimensional hybrid code that auxiliary gas injection may flatten radial plasma potential profile of a high-power Hall thruster discharge decreasing thruster wall erosion rate, which determines thruster lifetime [14]. The aim of this article is to revisit the auxiliary gas injection method on a low power Hall current accelerator, and examine it in more detail using a 2D3V (two-dimensional in space and three-dimensional in velocity) fully kinetic code. Fully kinetic codes are computationally expensive and time-consuming because electrons and ions are treated as discrete particles and velocity distribution function is computed self consistently; however, they can capture non-quasi-neutral and non-equilibrium plasma physics unlike the fluid and hybrid codes. Since the repetitive experiments with different accelerator configurations are laborious and costly, fully kinetic simulation can be an effective tool.
A 500 W class laboratory model Hall current accelerator of the University of Tokyo named UT-62 [15] is used as the testbed. Accelerator thrust performance, discharge channel wall erosion and plasma properties such as plasma number density and electron mobility are quantified for a fixed discharge channel geometry, magnetic flux density, channel material, discharge voltage (V d ) and anode mass flow rate (m). Thrust efficiency ( is decomposed into the product of energy efficiency (η E ), propellant efficiency (Φ P ) and beam efficiency (Ψ B ) according to the efficiency decomposition scheme proposed by Brown et al [16] to analyse thrust generation and loss mechanisms.
The remainder of this article is organized as follows. Section 2 gives an overview of the numerical model and simulation conditions; and section 3 presents the simulation results and discusses on the findings.

Description of the model
The two-dimensional (axial and radial) axisymmetric fully kinetic particle-in-cell (full-PIC) code used in this article is the same as the one described in previous publications [17,18]. The full-PIC code relies on an artificial mass ratio model to speed up computations and a semi-implicit method to derive electric field from Poisson's equation permitting stable simulations for larger time steps than electron plasma frequency (ω p ∆t ≫ 1) and larger grid spacing than Debye length (∆x/λ D ≫ 1). The computational time-step (∆t = 5 × 10 −10 s) was chosen such that the plasma oscillation and electron cyclotron motion were captured. The grid spacing and macroparticle size were decided according to the parametric study conducted previously for this specific accelerator [18]. A grid spacing of 0.2 mm was sufficient to accurately capture the sheath potential drop, even though it is larger than the Debye length (λ D = ∼0.05 mm), and a macroparticle size of 1 × 10 7 was a good trade-off between the numerical error from so-called discrete particle noise and the computing cost.
The Lorentz force equation of motion of charged particles is solved using the 4th-order Runge-Kutta method. Inter-particle collisions are modelled by the direct simulation Monte Carlo method. Neutral atoms, singly charged ions, doubly charged ions, primary electrons and electron induced secondary electrons; and corresponding ionization, excitation and elastic scattering collisions are included in the simulation. Triply charged ions and charge exchange collisions are neglected due to small collision cross-section [19].
Electrons are injected through plume boundary cells according to a Maxwellian distribution with an assumed temperature of 2 eV to start and sustain the discharge. Note that the full-PIC code does not self-consistently simulate the accelerator-cathode coupling process [20] as the beam neutralization process occurs outside the computation domain. Bohm diffusion is implemented in terms of virtual scattering collision between neutrals and electrons, v Bohm = 1/16(eB/m e ) [21], in order to compensate for the low electron mobility across the magnetic field due to the lack of azimuthal dimension.
The surface charge accumulation on the dielectric walls due to collisions with the charged particles is included in the electric field calculation. Xenon neutral atoms are injected from the anode with a thermal energy of 500 K, and temperature of the discharge channel walls are set at 700 K. Neutral atoms and ions are thermalized on the discharge channel wall surface with thermal accommodation coefficient (TAC), which defines the percentage of heat transfer between heavy particles (ions and neutrals) and solid wall surfaces. When ions and neutrals atoms strike the anode wall and discharge channel walls, they are reflected with an energy determined by a TAC of 0.9. The erosion rates of the boron nitride (BN) discharge channel walls are calculated according to a model developed by Yamamura and Tawara [22] assuming a threshold energy of 35 eV for ion sputtering.
The code was parallelized using the message passing interface and portable, extensible toolkit for scientific computation,

Simulation conditions
The UT-62 accelerator has an outer discharge channel diameter of 62 mm, an inner discharge channel diameter of 48 mm and a discharge channel length of 21 mm. It employs a single electromagnetic coil consisting of 180 turns of 1 mm diameter copper wire. The magnetic field was calculated using the software Finite Element Method Magnetics. Magnetic field strength varies between 0.005 and 0.03 Tesla along the discharge channel centreline (for a 2 A coil current), and its strength increases from the anode towards the thruster exit. The magnetic field induced by the plasma is considered to be negligible.
The simulated operating conditions of the accelerator are discharge voltage of 300 V, xenon mass flow rate of 1.36 mg s −1 and electromagnetic coil current of 2 A corresponding to a maximum magnetic flux density of 0.028 Tesla at the discharge channel centreline. Simulations were run for 0.4 million steps for each case, and the results were timeaveraged over 0.2 ms wall clock time every 100 simulations steps, which is approximately equal to the time a neutral particle takes to move from the anode to the plume boundary i.e. exit the computation domain. The convergence of the simulation was judged by the time-history of the total number of simulation particles.

Accelerator thrust performance
The UT-62 accelerator has a measured thrust (T), discharge current (I d ), thrust efficiency (η T ) and specific impulse (I sp = T/mg) of 16.4 mN, 1.0 A, %33 and 1226 s respectively for a measured vacuum chamber background pressure of ∼3 × 10 −5 Torr against a 1.36 mg s −1 anode and 0.27 mg s −1 hollow cathode xenon mass flow rates. The full-PIC code's prediction of the accelerator performance (16.0 mN, 1.1 A, %29 and 1201 s) matches very well with the measured values with a maximum discrepancy less 3% for the thrust and 10% for the discharge current. Note that the numerical discretization errors as well as the uncertainties in model physical parameters such as TAC and wall temperature may affect prediction accuracy of the fully kinetic codes [18]. Effects of vacuum chamber background pressure [23] and cathode coupling voltage are not included in the simulation, which may account for the small discrepancy between measurements and predictions.
The 2D plasma density plot (see figure 4 nom) shows an electron concentration in the vicinity of outer discharge channel wall because radial magnetic field strength along the inner wall surface is ∼1.5× higher than the outer wall surface. The plasma density reaches its maximum near the outer wall at an offset distance of 2 mm from the discharge channel centre line (r = 29.5 mm). Figure 1 presents normalized distributions of the main plasma properties at this location for a quantitative comparison instead of the thruster centreline. The peak electron temperature nearly overlaps that of the radial magnetic field. Plasma potential profile suggests a sharp increase in axial electric field and local maximum in the ion velocity near the discharge channel exit. These features are representative of the characteristics of Hall current plasma accelerators [24,25].
There is an interesting feature of plasma density in the near anode region seen in figure 1. The plasma density goes up (∼0.8 × 1018 m -3 ), and then goes down (∼0.4 × 1018 m -3 ), before going monotonically up. This feature is in fact difficult to confirm with experiments because plasma density measurements with Langmuir probe usually have an uncertainty on the order of ±50%. One possible explanation is that there is a slight plasma potential drop accompanied by electron heating in the near anode, which may improve ionization rate there. Ionization collision distribution plot (not given here), which roughly matches with the plasma density distribution, support this.

Auxiliary gas injection
In a conventional Hall current accelerator, the anode electrode also serves to inject neutral atoms to the discharge channel through an annular slot or azimuthally arranged orifices. In the auxiliary gas injection method, additional gas flow is injected from discharge channel inner and outer walls as illustrated in figure 2. Auxiliary gas injection into the discharge channel from the inner and outer walls was varied in the range of 0%-100% with 20% increments, while keeping the total xenon mass flow rate constant at 1.36 mg s −1 . The amount of gas flow rate shared equally between the inner and outer walls and the gas flow is injected from a circumferential slot of 2 mm width as this is easy-to-implement in practice. A smaller distributed area of injection is preferred as it boosts ionization efficiency due to larger neutral density.
Auxiliary gas was injected at a fixed location of z/L c = 0.8 where wall incident ion mean energy starts to cross the sputtering threshold of BN as shown in figure 3. Although ionization and acceleration regions overlap throughout the discharge channel, acceleration becomes dominant at this location because electron temperature begins to rise significantly, while plasma density decreases (see the dashed line in figure 1). The auxiliary gas injection at an upstream location (z/L c < 0.8) is not expected to have noteworthy influence.
The auxiliary gas injection at an downstream location (z/L c > 0.8) may have large influence on wall erosion rates as the exit plane (z/L c = 1.0) is the region of the maximum erosion. However, this might also affect adversely the thrust performance due to increased unionized neutral propellant current (I n ). Therefore, the location of z/L c > 0.8 represents a good compromise between thrust performance and wall erosion rate in line with Garrigues et al [14]. The C3 case, which has a thrust value of 15.0 mN and maximum erosion rate of 14.4 mm kh −1 , was also simulated with injection locations of z = 0.7L c and 0.9L c . It was found that thrust reduces to 15.4 and 14.6 mN respectively from its nominal value of 16.0 mN, and maximum erosion rate reduces to 17.4 and 11.0 mm kh −1 from its nominal value of 21.5 mm kh −1 . Reduction in wall erosion rates and thruster performance is larger for a downstream location verifying our rationale for the value of the auxiliary neutral injection location. Table 2 summarizes the simulation cases investigated in this study.   Table 3 gives thruster performance results from the simulation of the investigated cases and anode efficiency analysis. It can be seen that thrust (F), specific impulse (I sp ) and anode or thrust efficiency (η t ) decreases monotonically as the percentage of the neutral gas flow injected from the discharge channel walls increases; whereas discharge current oscillation amplitude (∆ = RMS/I d ) increases. Discharge current (I d ), singly charged ion beam current (I b1 ) and doubly charged ion beam current (I b2 ) also decrease, but they remain the same until at least half of the total gas flow is injected from the walls (m anode ⩾ m wall , where m wall = m inner _ wall + m outer _ wall ).

Thruster performance
Neutral propellant current, defined as the number of neutral particles leaving the thruster exit without getting ionized, starts to increase when the auxiliary gas injection is more than half of the total gas flow (m wall > m anode ). The ionization and acceleration processes are interrelated, it is therefore useful to separate thrust efficiency into three partial efficiencies, namely energy efficiency, propellant efficiency and beam efficiency in order to better understand the effect of the auxiliary gas flow. The energy efficiency can be further decomposed into voltage utilization efficiency and current utilization efficiency as follows: η E = η V η C . With increasing auxiliary gas injection ratio (m wall /m anode ), voltage utilization efficiency, also referred to as the acceleration efficiency, decreases substantially implying that the potential drop experienced by ions before leaving the discharge channel, i.e. electric field inside the channel that accelerates ions, becomes lower. Charge utilization efficiency increases because decrease in discharge current is sharper than of the thrust-generating ion beam current. However, decrease in energy deposited on the propellant, i.e. ion jet kinetic energy, is higher than decrease in input discharge power leading to decreased overall energy efficiency. Propellant efficiency, which is an indicator of the ionization level, also decreases, but for wall mass flow rate values exceeding anode mass flow rate. Interestingly, beam efficiency tends to increase as more propellant fed near the thruster exit plane.
For the cases C1, C2 and C3 in which m anode ⩾ m wall , the observed loss in thrust efficiency could be attributed to voltage utilization since the propellant efficiency is not affected by gas injection on the walls while current utilization and beam efficiency improve contrary to the voltage utilization. The theoretical thrust is proportional to ∼cosθI b √V b , where cosθ is beam divergence, I b is beam current and V b is beam voltage, which is the potential across which the ions are accelerated. It can thus be suggested that gradients in the axial profiles of plasma potential near the thruster exit plane decrease. As for the cases C4-C6 in which m wall > m anode , propellant utilization also contributes to the thrust efficiency loss.
In order to utilize maximum amount of discharge voltage for the acceleration of ions, it is important to sustain the local plasma potential in the upstream region of the discharge channel close to the anode electrode potential. As m anode decreases and m wall increases, plasma conductivity will be lower in the upstream region due to lower electron-neutral collision probability moving ionization region downstream, which can explain the decrease in voltage utilization.
Further understanding on the effect of auxiliary gas injection is obtained from the plasma properties presented in the next section. Figure 4 shows simulated plasma density distribution with ionstream traces. As it can be seen, plasma density is not concentrated in the vicinity of outer discharge channel. Radial magnetic field strength along the outer dielectric wall surface is lower compared to the inner wall surface because the UT-62 accelerator employs a single inner electromagnetic coil. This may lead to an increase in electron diffusion and non-uniform plasma density. In the nominal case, dense plasma arises in the near anode region where the neutral density is the highest, and extends across the discharge channel. As the injected gas flow rate from the discharge channel walls increases, the plasma column become shorter and smaller, and an ionization region  is shifted downstream. Plasma density in the discharge channel including the near anode region decreases by ∼20% from a maximum of 1.1 × 10 18 m -3 -0.9 × 10 18 m -3 . Figure 5 shows simulated plasma potential distribution with applied magnetic field lines. The equipotential lines follow the magnetic field lines more closely inside the discharge channel as the ratio of plasma potential to electron temperature is greater [26]. Approximately two third of the plasma potential drops outside the discharge channel in the nominal case. The anode and near-anode quasi-neutral plasma (z < 10 mm) are both at the same potential, indicating electron axial velocity is on the order of the electron thermal velocity in this region. However, an electron-attracting anode sheath 'positive anode fall' is established with increasing m wall /m anode and more than half of the discharge voltage drops inside the discharge channel in the cases of C4-C6. The need for enhanced ionization to increase the electron flux toward the anode might explain the formation of a positive fall [27]. Figure 6 shows simulated electron temperature distribution. Hot electrons extend further upstream on the inner wall surface. A possible explanation may be that electron-wall collision frequency is lower there due to higher magnetic field strength reducing electron energy loss. Electron temperature near the thruster exit plane decreases considerably with increasing m wall /m anode . In the near anode region, electron temperature increases due to the formation of the positive fall with the maximum electron temperatures of 4 and 12 eV for the nominal case and C6 respectively.

Plasma properties
The illustrated plasma properties support the findings from the anode efficiency analysis. The decrease in propellant efficiency with increasing m wall /m anode explains the lower plasma density. The increase in beam efficiency may be understood from the ion-stream traces. In the nominal case the first two ion-stream lines from the top diverges sharply in the radial direction outside the discharge channel than that of the case C6. The decrease in the voltage utilization efficiency may be seen in the plasma potential distribution, which may account for the observed decrease in the electron temperature. Figure 7 shows electron mobility perpendicular to magnetic field lines at an offset distance of 2 mm from the discharge channel centreline near the outer wall (r = 29.5 mm). It was calculated from the simulated plasma properties according to Koo and Boyd formulation [28] as follows: µ ⊥ = j e ⊥ / (E ⊥ + ∇ ⊥ p e /en e ) where j e ⊥ is electron current density  across the magnetic field lines and p e is electron thermal pressure (p e = kT e n e ).
Neutral atoms play a significant role in increasing the crossfield electron mobility through electron-neutral collisions, which is proportional to the neutral density. Thus, increasing the mass flow rate from the walls increases the electron mobility yielding to a peak occurring near the thruster exit plane. For the same reason, the electron mobility inside the discharge channel decreases with increasing m wall /m anode giving cause for the lower plasma conductivity in the upstream region mentioned in the anode efficiency analysis. However, in the near anode region (z/L c < 0.5), the electron mobility may be higher despite the decreasing mass flow rate, which may be related to the formation of the electron-attracting anode sheath. Ionization oscillations can be excited due to reduced electron transport and insufficient electron temperature (T e < T e,min , T e,min = 12.4 eV for xenon) [29]. This might explain the observed highly oscillatory behaviour of the discharge current (∆ > 0.2) in the low frequency range .
The choice of anomalous electron transport model remains to be open question as there are various discussions about which and how to employ the model (e.g. axially-varying or uniform) [30]. The fully kinetic code employs the simplest and acknowledged Bohm diffusion model, which includes a spatially uniform coefficient, namely 1/16B, to account for the enhanced mobility induced by the azimuthal electrostatic instabilities. The anomalous electron transport model may affect the plasma property profiles. For example, the use of a constant Bohm coefficient may result in the electric potential drop and electron temperature to peak outside the discharge channel exit as we reported in a previous sensitivity analysis of the fully kinetic code [18]. Experimental plasma properties for the UT-62 thruster were not available for a comparison with the simulated ones. However, measurement accuracy of the global quantities like thrust and discharge current is much higher than measurement of plasma properties due to intrinsic difficulties in probing flowing supersonic and magnetically confined plasmas. Additional simulations were run for the nom, C2, C4 and C6 cases with the Bohm diffusion model turned off to clarify the effect of the choice of electron transport model, and it was found that the overall conclusions do not qualitatively change. Further investigation on the Bohm diffusion model itself is beyond the scope of this work.

Discharge channel wall erosion
Flight-qualified low-power Hall current accelerators have a lifetime limited to a few 1000 h. Increasing lifetime of low power Hall current accelerators (<1000 W) such as for deep space missions remains an active area of research [31,32]. The discharge channel wall surfaces erode by high-energy ion bombardment, rate of which usually determines the accelerator lifetime. Erosion rate is a challenging parameter to obtain from experiments. Figure 8 shows the simulated wall incident ion flux distribution for the inner and outer walls. Ion flux into the outer wall greatly exceeds that of the inner wall due to non-uniformly distributed radial magnetic field.
Peaks of the ion flux occur at four distinct regions: near anode region, mid-channel, gas injection location and thruster exit plane, which also correspond to the peaks of the plasma density distribution. With increasing m wall /m anode , ion flux into the walls decreases except for the case of C2. High plasma potential region recedes toward the anode, and the equipotential lines, which coincide with the shape of the magnetic field lines, are slanted towards the outer channel as shown in figure 5. This might lead to acceleration of more ions from the bulk plasma toward the outer wall in the cases of C2, C4 and C6. Nevertheless, anode mass flow rate and, as a result, neutral density in the middle of the channel is much higher for C2, and this explains why the ion flux into the outer wall is larger than the others. Figure 9 shows the simulated wall incident ion mean energy distribution for the inner and outer walls. Ion mean energy decreases at the gas injection location as the average  peak electron temperature decreases with more frequent electron-neutral collisions. Ion mean energies are greatest at the thruster exit where the axial electric field also peaks. With increasing m wall /m anode , ion mean energy at the upstream of the discharge channel increases, indicating increase in the radial electric fields causing radial expansion of ion production. Figure 10 shows the simulated discharge channel wall erosion rates for the inner and outer walls calculated according to the wall-incident ion flux and ion mean energy. In the nominal case, the channel upstream is essentially erosion-free for the nominal case since the incident ion mean energy is lower than the sputtering threshold energy of the ceramic walls. The maximum erosion happens at the discharge channel exit where the ion flux and energy are also the highest with rates of 5.6 and 21.5 mm kh −1 respectively for the inner and outer walls. Ions bombarding the inner wall surface has lower mean energy than the outer wall surface, but this is not enough to explain the large variation in the observed erosion rates. Both ion flux and energy should contribute to the channel wall erosion rate.
The wall erosion rate increases exponentially along the discharge channel, reaching a maximum at the channel exit for all of the cases. As expected from the ion flux and ion mean energy plots, the maximum erosion rates of both inner and outer wall decrease with increasing m wall /m anode . Yet, wall erosion increases at the channel upstream as a result of the increased ion mean energy.
Injecting a percentage of the neutral atoms from the discharge channel walls alleviates erosion, and thus increases the accelerator lifetime, but aggravates the thrust performance concurs with the previous two-dimensional hybrid model simulation results of a SPT-100 Hall accelerator reported by Garrigues et al [14]. Nevertheless, the maximum erosion rate is reduced merely 33% for 50% of the xenon injection at z/L c = 0.8, not by a factor of 10 for the same injection percentage and location despite the fact that the maximum erosion rate of the UT-62 accelerator is much higher with respect to the SPT-100 accelerator in the nominal case (5.2 mm kh −1 ). This may be because plasma-wall sheath is not self-consistently modelled in the hybrid model unlike the full-PIC code, leading to significant overestimation of the mean energy of the incoming ions on the walls.
Lower wall erosion rate enables longer lifetime, but lower thrust translates into longer acceleration times. This suggests there should be a trade-off between thruster lifetime extension and thrust decrease. In the case of C1, which represents the lower bound in terms of auxiliary gas injection percentage, the decrease in thrust is less than 1%, while the maximum erosion rate decreases as large as 7%. In the case of C6, which represents the upper bound, the decrease in thrust is 59%, while the maximum erosion rate decreases as large as 90%. The optimum neutral injection percentage maximizing the erosion reduction and minimizing the thrust performance loss seems to be 40% i.e. in the case of C2, where the thrust and the maximum erosion rate decrease by 2% and 17%, respectively. The state-of-the-art low power Hall accelerators have a lifetime of a few 1000 h (<4000) at most, a wall erosion reduction on the order of 20% may extend the accelerator up to 800 h with the performance is almost unaffected. Therefore, a low power Hall thruster with auxiliary injection is relevant for a practical spacecraft purpose, and this method requires only small modifications to the thruster structure.

Conclusion
Additional gas flow injection through the discharge channel sidewalls of a low power Hall current accelerator was studied by keeping total mass flow rate constant using a 2D full-PIC code. Increasing the amount of injected neutral atoms from the discharge channel walls improves the accelerator lifetime, but at the expense of thrust performance. Voltage utilization appears to be the main mechanism responsible for the thrust performance loss. With increasing wall mass flow rate, the plasma potential at the upstream discharge channel cannot be sustained close to the anode potential. A shorter discharge channel in which the anode electrode is closer to the channel exit might improve the voltage utilization.
Effective electron mobilities calculated at the discharge channel centreline provide evidence that plasma conductivity reduction due to decay of the electron-neutral collision frequency may explain the decreasing voltage utilization. Discharge current becomes highly oscillatory possibly because of the decreasing voltage utilization, and thus electron temperature or equivalently ionization rate coefficient. Maximum electron temperature decreases sharply compared to the plasma density. Anode sheath fall becomes positive to compensate for the decreasing ionization rate coefficient, and increase the electron flux toward the anode.
Maximum erosion rate of the sidewalls at the discharge channel exit plane decreases, but wall erosion increases at the channel upstream owing to increasing radial electric fields and radial expansion of ion production. Even though improvement of the accelerator lifetime by the auxiliary gas injection method is not as large as that reported in [14], extension of lifetime by ∼20% appears to be possible with little impact (2%) on the thrust. There are also unexplored operating regimes that might lead to better performance. For example, the accelerator may use two separate propellants supplying the one with higher ionization energy from the discharge channel side walls. Due to their larger surface-to-volume ratio, this method may work better with very low power Hall accelerators (<100 W) for which there is a surge of commercial interest. Finally, mechanisms causing erosion reduction may be a result of the interplay between magnetic field and local neutral density, so one may pay attention to magnetic field distribution.

Data availability statement
The data that support the findings of this study are available upon reasonable request from the authors.