Numerical investigation of reversed gas feed configurations for Hall thrusters

by   Stefano Boccelli, et al.

A reversed gas feed configuration for Hall thrusters is proposed and studied numerically. As an alternative to standard direct injection, we investigate the effect of injecting the propellant near the channel exit and directing it backwards, towards the anode. The resulting neutral density and average velocity fields are studied with the Direct Simulation Monte Carlo method for both cold and warm anode conditions. The residence time of neutral particles inside the channel and in the ionization region is computed using a test-particle Monte Carlo method. The computations indicate that the reversed injection allows for increasing the mass utilization efficiency of a standard feed configuration from 2-5 efficiency.



There are no comments yet.


page 6

page 10

page 18


A higher order moment preserving reduction scheme for the Stochastic Weighted Particle Method

The Stochastic Weighted Particle Method (SWPM) is a Monte Carlo techniqu...

Discovering Quantum Phase Transitions with Fermionic Neural Networks

Deep neural networks have been extremely successful as highly accurate w...

A meshfree arbitrary Lagrangian-Eulerian method for the BGK model of the Boltzmann equation with moving boundaries

In this paper we present a novel technique for the simulation of moving ...

Determinant-free fermionic wave function using feed-forward neural networks

We propose a general framework for finding the ground state of many-body...

Accelerating the estimation of energetic particle confinement statistics in stellarators using multifidelity Monte Carlo

In the design of stellarators, energetic particle confinement is a criti...

A fast GPU Monte Carlo Radiative Heat Transfer Implementation for Coupling with Direct Numerical Simulation

We implemented a fast Reciprocal Monte Carlo algorithm, to accurately so...

Classifying Multi-Gas Spectrums using Monte Carlo KNN and Multi-Resolution CNN

A Monte Carlo k-nearest neighbours (KNN) and a multi-resolution convolut...
This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

1 Introduction

Hall effect thrusters (HET), also known as stationary plasma thrusters (SPT) [1, 2], are effective devices for attitude control and maneuvering of spacecraft, and are envisaged as propulsion systems for future manned missions inside the solar system. The benefit of employing electric propulsion engines, as opposed to chemical thrusters, lies in the high specific impulse, that exceeds 1000 s and can reach 5000 s or more for high specific impulse HETs [3, 4]. Among the figures of merit, the mass utilization efficiency (or propellant efficiency) plays a crucial role [5]:


with the mass flux of ions ejected by the thruster and the injected flow of propellant. If a hollow cathode is used to neutralize the plume, then includes an additional cathode loss contribution. In this work, we neglect this for simplicity, assuming that no propellant is lost through the cathode. This would actually be the case if a filament cathode is employed, as frequently done in experiments [6].

The mass utilization efficiency varies widely with the thruster size and operating conditions [7, 8, 9], and can range from below 0.5 for micro-thrusters up to above 0.9. The gas injector design plays a crucial role in the mass utilization efficiency, as it influences the residence time of neutral particles inside the chamber, and thus the possibility that neutrals are ionized by electron impact collisions. A large number of gas feed configurations have been tested throughout the years [10], including injection ports located at the middle of the channel [11] and various anode configurations.

With possible exceptions inside the gas injector, the neutral gas in the Hall thruster channel is characterized elsewhere by strong degrees of rarefaction, as expressed by the Knudsen number , where is the particle mean free path and the characteristic size of the thruster channel. Besides the channel size, the value of for Hall thrusters also strongly depends on the mass flow rate. In this work we consider the SPT-100 thruster [12], where denotes the outer diameter in mm and with an inner diameter of 70 mm. For such thruster, the Knudsen number based on the channel width lies roughly in the range of 0.1 to 1, as will be discussed in Appendix A. The neutral flow is therefore transitional and shall be analyzed by use of the kinetic theory of gases [13], employing for example the Direct Simulation Monte Carlo (DSMC) method [14]. Such method can treat equally accurately both collisional and collisionless regimes. Other collisionless approaches such as view-factor methods (see for example [15]) could prove useful for low mass flow rates or smaller devices such as the SPT-50 or SPT-20 thrusters [6], characterized by a larger surface-to-volume ratio.

Aim and structure of this work

In this work, we investigate via numerical simulation the possibility to employ a reversed gas feed configuration: the propellant is injected near the outer end of the channel and directed inwards, towards the anode (see Fig. 1). Such strategy clearly presents a number of additional technical challenges in the design of the gas feed system, if compared to a direct injection. The present study does not consider such technical difficulties, but merely investigates the resulting neutral gas dynamics, once injected into the chamber. For simplicity, all plasma effects will be neglected, except for the increased temperature of the anode [16, 17]. Electron-neutral and ion-neutral collisions will be neglected, although these are known to affect the neutral temperature and momentum. Despite the heavy simplifications, results will provide useful insights on the neutral behavior for the considered conditions.

This work is structured as follows. At first, a heavily simplified analysis of the problem is discussed in Section 2

