Influence of a power supply model on simulated Hall thruster discharge voltage oscillations

Hall effect thrusters are spacecraft propulsion devices that use electricity to ionize and accelerate neutral propellant gas. The thruster consists of an annular or a cylindrical cavity with one end open to the ambient environment. The capped upstream face houses the anode which also serves as the propellant inlet. An externally mounted cathode provides a source of electrons, some of which flow towards the anode, completing the circuit. Magnets placed along the thruster circumference generate a predominantly radial magnetic field to retard the axial motion of electrons. The ions created in this region of increased electron density are accelerated out of the device by the potential gradient established by the two electrodes, generating thrust. One interesting aspect of a Hall thruster operation is that despite being powered by a direct current (DC) power supply, the discharge current exhibits highly oscillatory behavior. While the dominant low-frequency term, the so-called “breathing mode,” arises from the prey-predator nature of the ionization/convection process, many other contributions exist in the complex topology defined by the imposed electric and magnetic fields. A thorough review of plasma oscillations present in a Hall thruster can be found in Choueiri.¹ 

To date, virtually all Hall thruster simulations have utilized a simplistic model of the power processing unit (PPU). This model assumes that the voltage difference between the anode and the cathode remains constant. Unfortunately, this simplification fails to take into account the impedance effects introduced by the finite length electrical harness used to power the thruster. Given the nonlinear nature of Hall thruster plasma, there is an obvious capacity for interaction between the dynamics of the device with the dynamics of the electrical circuit powering these devices. This impact can be significant in ground testing since non-vacuum rated PPUs must be placed outside the vacuum chamber. For large chambers, the connecting electrical harness can easily exceed ten meters, leading to significant discrepancies between ground testing and flight operations. These discrepancies were investigated recently by Pinero,² who studied the impact of finite harness impedance on the discharge voltage ripple when driven by an predefined, experimentally-derived HET current source.

The purpose of this paper is to continue this investigation by studying finite impedance harness effects on a physics-based model of an HET in a self-consistent manner. This is accomplished by developing a dynamic power supply model in the industry-standard Hall thruster simulation code HPHall.³ The code is applied to a numerical model of the SPT-100 thruster.⁴ Several simulations are run from identical initial conditions with the only varying factor being the length of the simulated electrical harness. Simulations clearly identify an impedance threshold where HPHall begins to interact significantly with the harness fluctuations.

Before continuing, it is important to summarize the mathematical formulation implemented in HPHall. HPHall is a hybrid Particle in Cell (PIC)⁵ code in which ions and neutrals are represented by simulation particles, but electrons are treated as fluid. Positions and velocities of the heavy species are updated from the equations of motion. Ions are not magnetized and thus the reduced form of the Lorentz force $F→=qE→$ is used. The electric field is obtained from the electrostatic assumption, $E→=−∇ϕ$⁠. The method used for computing plasma potential ϕ is problem specific.

Inside the device, we can write the current balance I_d = I_i + I_e where the terms correspond to the discharge current reaching the anode, and the ion and electron currents in the beam. The ion current term is computed directly from particles, $j→i≡I→i/A=eZiniu→i$⁠, where Z_i is the average ionization state. On the other hand, electron motion is is governed by the Generalized Ohm’s Law,

with σ_e being the electron conductivity, and p_e = n_ekT_e is the electron partial pressure. Because electrons are magnetized, their dynamics can be decoupled into the motion along and across the magnetic field line. Ignoring magnetic gradient effects, σ_e,∥ ≫ σ_e,⊥. Electrons therefore respond to tangential disturbances almost instantaneously. This leads to two important simplifications. First, Equation 1 reduces to $ene∂ϕ/∂t^=∂pe/∂t^$ in the tangential direction, and second, it allows us to assume constant electron temperature along each field line. This simplification ignores the localized cooling associated with secondary electron emission. Therefore, as originally outlined by Morozov,⁶ the electric potential along a field line can be written in terms of a “thermalized potential,” ϕ = ϕ+ (kT_e/e) ln(n_e). With quasineutrality, n_e = n_i, the task of computing potential within the discharge reduces to finding ϕ^* and T_e, both of which are constant along each magnetic field streamline.