, and the ionization probability is written in terms of the neutral residence time in the ionization region. In Section 

3 the neutral flow field is investigated more in detail for direct and reversed injection configurations using DSMC simulations. In Section 4 a test-particle Monte Carlo algorithm is implemented and parallelized on a GPU. The algorithm follows single particles in the thruster channel, moving in a pre-computed DSMC solution (for collisional conditions) or in a vacuum (for collisionless cases). Such algorithm is employed to study the average residence time of particles in different configurations in Section 5. Distributions and average values of residence times are discussed. Finally, the obtained residence times are employed to study the increased mass utilization efficiency for the reversed injection against a direct injection efficiency. Conclusions are drawn in Section 6.

2 Direct vs reversed injection: simplified analysis

Figure 1: Single particle trajectory with specular wall scattering. Left: direct injection from the anode. Right: reversed injection. Ionization region is labeled “iz”.

A strongly simplified view of the effect of reversed injection can be obtained via a preliminary model, where we assume (i) monoenergetic particles (no thermal velocity), (ii) free molecular conditions and (iii) specular wall collisions. More representative conditions will be addressed in the following sections. The problem is depicted in Fig. 1. Assuming monoenergetic neutrals implies that all particles have the same initial axial velocity , and therefore all neutrals injected inwards will enter the channel. In free molecular flow conditions, gas phase collisions are negligible: every neutral moves independently from the others, and its residence time inside the thruster channel depends on its velocity and on the interaction with the walls. If we additionally assume a specular reflection wall-collision model (zero energy and momentum accommodation), the time of permanence inside the ionization region for neutral atoms injected from the anode (Fig. 1-(a)) simply reads , where is the length of the ionization region. On the other hand, by injecting a neutral from the open end of the channel (with the same value of the axial velocity, but opposite sign) such particle would cross the ionization region twice, and the residence time would double: . Since electron-neutral ionizing collisions can be assumed to be uncorrelated events, the probability that a given neutral particle ionizes during the time interval is


where is the electron-neutral ionization frequency, that depends on the local electron density, temperature and in general on their distribution function [18]. For the goals of this work, it is not necessary to specify further. The ionization probability in the reversed injection case can thus be obtained from the direct injection case as


In this simplified model, we have , such that for all conditions, with the largest improvements for smaller ionization probabilities.

Clearly, this model is heavily simplified. In a more realistic scenario, the non-zero temperature of neutrals is such that the residence times are characterized by a distribution that depends on the geometry. Each different particle has a different ionization probability, slow particles being naturally more likely to ionize since they spend more time inside the ionization region. Moreover, due to the finite temperature, a fraction of particles injected near the channel end will be lost in the plume prematurely. Such particles will constitute in a loss in terms of residence time.

The effect of wall collisions also strongly affects the dynamics of the problem, according to the value of the accommodation coefficient [13]. This may significantly depart from 0 (specular reflection) and approach unity (diffuse collisions), limiting the mobility along the channel. Finally, the presence of gas-phase collisions would have the effect of pushing atoms towards the exit. The relative importance of these phenomena depends on the Knudsen number, and thus on the thruster size and mass flow rate. All mentioned effects will be considered in the remaining of this work for two direct injection cases and for various reversed injection configurations.

3 DSMC simulations

In this section, we perform DSMC simulations of direct and reversed injection strategies, as to obtain density, velocity and temperature fields of the neutral gas. These simulations are complementary to the simplified analysis of Section 2, and directly include the effect of finite thermal velocity, diffuse wall collisions and gas-phase collisions.

The simulations are performed with the SPARTA DSMC code [19] on a 2D computational domain with axial symmetry, and represents a slice in the radial-axial plane , as shown in Fig. 2. The considered geometry is characteristic of the SPT-100 thruster [12], and the propellant mass flow rate is fixed at mg/s, injected either uniformly from the whole anode (case A1), or from a small slit in the anode (case A2). In cases “B” the reversed injection is considered, with angle with respect to the axial direction ( is a radial injection), and a slit centered at axial position . Cases B30_1, B60_1 and B90_1 consider different injection angles of and respectively, for a given slit position located close to the exit plane. In cases B30_2, B30_3 and B30_4, we keep the angle and progressively move the slit towards the anode.

= 22.5 mm 18.5 15.5 12.5
B30_1 B30_2 B30_3 B30_4
B60_1 - - -
B90_1 - - -
Table 1: Test cases for the reversed injection configuration.

Xenon gas is considered. All cases are run for two different conditions of wall temperature and anode temperature :

  • Cold anode: K;

  • Warm anode: K, K.

Figure 2: Computational domain and boundary conditions for the SPT-100 DSMC simulations. Axial-radial plane.

Actual wall and anode temperatures for real thrusters strongly depend on the operating conditions and require some hours of operation before reaching a steady state [20, 21]. In cases A1 for cold and warm anode, atoms are injected from radial position to mm, assuming a Maxwellian distribution with zero average velocity and the anode temperature. The inlet flux is thus completely thermal, and the number density is set as to result in the required mass flux. The slit size is 1 mm and for the B cases is located on the inner wall. The boundary condition to reproduce a physically accurate orifice expansion requires particular care (see for example [22]) and ultimately depends on the detailed structure of the gas feed. For simplicity, we assume sonic injection, at the anode temperature in case A2 and at wall temperature for cases B. Particles injection is therefore performed by sampling a Maxwellian distribution with sonic average velocity, equal to m/s for a temperature of K and m/s for injection temperature of K. For case A2 the average injection velocity is axial, while for cases B its components form an angle with respect to the side walls. The number density is finally chosen as to give the required mass flux.

Figure 3: DSMC simulations. Left and Center: density and velocity at the channel centerline. A1, A2: direct injection. B: reversed injection. Right column: number density contours for cases A1, B30_1 and B30_4, cold anode.

The DSMC cell dimension is chosen as to resolve the local mean free path, and a discretization of 280280 cells in the whole domain proved adequate. The time step is set to 1 s as to properly resolve the mean free time. The number of simulated particles in steady state is in the order of , providing sufficient resolution inside the channel and in the first part of the plume. A variable soft sphere (VSS) collision model is employed [14], where the cross section is , and the diameter depends on the relative velocity ,


with relative mass . The deflection angle reads


where is the impact parameter. In this work, we use a reference diameter at the reference temperature , viscosity index and scattering parameter (see [14]). Each simulation takes roughly 15 minutes on a laptop (Intel i5 processor, single core computation). The resulting density and velocity fields are shown in Fig. 3 for the cases A1, A2, and for the angled injection cases of Table 1. In all cases, the effect of the anode temperature is seen to decrease the density and increase the average velocity, as expected. The reversed injection strategy appears to have a strong impact on the averaged fields. All three injection angles level out the density inside the chamber and bring the average velocity to zero, with the injection at providing the highest density. Referring to the center column in Fig. 3, a low (zero in this case) average velocity can be obtained by the reversed injection approach, suggesting that the average residence time of neutrals can be increased. In terms of flux balance, the reversed injection configuration mimics the flow in a closed end channel, where the global mass flux is zero. Still, individual particles reach the anode via thermal motion, eventually being reflected back and leaving the channel. The rarefied conditions are such that the average fields only show a marginal picture of the problem, and a more in-depth analysis of residence times should be performed (see sections 4 and 5).

Observing the density contour in the right column of Fig. 3, it appears that the actual direction of injection is much closer to . This can be explained by considering that the velocity field is an average of (1) the population of injected particles at and (2) the population of scattered particles that are leaving the domain. For injections at very small angles, the second effect can be expected to increase, as the fraction of particles colliding with the inner wall increases.

The effect of shifting the slit towards the anode is shown in Fig. 3-Bottom row. As could be expected, the density and velocity profiles approach the uniform injection case A1, as the slit position approaches the anode. For such cases, only the warm anode case is shown. Results for the cold anode case are qualitatively analogous and will not be shown here.

3.1 Effect of ionization and ion-neutral interactions on the DSMC fields

As mentioned, ionization has been neglected in this study. For this reason, we can expect the density profiles of Fig. 3 to be somehow over-predicted, especially after neutrals have crossed the ionization region. The accuracy of our simulations ultimately depends on the operating conditions, in terms of magnetic field and applied potential. For example, Koo and Boyd [23] report the neutral density profiles for an SPT-100 thruster in presence of ionization, for a mass flow rate comparable to this work. In their simulations, the neutrals density appears to be strongly affected by ionization, differing by up to a factor 10 from our neutral results, in some parts of the domain.

It should be noted that, if one considers completely collisionless operating conditions, these potentially large discrepancies would have no effect in terms of neutral residence times, since collisionless conditions make an accurate prediction of the baseline fields completely unnecessary. In the situations considered in this work, collisionality appears not to be the predominant factor, as also discussed in Appendix A. Although it does play a significant role (as will be further discussed in Section 5.2), we could expect our results to be anyway indicative of the general behavior. In any case, more complete simulations are suggested as a future work to investigate this to a full extent.

Secondly, ion-neutral interactions (either charge exchange – CEX – or momentum exchange – MEX) can be expected to modify both the baseline fields and the average residence time of neutals in the ionization region. This effect becomes more important for higher ion densities and higher acceleration potentials. To get an estimation of the importance of such effect, we propose here a simple analysis of the Knudsen number based on the ion-neutral collisionality. We denote by

the frequency at which a given neutral particle is hit by a fast-moving ion. We consider an ion density , a cross-section and an ion velocity (ions are created at a low velocity and typically accelerate up to, say, , see for example [24]), giving . The residence time of a neutral particle inside the channel can be roughly estimated as , with and , such that . From this simplified analysis, the mean ion-neutral collision time and the neutral transit time appear comparable . Therefore, ion-neutral collisions are not the predominant contribution to the problem, but completely neglecting them, as we do in the present work, is a strong assumption, that should be relaxed in a future work.