Following Fife, Equation 1 can be rewritten for the normal direction as

where μ_e,⊥ ≡ σ_e,⊥/(en_e) is the cross-field electron mobility. This term is generally computed by considering the classical, collision-driven transport, and overlaying it with an “anomalous” Bohm term. The optimization of Hall thruster simulations generally involves tuning the relative strength of the Bohm term, which may be location-specific.⁷ The cross-field electron velocity $ue,n^$ is one of the terms feeding into the energy equation,

Here S is a volumetric rate of change of energy density due to elastic and inelastic processes such as heavy-species collisions and ionization. Solution of this equation provides the variation of electron temperature. The challenge in using Equation 2 rests in that it contains two unknowns: $ue,n^$ and ϕ^*. In order to close the system, we use the current balance and let $ue,n^=(Ii−Id)/(eneA)$⁠. Equation 2 thus governs the spatial evolution of ∂ϕ^*/∂n. This formulation however makes the discharge current an input variable. HPHall’s potential solver begins with an initial guess for I_d to compute T_e and the corresponding $∂ϕ*/∂n^$⁠. This derivative is then numerically integrated to obtain the potential on the anode, taking into account a model-based correction for the anode sheath. The resulting discharge voltage, V_d = V_a − V_c is compared to the desired value. A Newton-Raphson solver subsequently produces a new guess for I_d, and the solver iterates until the desired discharge voltage is satisfied within some tolerance, ϕ_tol = 0.2 V. In other words, it is possible to treat the entire HPHall potential solver as a black box which computes the discharge current given some input discharge voltage, $Id=HVd,X$⁠, where X represents the instantaneous plasma configuration of the thruster. It is important to realize that even if V_d remains constant, I_d varies based on the temporal evolution of the plasma.

Instead of assuming constant discharge voltage, we can consider the circuit sketched in Figure 1 and solve for the time-dependent variation in discharge voltage due to the finite harness impedance. This simple circuit models the properties of the harness using a lumped element model suggested by Pinero.² Note that this approach inherently assumes that the entire atmospheric leg of the electrical system, from the PPU to the breakout box on the outside of the chamber, provides an idealized constant voltage power supply.

In our circuit, the cathode potential is zero, and hence V_d = V_A. R, L, and C are the total resistance, inductance, and capacitance of the line. However, it is more convenient to specify these properties as values per unit length. We can then let R = rx, L = lx, and C = cx, where x is the harness length, allowing us to easily study the effect of the harness length by changing a single variable. Assuming that the power supply maintains constant voltage V_s, the anode voltage is given by

Current balance at node 1 in Figure 1 leads to

The current through the capacitor is given by I_c = CdV_a/dt. The discharge current I_d is modeled as an arbitrary, time varying function of discharge voltage, I_d = H(V_a). We next use the standard finite difference discretization to rewrite the derivative terms. We use first-order differencing here, but a higher order scheme could be easily implemented as part of future work. We label values at the present and the prior time step with superscript k and k − 1, respectively. Equation 4 can then be rewritten as

In this formulation, all known terms have been moved to the right side, and equation 5 was used to rewrite the power supply current in terms of I_c and I_d. The power supply logic utilizes the same Δt as HPHall.

Similarly, we have

The resulting system can be written in matrix notation $Au→=b→(u→)$ with

and

This non-linear system is easily solved using the Newton-Raphson (NR) method, $v→m+1=v→m+J−1F→v→m$ where $F→v→m=Av→m−b→v→m$ and J_i,j = ∂F_i/∂v_j = A_i,j − ∂b_i/∂v_j is the Jacobian of $F→$⁠. The Jacobian can also be written as J = A − P, with