4 Test-particle Monte Carlo simulations

Obtaining the residence time of neutral particles inside a region could be easily done by the DSMC method. However, such feature is usually not included in standard software, that would need in most cases some ad-hoc modifications, such as attaching one additional field –namely, the residence time– to the particle object, and updating it during the particle motion. In the present work we pursue a different strategy, referred to as test-particle Monte Carlo (TPMC), where a single test particle is moved in a background of gas, whose properties are given (known from the previous DSMC simulations).

The procedure goes as follows. First, a target particle is injected from the gas feed, with velocity sampled from a Maxwellian distribution at the injection conditions. Then, the particle is moved ballistically for a time step, detecting wall collisions, if any. The test particle collides with the background (assumed Maxwellian) with a probability that depends on the local density and temperature fields. Finally, when the particle leaves the channel, the residence time is saved and a new target particle is processed.

This strategy is ultimately a simplified version of the DSMC method, where particle collisions happen with a given background and not between simulated particles. As a result, the scheme implementation is trivial and heavily parallelizable. This procedure is particularly suited for describing free molecular flows, where one simply has to remove the gas-collision phase. This can be useful for lower propellant flow rates and/or thrusters with a larger surface-to-volume ratio, when the Knudsen number is large (see Appendix A). No DSMC computations for the background are needed in this case.

Figure 4: Cartesian and cylindrical reference systems for the test-particle Monte Carlo simulations.

The reference system for the simulations is shown in Fig. 4. Advection is performed in the three directions , and after each time step, the particle position is rotated back into the plane, as to impose axial symmetry. First, cylindrical coordinates are found by


The particle position vector is rotated back to the symmetry plane by the angle



and the same rotation is performed on the particle velocity vector :


When a particle crosses the walls during the advection phase, a wall-collision is detected: the particle is reflected into the domain, sampling its new velocity from a half-Maxwellian, characterized by the wall temperature. This coincides with assuming a diffusive scattering kernel, with unitary accommodation coefficient. The new particle position after the wall scattering event is computed accordingly.

The algorithm is finished if one wishes to consider a collisionless gas and only particle-surface interactions. In the collisionless case, the resulting residence times would scale linearly with the thruster size, and results obtained for a given configuration could easily be extrapolated to scaled geometries. Residence times can also be expected to scale with the square root of the surface temperatures, with higher temperatures causing higher thermal velocities and lower residence times. However, having a different anode and wall temperatures make the scaling less trivial.

Collisions of the target gas with the background are performed as follows. First, density, velocity and temperature are imported from the DSMC simulation (see for example some fields in Fig. 5). The values are provided as tabulated on the

plane as a 2D matrix and are sampled at the current particle location via bilinear interpolation (see for example

[25]). This gives the local background gas density , average velocity and temperature . The relative velocity between the particle and the background is .

Figure 5: Density and temperature profiles for case B30_1, warm anode, extracted from a strip of the DSMC simulation.

Collisions are then evaluated assuming a hard-sphere (HS) collision model. For a Maxwellian background and HS collisions, the velocity-dependent collision frequency is easily derived [26] and reads


with the xenon hard-sphere cross-section from [14], the Boltzmann constant and the mass of xenon background particles. The function reads


For consistency with the DSMC computations, one could employ the variable-hard-sphere (VHS) or the variable-soft-sphere (VSS) collision models instead of the simpler HS. In such cases, expressions of the velocity-dependent collision frequency are not trivial and one should instead compute it numerically. For the considered conditions, we have verified by numerical simulation that employing such models does not change considerably our results. Therefore, we only consider the HS model here. From the collision frequency, the probability that the test particle experiences a collision during the time step is then computed by (see for example [27])


It should be stressed that the time step needs to ensure a collision probability reasonably smaller than 1. Since the collision frequency is position-dependent, one possibility consists in updating the time step along the particle trajectory. A collision is performed if a random number satisfies the condition: . If it is the case, then a neutral colliding partner is sampled from the background Maxwellian distribution and a collision is performed by rotating the relative velocity vector by a random angle [28].

The discussed algorithm was implemented in two versions, one serial and one exploiting GPU acceleration with CUDA, available at [29]. As all particles are independent, the algorithm is heavily parallelizable. Running the CUDA version on an NVidia GTX 760 GPU, it was possible to obtain speedups of over 1600 times, for simulated particles, with respect to the serial program, running on an Intel CORE i5 processor with 12 GB of RAM. The wall time for simulating 1 million particles on a background field of grid points is in the order of a few seconds for the parallel code.

5 Results of test-particle Monte Carlo simulations

The scheme sketched in Section 4 is here employed to analyze the residence time of neutrals. The computational domain only includes the SPT-100 channel, such that particles that reach the plume cannot be scattered back and are effectively lost. Inside the channel, ionization typically happens in a region located before the exit plane. For the considered geometry, we assume the ionization region is located between and (see [30]). The computational domain and the ionization region are shown in Fig. 6. For a better clarity, the distribution is shown for the logarithm of residence times .

Figure 6: Computational domain for the test-particle Monte Carlo simulations (blue dashed box) and ionization region “iz” (red dashed box). Left: direct injection. Right: reversed injection.

First, we consider in Section 5.1 the free molecular limit, and run the test-particle Monte Carlo neglecting the gas-phase collision step. This step allows to observe some physical features that would be difficult to notice in the collisional case. Collisions are then introduced in Section 5.2, providing a more reliable figure of merit and the average residence time inside the ionization region is computed. Finally, this data is used in Section 5.3 to estimate the effect of reversed injection on the mass utilization efficiency. In all simulations, the particles injection is performed as described for the DSMC case (Section 3). Computations are performed on test particles.

5.1 Free-molecular limit

We start by neglecting all gas-phase collisions and consider 8 configurations: the direct injection cases A1, A2, and the reversed injection cases of Table 1. For every particle injected in the domain, the average residence time inside the whole channel is tracked. The resulting distribution of residence times is shown in Fig. 7.

Figure 7: Normalized distribution for the logarithm of residence times in the whole channel, in the free molecular limit. Effect of the injection angle (Left) and slit position (Right).

The distribution for the uniform injection appears rather simple. The slit injection, on the other hand, shows a larger and bi-modal shape: the lower residence-time peak is constituted by particles that, being injected with a higher velocity at the anode slit, reach the exit without colliding with the walls. This confirms the intuitive statement that the average velocity in the channel should be kept as low as possible for the sake of the residence time. The very same effect is visible in the warm anode case (Fig. 7-Bottom), where the high anode temperature further increases the low-residence-time bump, by increasing the amount of particles that leave the domain without wall collisions. Indeed, considering the cold anode case, a particle injected at the (sonic) average velocity of m/s takes s to travel a distance m, while for the warm anode case, the velocity is m/s, reducing the time to s. This matches the position of the low-residence-times peaks for the slit injection. Between the cold and the warm-anode cases, the high-residence-times tail is unchanged, these particles being scattered and thermally accommodated to the walls, which are set at K for both simulations.

Considering the reversed injection cases B30_1, B60_1 and B90_1, one can also observe a bi-modal shape. This time, the low-residence-time tail is constituted by the small fraction of particles with a positive axial velocity at injection due to the non-zero injection temperature. Such particles leave the domain prematurely and are effectively wasted. Indeed this peak: (1) decreases for more inclined injections and (2) is not affected by changes in the anode temperature.

From Fig. 7-Left, it appears that all cases perform worse than the uniform injection case A1. However, the injection at appears the best candidate in this pool. The performance can be increased by shifting the injection slit inwards. The effect of different slit axial positions is shown in Fig. 7-Right, where a fixed injection angle is chosen. As one could expect, moving the slit inward drastically reduces the amount of lost particles and thus increases the residence time inside the channel. In the limit of , one can expect to retrieve a behavior very similar to the uniform injection: most injected particles would collide on the anode and be scattered back in the chamber.

It should be stressed that, while we considered here the geometry of the SPT-100 thruster, our collisionless results map linearly to scaled geometries such as the SPT-50 thruster. Note that the mass flow rate does not have any influence in the present collisionless analysis, where all particles behave independently. Nonetheless, an analysis of the Knudsen number suggests that the collisionless assumptions are not far from reality, even in the presently considered conditions (see Appendix A). The analysis performed in this section gives useful insights in the behavior of neutrals in the different configurations. Identifying the best configuration requires further analysis and considering the finite width of the ionization region. This is done in the next section.

5.2 Collisional case and residence time in the ionization region

In this section, we address the full problem. Gas-phase collisions are introduced in the test-particle Monte Carlo algorithm, and the average residence time inside the ionization region is analyzed for a number of configurations. The previous collisionless analysis had shown the injection cases as the best candidates, as they minimize the fraction of prematurely lost particles. In this section we only consider such injection angle.

Differently from what we did in Section 5.1, the distribution of residence times will not be shown for the cases of this section. Indeed, since we only consider the ionization region and not the whole domain, all particles leaving the channel prematurely wouldn’t be accounted for, and this would make a comparison of the results much less intuitive. Only the average residence times will be considered, to which all lost particles simply do not contribute.

First, we perform a parametric analysis of the residence time at different slit positions, for the collisionless case. Then, we repeat the simulations for a selected number of configurations by introducing the background gas, obtained from additional DSMC simulations. The “residence time” variable attached to each simulated particle is updated only when the particle is inside the ionization region , and does not account for the time spent by the particle either before (closer to the anode) or after it.