Here $vm≡ukm$ is the solution vector at the m^th iteration of the NR solver and simulation time step k. The solver iterates until ||v_m − v_m−1|| < ϵ_tol. The term ∂H(V_a)/∂V_a is computed from [H(V_a + ΔV_a) − H(V_a)]/ΔV_a with ΔV_a = 2ϕ_tol. Therefore, each iteration of the PPU NR solver requires two calls of the HPHall potential solver.

The above model was next applied to a simulation of the SPT-100 thruster. This 1.35 kW thruster, initially developed in the Soviet Union, is one of the earliest prototypes of the Hall thruster (also known as Stationary Plasma Thruster) technology. Due to its relatively simple design and the availability of magnetic circuit topology, this thruster is a common test-bed for Hall thruster simulation codes. The simulations were performed using the nominal 300 V power supply voltage, which, in the native implementation, translates to V_d = 300 V. The simulations were run for 100,000 Δt = 5 × 10⁻⁸ s time steps. Figure 2 shows a snapshot of the instantaneous results. The horizontal axis is the thrust direction and the vertical axis is the radius. These plots demonstrate the oscillatory nature of a Hall thruster by plotting the instantaneous ion density and plasma potential at two different simulation time steps. Note that for this entire investigation, the magnetic field is held constant, reflecting a non-flight-like configuration in which the magnetic circuit is powered separately from the discharge circuit.

Figure 3 compares the typical discharge voltage oscillations. These cases used the values given by Pinero: r = 5 mΩ/m, l = 380 nH/m, and c = 59 pF/m. Each section plots the discharge voltage over 2000 time steps. As can be seen, the magnitude of the oscillations increases with harness length x. The non-zero oscillations with x=0 are due to the finite ϕ_tol = 0.1 V. The magnitude of these oscillations increases proportionally with harness inductance. In this plot we also show a case in which the inductance was decreased by 50%. The resulting ripples obtained with 10 m, 0.5L is comparable to the 5 m case with 1L. The discharge voltage oscillations can also be compared by computing the root mean square (RMS) response for each harness length. This is done in Figure 4. This was done for both the full and the halved value of harness inductance. The figure also plots the data from Figure 12 in Pinero, which was scaled by a factor of three. The need for some scaling was expected, since the Pinero study was made for the 600 V, 12.5 kW HERMeS thruster. Our self-consistent model reproduces the same linear trend that was reported previously. However, in both self-consistent simulations, the discharge voltage ripple begins to exhibit nonlinear behavior for sufficiently high harness impedance.

We can similarly compare harness effects on discharge current I_d in Figure 5 and on the power spectrum of the AC component of the discharge current in Figure 6. The major observation is that the dominant low frequency peak at 20 kHz is not strongly affected by the harness length - indicating that the reactance of the harness has not fundamentally changed the breathing mode behavior of the numerical HET discharge - despite the obvious changes to the ripple voltage. Second, there is some dampening of the high frequency peak at 2.5 MHz but this is consistent with the low-pass filter (cutoff frequencies between 1.3 Mz/25 m harness and 16.8 MHz/2 m harness) nature of the harness model.

This paper demonstrated a self-consistent electrical harness model coupled to a physics-based Hall thruster code. It was used to confirm a prior experiment-based study that found that discharge voltage ripple increases with harness length. We also verified that decreasing the inductance of the longest harness studied in the prior study has the effect of decreasing rms voltage ripple. We were able to extend the analysis to study non-linear coupling effects of extremely long harness lengths. Beyond 10m, the effect of increased harness lengths gradually diminished and eventually reversed in sign to effectively decrease the voltage ripple.

This effort was funded in part by AFRL under an ERC, Inc. subcontract FA9300-15-C-0004. LB would like to acknowledge helpful discussions with Daniel Gerlach on performing circuit analysis.