Figure 8: Average residence time of neutrals in the ionization region. Vertical lines: assumed ionization region.

The average times in the ionization region are shown in Fig. 8. The effect of the anode temperature appears stronger for the collisionless cases, but much weaker when collisions are enabled. Considering large values of (slit close to the exit plane), a number of neutrals are lost prematurely, as previously observed in Section 5.1. Such particles never enter the ionization region, and the residence time is correspondingly low. On the other hand, positioning the slit between the ionization region and the anode does not bring substantial improvements. In the limit of , the residence time curve flattens and reaches the value obtained by direct injection. Indeed, all particles injected too early in the channel would bounce on the anode and cross the ionization region only once, when leaving the channel. The best improvements are obtained by centering the slit around , in the second half of the ionization region. In this configuration, the maximum gain in the average residence time is for the collisional case and about 40% in the collisionless limit (representative of lower and/or smaller thrusters).

Fig. 9 shows two sample trajectories obtained with the test-particle Monte Carlo method. The previous considerations on the value of the Knudsen number appear evident, and one can see that the mean free path is of the same order (or slightly smaller) than the channel dimension for the considered conditions. The free molecular regime is quickly approached for lower mass flow rates (resulting in a smaller density, and thus a larger mean free path), or for smaller thrusters (reducing the characteristic length); also see Appendix A.

Figure 9: Dump of two sample trajectories, for case B30_3. Gas feed is shown in red. Some trajectories appear curved due to the azimuthal motion out of the plane.

5.3 Estimation of the mass utilization efficiency

Although electron-neutral collisions have been neglected in terms of neutral dynamics (in view of the large mass ratio), we can obtain nonetheless an estimate of the mass utilization efficiency. This is done by assuming that, while inside the ionization region, neutral particles ionize by collisions with a cloud of background electrons with a given ionization frequency . Each simulated particle ionizes with a probability , that depends on the residence tau spent in the ionization region, as discussed in Eq. (2). The mass utilization efficiency represents the fraction of ionized over injected particles and can thus be found from


Our goal is now to estimate how much improve can be obtained by use of the reversed injection strategy. Let us first consider a thruster showing an efficiency when operated with a standard gas feed. From a simulation of such configuration, it is possible to obtain the array of residence times . A value of can be then determined for such value of , by solving iteratively Eq. (12). The value just found represents an average of the the plasma state in the ionization region. Then, a simulation of the reversed configuration is considered, and a new array of values is obtained. By plugging such values into Eq. (12), together with the obtained value of , gives the efficiency of the reversed configuration, for the same plasma conditions. The percent gain is then computed as . Maps of the percent efficiency gain are shown in Fig. 10 for both the collisionless and collisional cases, at different positions of the side slit. Only the reversed injection angle of is considered. In such maps, the collisional cases are representative of the simulated SPT-100 configuration with . The collisionless maps on the other hand can be directly applied to scaled versions of the SPT-100 geometry, provided that the ionization region is also in similarity.

Considering thrusters with a high initial efficiency , our computations suggest that the gain to be expected lies below the 15% in the best case scenario. On the other hand, for low-efficiency thrusters (say, ) the efficiency improve can reach the 20–30%.

Figure 10: Percent gain for different positions of the slit center and different values of . Injection angle .

6 Conclusions

In this work, we have investigated numerically the possibility to employ a reversed injection strategy as to increase the mass utilization efficiency of Hall Effect Thrusters. In place of injecting neutrals from the anode, the injector is located along the inner wall, closer to the channel exit, and injection proceeds backwards. Only the neutral gas flow is considered and interaction of neutral particles with electrons and ions are neglected. Such interactions can cause an increase in neutral temperature and average velocity, but such analysis is beyond the scope of this work. The only plasma effect considered here is the increased anode temperature, that strongly affects the neutral dynamics.

A preliminary analysis of the reversed injection configuration is first considered, assuming heavily idealized conditions (monoenergetic particles, specular wall reflection conditions and no gas-phase collisions). DSMC simulations were then performed assuming a more realistic diffuse wall scattering, finite temperature and enabling gas-phase collisions. The considered geometry is that of a SPT-100 thruster, with propellant mass flow rate mg/s. A baseline direct injection case with uniform injection, a slit injection and various reversed injection configurations are considered, for two different anode temperatures as to simulate more realistic operation conditions. Results confirm a strong influence of the injection strategy on the density and average velocity fields.

In order to study more into detail the residence time of neutral particles in the channel, a test-particle Monte Carlo scheme was also employed. From collisionless simulations, we have obtained the distributions of residence times of neutral particles inside the thruster channel. From this analysis, a reversed injection at an angle of appears to be the best candidate among the considered configurations, as it reduces the fraction of neutrals prematurely lost in the plume.

A parametric analysis was then performed for both collisionless and collisional cases, allowing to identify the best configuration. A slit located at appears to bring the best improvements, namely an increase of the residence time from 29% (collisional) to 40% (collisionless simulations) over the baseline direct injection case. Finally, the change in mass utilization efficiency is analized for different slit locations, allowing to identify the optimal configuration for a given initial thruster efficiency. Increases in the order of 5–15% are to be expected for thrusters with an efficiency of according to the operating conditions. A higher gain up to 20–30% could be expected for thrusters with a lower efficiency of .

As a concluding note, it should be stressed that the real performance of a gas feed configuration depends on the actual plasma properties inside the thruster, which strongly depend on the neutral flow itself. Also, besides the beneficial increase in neutral residence time, one should carefully consider for instance the effect of a higher neutral density on electron transport [16].


We thank Mr. Willca Villafana (CERFACS, Toulouse), Mr. Maciej Jakubczak (Institute of Plasma Physics and Laser Microfusion, Warsaw) and Dr. Federico Bariselli (VKI, Sint-Genesius-Rode) for the interesting discussions about this work.

Appendix A Analysis of the Knudsen number

The Knudsen number gives an indication of the rarefaction regime. Quantities and are usually global values for the problem. However, it is possible to compute maps of the local Knudsen number, by considering for example the local mean free path. This is shown in Fig. 11 for the uniform injection case, where the computed Knudsen number is based on the channel width , and the local mean free path is computed as , with the local number density and the HS cross-section for xenon. Fig. 11-Left shows the cold anode case, where all surfaces have the same temperature. In such conditions, we see a uniform and continuous decrease in the local Knudsen number, that directly reflects the progressive gas expansion. The anode temperature is seen to increase the local Knudsen number, due to the lower associated density. At the anode position, the warm emitted particles accommodate quickly to the colder walls, thus creating the curved contours visible in Fig. 11-Right, that reflect a higher density near the walls. It should be stressed that the present analysis of the Knudsen number is only indicative in terms of orders of magnitude.

The Knudsen number contours of Fig. 11 can be extrapolated to scaled thruster geometries and to different mass flow rates. In particular, for a given density, the Knudsen numbers are inverstly proportional to the thruster size, and would double by halving the chamber size. Also, at the considered rarefaction conditions, the chamber density is roughly linear with the mass flow rate , for a selected injection configuration. Hence, the mean free path is proportional to , and we expect the Knudsen number to scale linearly, as . This linear scaling with the mass flow rate is expected to eventually break down as low Knudsen numbers are approached, in the transitional or continuum regime, but this is likely to be far from the operating conditions of Hall thrusters.

Figure 11: Contours of the local Knudsen number from the DSMC simulations of Section 3, with and for the direct uniform injection configuration. Left: cold anode. Right: warm anode.


  • Zhurin et al. [1999] Zhurin, V. V., Kaufman, H. R., and Robinson, R. S., “Physics of Closed Drift Thrusters,” Plasma Sources Science and Technology, Vol. 8, No. 1, 1999, p. R1., URL
  • Morozov and Savelyev [2000] Morozov, A., and Savelyev, V., “Fundamentals of Stationary Plasma Thruster Theory,” Reviews of plasma physics, Springer, 2000, pp. 203–391.
  • Goebel and Katz [2008] Goebel, D. M., and Katz, I., Fundamentals of Electric Propulsion: Ion and Hall Thrusters, Vol. 1, John Wiley & Sons, 2008.
  • Hofer et al. [2006] Hofer, R. R., Jankovsky, R. S., and Gallimore, A. D., “High-Specific Impulse Hall Thrusters, Part 1: Influence of Current Density and Magnetic Field,” Journal of Propulsion and Power, Vol. 22, No. 4, 2006, pp. 721–731.
  • Brown et al. [2009] Brown, D. L., Larson, C. W., Beal, B. E., and Gallimore, A. D., “Methodology and Historical Perspective of a Hall Thruster Efficiency Analysis,” Journal of Propulsion and Power, Vol. 25, No. 6, 2009, pp. 1163–1177.
  • Guerrini et al. [1997] Guerrini, G., Michaut, C., Dudeck, M., and Bacal, M., “Parameter Analysis of Three Ion Thrusters (SPT-20, SPT-50 and A-3),” European Spacecraft Propulsion Conference, Vol. 398, 1997, p. 441.
  • Raitses et al. [1998] Raitses, Y., Ashkenazy, J., and Guelman, M., “Propellant Utilization in Hall Thrusters,” Journal of Propulsion and Power, Vol. 14, No. 2, 1998, pp. 247–253.
  • Komurasaki and Arakawa [1996] Komurasaki, K., and Arakawa, Y., “Performance Calculation of Hall Thrusters,” Acta astronautica, Vol. 38, No. 3, 1996, pp. 185–192.
  • Kurzyna et al. [2019] Kurzyna, J., Jakubczak, M., Szelecka, A., and Riazantsev, A., “Preliminary Tests of HIKHET Laboratory Model at IFPiLM,” 36th International Electric Propulsion Conference, 2019, pp. 2019–591.
  • Reid and Gallimore [2007] Reid, B. M., and Gallimore, A. D., “Review of Hall Thruster Neutral Flow Dynamics,” 30th International Electric Propulsion Conference, 2007, pp. 2007–038.
  • Garrigues et al. [2004] Garrigues, L., Hagelaar, G., Boniface, C., and Boeuf, J., “Optimized Atom Injection in a Hall Effect Thruster,” Applied physics letters, Vol. 85, No. 22, 2004, pp. 5460–5462.
  • Mitrofanova et al. [2011] Mitrofanova, O., Gnizdor, R. Y., Murashko, V., Koryakin, A., and Nesterenko, A., “New Generation of SPT-100,” 32nd International Electric Propulsion Conference, Wiesbaden, Germany, 2011.
  • Kogan [1969] Kogan, M. N., Rarefied Gas Dynamics, Plenum Press, New York, 1969.
  • Bird and Brady [1994] Bird, G. A., and Brady, J., Molecular Gas Dynamics and the Direct Simulation of Gas Flows, Vol. 42, Clarendon press Oxford, 1994.
  • Araki [2020] Araki, S. J., “Radiosity View Factor Model for Sources with General Distribution,” Journal of Computational Physics, Vol. 406, 2020, p. 109146.
  • Book and Walker [2010] Book, C. F., and Walker, M. L., “Effect of Anode Temperature on Hall Thruster Performance,” Journal of Propulsion and Power, Vol. 26, No. 5, 2010, pp. 1036–1044.
  • Martinez and Walker [2013] Martinez, R. A., and Walker, M. L., “Propellant Thermal Management Effect on Neutral Residence Time in Low-Voltage Hall Thrusters,” Journal of Propulsion and Power, Vol. 29, No. 3, 2013, pp. 528–539.
  • Boccelli et al. [2020a]

    Boccelli, S., Giroux, F., Magin, T., Groth, C., and McDonald, J., “A 14-Moment Maximum-Entropy Description of Electrons in Crossed Electric and Magnetic Fields,”

    Physics of Plasmas, Vol. 27, No. 12, 2020a, p. 123506.
  • Plimpton et al. [2019] Plimpton, S., Moore, S., Borner, A., Stagg, A., Koehler, T., Torczynski, J., and Gallis, M., “Direct Simulation Monte Carlo on Petaflop Supercomputers and Beyond,” Physics of Fluids, Vol. 31, No. 8, 2019, p. 086101.
  • Mazouffre et al. [2006] Mazouffre, S., Echegut, P., and Dudeck, M., “A Calibrated Infrared Imaging Study on the Steady State Thermal Behaviour of Hall Effect Thrusters,” Plasma Sources Science and Technology, Vol. 16, No. 1, 2006, p. 13.
  • Huang et al. [2011] Huang, W., Gallimore, A. D., and Hofer, R. R., “Neutral Flow Evolution in a Six-Kilowatt Hall Thruster,” Journal of Propulsion and Power, Vol. 27, No. 3, 2011, pp. 553–563.
  • Sharipov [2004] Sharipov, F., “Numerical Simulation of Rarefied Gas Flow Through a Thin Orifice,” Journal of Fluid Mechanics, Vol. 518, 2004, p. 35.
  • Koo and Boyd [2003] Koo, J. W., and Boyd, I. D., “Computational Model of an SPT-100 Thruster,” Proc. of 28th International Electric Propulsion Conference, 2003.
  • Boccelli et al. [2020b] Boccelli, S., Charoy, T., Alvarez Laguna, A., Chabert, P., Bourdon, A., and Magin, T. E., “Collisionless Ion Modeling in Hall Thrusters: Analytical Axial Velocity Distribution Function and Heat Flux Closures,” Physics of Plasmas, Vol. 27, No. 7, 2020b, p. 073506.
  • William [2007] William, H. P., “Numerical Recipes,” The art of scientific computing, Vol. 745, 2007.
  • Cercignani [1988] Cercignani, C., The Boltzmann Equation and its Applications, Springer, New York, 1988.
  • Vahedi and Surendra [1995] Vahedi, V., and Surendra, M., “A Monte Carlo Collision Model for the Particle-in-Cell Method: Applications to Argon and Oxygen Discharges,” Computer Physics Communications, Vol. 87, No. 1-2, 1995, pp. 179–198.
  • Garcia [2000] Garcia, A. L., Numerical Methods for Physics, Prentice Hall, Englewood Cliffs, NJ, 2000.
  • Boccelli [2020] Boccelli, S., “TPMC-2D,”, 2020.
  • Ahedo et al. [2001] Ahedo, E., Martınez-Cerezo, P., and Martınez-Sánchez, M., “One-Dimensional Model of the Plasma Flow in a Hall Thruster,” Physics of Plasmas, Vol. 8, No. 6, 2001, pp. 3058–3068.