Challenges with the self-consistent implementation of closure models for anomalous electron transport in fluid simulations of Hall thrusters

The performance of closure models for the anomalous electron transport when self-consistently implemented in a fluid model for a Hall effect thruster is investigated. This cross-field transport, which is orders of magnitude higher than classical collisional transport, is represented as an effective collision frequency. The proposed closure models relate this transport coefficient to local fluid properties of the plasma. Before implementation, the models are calibrated against values of the collision frequency inferred empirically from a 9 kW Hall thruster at 300 V and 15 A. It is found that even though closure models match the empirical collision frequency values, they diverge from these values when implemented self-consistently in a Hall thruster code. Possible drivers of this behavior are examined, including the role of non-linearity in the governing equations of the Hall thruster fluid model, artifacts from using time-averaged calibration data, and the non-uniqueness of the empirically-inferred collision frequencies. These results are discussed in the context of their implications for discovering and validating new closures necessary for enabling fully-predictive Hall thruster models.


Introduction
The ability to model Hall thrusters predictively is a key enabling capability for future developments in this technology [1,2]. These devices are a common type of electric propulsion * Author to whom any correspondence should be addressed.
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. employed for station-keeping and orbit-raising of spacecraft, and they are increasingly being baselined for deep space missions [3,4]. While this technology is already widely used, new mission applications are calling for Hall thrusters that operate in expanded power regimes and for unprecedentedly long lifetimes. Modeling and simulation tools are correspondingly becoming increasingly important to guide technology development in these new operating regimes. Despite their widespread operational use, however, there remain unresolved questions about key aspects of the governing physics. These questions have largely precluded the development of fully predictive engineering models.
The existence of so-called 'anomalous' cross-field electron transport has historically been the most critical of these poorly-understood processes [5]. Hall thrusters employ a radial magnetic field to confine electrons which is applied perpendicular to an axial electric field. While the ions are unmagnetized and therefore accelerated by the electric field, the electrons are magnetized and subsequently trapped in an azimuthal E×B drift. This confinement is critical to the efficient operation of the thruster. Classically, collisions with heavier species can interrupt the azimuthal drift, inducing electrons to move in the cross-field direction. However, the amount of observed electron transport is at least an order of magnitude more than that predicted by classical collisional theory [5]. This suggests there is an additional anomalous effect driving electrons across the magnetic field lines and increasing their mobility. Given the importance of the electron dynamics to the operation of Hall thrusters, including this effect is necessary to model these devices predictively.
Many researchers have therefore attempted to account for anomalous electron transport in simulations. The most direct approach has been to model the electron dynamics kinetically using particle-in-cell [6,7] or Eulerian direct kinetic [8][9][10] methods. These efforts have led to new insights about the nature and scaling of this enhanced transport. Most notably, results from these high-fidelity simulations have contributed to a growing consensus that microscale turbulence is a dominant driving mechanism for transport [5,[11][12][13]. With that said, kinetic methods are computationally expensive to a degree that makes them impractical for many engineering applications.
Modeling the electrons as a fluid provides an alternative method for simulating the Hall thruster at a reduced computational cost. The drawback of this approach that the classical form of the plasma fluid equations does not give rise to anomalous electron transport. Instead, it must be approximated using a closure model, i.e. a model that expresses the transport as a function of the bulk propertiesdensity, temperature, etc-that are already known from the fluid equations. The challenge with this approach is that presently, there is not sufficient knowledge of the physics to inform a predictive, reduced-fidelity approximation of the transport.
To address this shortcoming, many proposed closure models have been proposed to date based on varying hypotheses about the mechanisms driving the electron dynamics [14][15][16][17]. These models have considered effects ranging from near wall conductivity to the onset of microturbulence. More recently, we have also been exploring the use of data-driven (DD) methods to try to regress measurements of the anomalous transport and identify new closure models. In light of the wide range of proposed models to date, this invites the question as to how best to compare these models and evaluate their fidelity.
To this end, it is common to use 'empirical' profiles of the anomalous transport that have been calibrated against experimental data. These profiles are usually found by assuming the transport is static in time but varies spatially in the thruster. The profiles are then adjusted until key predictions of the code agree with experimental measurement [18][19][20][21][22][23]. This process yields useful insight into the state of the Hall thruster plasma in regions of the discharge that are not accessible by nonperturbative diagnostics, i.e. within the discharge channel. The resulting empirical profiles in turn are commonly used as references for evaluating the quality of proposed closure models [14][15][16][17].
Despite this practice, these empirically-determined profiles may not be sufficiently representative of the electron dynamics to be used as reference for benchmarking closures. Experimental measurements, for example, have shown that the transport is time-dependent, which is a departure from the static assumption of these profiles [24]. Additionally, it has been shown that these profiles are non-unique in that the empirical transport may vary by an order of magnitude or more in certain regions of the discharge without significantly impacting the predicted operations of the thruster [22]. This makes the validation of closure models challenging, as two models may agree well in regions of the discharge most impacted by anomalous transport, but diverge in regions less sensitive to the electron dynamics. These previous observations may suggest that closure models that accurately predict empirical profiles may not do so when self-consistently implemented in a simulation. In light of this possibility-which would significantly impact efforts to develop predictive Hall thruster models-there is an apparent need to evaluate the fidelity of closure models trained on empirical data.
The goal of this paper is to determine if empirical profiles are appropriate for validating closure models of electron transport in Hall thrusters. This paper is thus organized in the following way. In section 2 we review the problem of anomalous transport in Hall thrusters and define the closure problem in a fluid framework. We present examples of two closure models-one physics-based and one DD-which have shown to agree with empirical profiles in previous studies [16,17]. In section 3, we describe our methodology for incorporating these closures self-consistently into a fluid-based Hall thruster model. We also describe the experimental data we use for model validation. In section 4, we presents the results of our simulations compared to experimental measurements. In section 5, we discuss our findings and their implications for the evaluation of closure models.

Closure modeling in Hall effect thrusters
We present in this section an overview of modeling anomalous transport in Hall thrusters in terms of an enhanced effective collision frequency. We then describe a commonly-employed method to approximate this collision frequency empirically. Finally, we give an overview of previous attempts at closure models for this collision frequency, highlighting two of the most successful expressions that we investigate in subsequent sections.  Figure 1 depicts the basic principles of operation of a Hall thruster. In these axisymmetric devices, a DC discharge voltage V d is applied between the upstream anode and a downstream cathode, setting up an axial electric field ⃗ E. There is also a magnetic field ⃗ B applied in the radial direction across the channel, the strength of which is tailored such that only the electrons are magnetized. This crossed-field configuration gives rise to an electron current density ⃗ j de in the Hall (i.e. azimuthal,θ) direction. The electrons in the Hall current collide with neutral atoms injected through the anode and ionize them. The resulting ions follow the electric field and are accelerated out of the channel, producing thrust.

Anomalous transport as a diffusive process
In an ideal collisionless plasma, the crossed fields trap all electrons in the Hall drift. In practice, however, collisions can induce diffusion in the direction perpendicular to the magnetic field. If we neglect the inertia of the electrons (a common assumption for Hall thruster plasmas), we can represent this axial current density with a generalized Ohm's law (see [17]): Here j e⊥ is the component of the electron current density vector in the field-perpendicular direction, q is the electron charge, n e is the electron number density, m e is the electron mass, ω ce = q| ⃗ B|/m e is the electron cyclotron frequency in radians per second, ν e is the total electron collision frequency, E ⊥ is the perpendicular component of the electric field vector ⃗ E, and T e is the electron temperature in eV. In Hall thruster operation, the electrons are strongly magnetized, i.e. ν e ≪ ω ce . We therefore can simplify equation (1) to This result shows that the cross-field electron current scales with the electron collision frequency. Physically, this relationship reflects the fact that collisions act as an effective drag force in the azimuthal direction. This force, crossed with the applied magnetic field, drives the electrons to drift in the axial direction. Higher collision frequency translates to increased drag force and therefore enhanced crossed field current.
For equilibrium plasmas, there are classical expressions for the electron-ion and electron-neutral collision frequencies, ν c . However, incorporating these classical collisions into equation (2) yields electron currents that are orders of magnitude lower than what is observed in experiment. In order to resolve this discrepancy and bring simulations in line with experiment, it is common to introduce an additional 'anomalous' collision frequency ν an (or equivalently an anomalous mobility [15]) such that ν e → ν c + ν AN . By varying this term, it is possible to achieve a cross-field current commensurate with measurement. It is common practice to employ empirical methods to determine appropriate values for ν AN [22].

Empirically determining the anomalous collision frequency
The most widely-used method of empirically determining the anomalous collision frequency is to assume that it is constant in time but may vary spatially. The spatial variation can in turn be adjusted until key aspects of the simulation predictions, e.g. the acceleration of ions and discharge current, match experiment. This has resulted in a canonical 'shape' for ν AN . Figure 2(a) depicts an example of this typical shape in a Hall thruster plotted as a function of distance from the anode along channel centerline. We generated this result from numerical simulations of a magnetically shielded Hall thruster operating with xenon at 300 V and 15 A (described in more detail in section 3.1). We also show in this result the classical collision frequency ν c and electron cyclotron frequency ω ce . Figure 2 shows key plasma properties like the electric field strength, electron temperature, and plasma density.
The properties of the anomalous profile can be understood in the context of the trends in the plasma properties and background magnetic field. In the near-anode region, the neutral density and thus the electron-neutral collision frequency are high such that ν c ∼ ν AN . Inside the discharge channel, the anomalous collision frequency is typically an order of magnitude greater than the classical frequency. However, near the channel exit plane, in the region of peak magnetic field strength, the collision frequency decreases by an order of magnitude. The corresponding increase in electron resistivity in this location results in strong electric fields and enhanced Ohmic heating. This is reflected by the profiles exhibited in figure 2(b). Without this transport barrier, the predicted gradients in plasma properties from simulation would be much more shallow than observed experimentally. Downstream of the exit plane, the anomalous collision frequency increases by up to two orders of magnitude to reach a maximum value in the near-plume region. It then declines with the magnetic field strength with increasing distance from the thruster exit plane.
Adjusting the spatial variation of ν AN to reflect these typical features has allowed fluid-based Hall thruster simulations to match experiment to a high degree of fidelity. These types of calibrated simulations have been extensively leveraged for the design and qualification of Hall thrusters [25]. However, despite the success of these empirical profiles and the commonalities in their shapes, they are not extensible. An empirical profile calibrated to work on one thruster will not work with a different thruster or even a different operating condition on the same thruster. For each new thruster or operational state, we instead must infer a new empirical profile from experimental measurements of the thruster [20,23]. This limits our ability to leverage empirical profiles for fully predictive simulations.

Closure models for anomalous collision frequency
As an alternative to empirically-inferred profiles for the anomalous collision frequency, we ideally would find a closure model for this parameter. Closure in this context refers to a functional form for that depends on the local plasma properties. The nomenclature stems from the fact that by introducing this new collision frequency, we have opened the governing fluid equations, i.e. there are an equal number of equations and unknowns. Identifying an equation for ν AN that relates to the classical fluid properties e.g. ν AN (T e , n e , etc), closes the governing equations again.
Empirical profiles can be viewed as a type of closure model in this framework. Instead of the collision frequency depending on plasma properties, the assumption is that the collision frequency depends only on position. As previously established, this approach to closure is not predictive. Instead, it is common to propose physics-inspired closure models by making assumptions about the mechanisms driving the enhanced transport. These approaches often lead to analytical expressions for ν AN as a function of plasma properties [14][15][16].
In order to evaluate these closure models, it is common to compare them to empirically-inferred profiles [14,16,20]. More specifically, one can evaluate the closure model by using the time-averaged outputs from a simulation run with the empirical profile as inputs. The quality of the model is then evaluated based on how well the two curves match. As discussed in [17], following this metric, two of the most effective closure models are the first-principles (FP) model from Lafleur et al [16] and the DD model by Jorns [17].
In the FP model, the anomalous collision frequency is assumed to result from wave-induced drag due to an electron drift instability. By making assumptions about the growth and saturation of this instability, the authors proposed the model Here, ⃗ u i is the ion velocity vector, c s is the ion sound speed, v de = | ⃗ E|/| ⃗ B| is the electron ⃗ E × ⃗ B drift speed, and c 1 is a unitless constant that can be adjusted. In contrast, The DD model stems from an investigation by Jorns [17]. This expression was found by first compiling a dataset of empirically-inferred collision frequency profiles (per section 2.2) as a function of simulated plasma properties. They then applied symbolic regression to this dataset to propose functional forms for the anomalous collision frequency. This yielded several DD closure models that not only matched the data used to discover the models, but were also able to predict new empirical profiles  (3)) and data-driven [17] models (equation (4)) to the reference empirical profile depicted in figure 2(a).
with a higher degree of accuracy than existing physics-based closures. One of the most promising and simple models from this work is given by where c 1 and c 2 are adjustable unitless constants.
In figure 3, we compare these two closures to the calibrated empirical anomalous collision frequency profile from figure 2(a). For both models, we have used the time-averaged output from the resulting simulation to determine the plasma property values in equations (3) and (4). As can be seen, the FP model (with c 1 = 0.1245) and the DD model (with c 1 = 2.39 and c 2 = 3.32) both capture the major features of the empirically-inferred anomalous collision frequency profile and are able to correctly predict the location of the minimum in this parameter. There are discrepancies between the models both upstream and downstream of the minimum. However, as shown in [22], the anomalous collision frequency may vary by up to an order of magnitude in these regions without substantially impacting predictions for performance and plasma profiles.
Despite how well these models are able to predict the empirical anomalous collision frequency profile, it remains to be seen how they might perform when coupled with the rest of the plasma. This is the most rigorous measure of the predictive power of these models. To this end, we now turn to integrating these models into a full Hall thruster simulation.

Methodology
In this section, we overview the fluid model we employ as well as the details of the thruster we simulated. We then define several key metrics for evaluating the performance of our simulations.

Thruster
For the simulations in this work, we used the geometry and operating conditions of the H9 Hall thruster ( figure 4(a)). This is a 9 kW-class system developed in a collaboration between the University of Michigan, the Jet Propulsion Laboratory, and the Air Force Research Laboratory [26,27]. It is magnetically shielded and has a center-mounted cathode that is typically electrically tied to the thruster body. The thruster's nominal operating conditions range from 300 to 600 V and 10 to 20 A with a cathode flow fraction of 7% of the anode flow rate. We compare in this work results from simulations applied to the domain shown in figure 4(a) to experimental measurements of performance, thrust, efficiency, and specific impulse. Additionally, we considered the spatially-resolved mean ion velocity along channel centerline. These data were obtained via the method of laser-induced fluorescence (LIF) [28]. This type of local experimental plasma data is preferred for model validation as it is collected non-invasively. In contrast, probe-based methods to perform local measurements have shown in some cases to perturb the plasma properties [29]. All experimental measurements were collected while the H9 was operating at 300 V and 15 A on xenon [23,30].

Fluid-based model
We performed all simulations in this work with Hall2De, a state-of-the-art axisymmetric multi-fluid/particle-in-cell Hall thruster code developed at the Jet Propulsion Laboratory [21], which has been leveraged extensively for simulating Hall thrusters [2,21,22,[31][32][33]. In the version of the code we use, each ion species and the electrons are treated as fluids, while the collisionless neutral flow is treated using a line-of-sight view factor algorithm. The governing equations include continuity for both ions and electrons, a momentum equation for ions, and a generalized Ohm's law (section 2) for electrons. In addition to anomalous collisionality in this Ohm's law, the code simulates classical electron-neutral, electron-ion, and ion-neutral charge-exchange collisions. Ions are assumed to be isothermal with a user-provided temperature, but a governing equation is solved for the electron energy. In this equation, we do not explicitly add any additional energy loss or additions mechanisms due to wave-driven effects. However, the anomalous electron collision frequency does modify the energy equation by altering the collision frequency used to evaluate the cross-field thermal and electrical conductivities. The latter effect enhances Ohmic heating. The outputs of the code include two-dimensional (z-r, axial-radial) spatial distributions of ion and electron densities, electron current density, electron temperature, plasma potential, electric field, and ion velocity [21,22,31,32,[32][33][34][35].
We show a detailed picture of the simulation domain with labelled boundaries in figure 4(b). As can be seen, the region of interest extends eight thruster channel-lengths downstream from the anode in the axial direction and two thruster radii from the symmetry axis in the radial direction. We applied electrically-insulating boundary conditions on the channel walls and pole covers and zero-gradient boundary conditions at the outflow boundaries. At the cathode, we fixed the neutral flux, ionization fraction, and electron temperature.
We summarize in table 1 the simulation parameters that we employed for this work. The discharge voltage, mass flow rate, background pressure, and cathode flow fraction were chosen to match the experimental operating conditions at 300 V and 15 A for the H9. While Hall2De can accommodate up to four ion fluids, we restricted our simulation to two fluids-one population of ions emanating from the cathode and the other from the main beam. We allowed for four xenon charge states (neutral through triply-charged) for each ion fluid. We did not use the particle-in-cell feature of the code. We fixed the cathode electron temperature and ionization fraction to 3 eV and 5%, respectively, which are consistent with experimental measurements of the H9 cathode [36]. Hall2De employs a computational mesh which is aligned with the applied magnetic field in order to reduce numerical diffusion. Following the grid convergence study detailed in appendix, we selected a mesh containing 3925 cells. Small cell sizes near the magnetic poles of the field aligned mesh restricted the timestep to 9 ns in order to maintain numerical stability. We ran the simulations for two milliseconds of simulated time each, enough time to allow the plasma to settle into a stable oscillation and to minimize the effects of transients on time-averaged quantities. This amount of simulation time required 28 h of wall time when using 8 CPU cores.

Metrics for comparing models
In order to quantitatively compare the performance of the investigated models, we consider five metrics. The first metrics is the thrust, which is the amount of propulsive force generated by the device. To compute the thrust from the simulation results, we integrate the flux of axial momentum over the outflow boundaries of the simulation: In the above expression, ⃗ u f,j is the velocity vector of particles belonging to fluid f with charge j, u z,j,f is the component of that vector parallel to theẑ axis, n f,j is the number density of fluid f with charge j, dS is the differential surface area,n is the surface normal vector, and M is the mass of a xenon atom. The second metric we consider is the discharge current I d , which is the current carried by ions and electrons from anode to cathode. We compute I d by integrating the sum of the ion and electron currents over the anode boundary surface: Here, q is the fundamental charge and ⃗ j e is the electron current density vector. The third metric for comparison is the anode efficiency η a , which measures the fraction of the discharge power converted into useful thruster power. It is defined as whereṁ a is the mass flow rate injected through the anode and V d is the discharge voltage. The fourth metric we employ is the integrated velocity error (IVE). This measures how well the simulation predicts the ion velocity along the discharge channel centerline. We define this as In the above, z 0 and z N are the axial locations of the first and last ion velocity data-points, respectively, and u i,1 and u i,2 are the axial components of the ion velocity obtained from two different measurements or simulations. Higher values of the IVE correspond to worse agreement between simulation and experiment. For example, if the ion velocity curves differ by 10% on average, then the IVE should be approximately 0.1. The final metric we use in this work is the integrated anomalous collision frequency error (IAE). This is a measure of how well the anomalous collision frequency profiles from two simulations agree with each other. This has the same functional form as the IVE: In summary, these metrics allow us to evaluate the efficacy of our simulations. Armed with these formulations, we turn in the next section to comparing our simulation results to experiment and to an empirically-calibrated reference simulation.

Results
In this section, we first describe the results of our calibrated reference simulation. We then compare the performance of the closure models against that of the reference simulation.

Reference simulation
To assess the performance of our chosen closure models, we compare them to a calibrated reference simulation of the H9 at 300 V and 15 A. We choose to make comparisons to a reference simulation instead of directly to experiment in order to  directly assess the impact of the anomalous transport on the simulation. Since all factors except for the anomalous transport model have been held equal between the reference simulation and the two self-consistent simulations, any deviations are due to the transport model and not due to physics that may differ from experiment. The result presented here is the same as we used to motivate our discussion in section 2.2. We briefly discuss here how this simulation was generated with Hall2De. In particular, we obtained the empirical anomalous collision frequency profile by varying its spatial dependence iteratively until the discharge current was within 0.5 A of the experimental value (table 2) and the ion velocity along the channel centerline matched LIF measurements to within IVE < 0.1. This approach to determining the empirical profile is common practice for Hall2De simulations [22,32]. Achieving this match in our case required 29 iterations. We show the results of the reference simulation as a solid line in figure 6(b). Experimental data is depicted as discrete markers. We compare in table 2 performance metrics from our simulations to experimental data. As can be seen, the current and IVE of the reference simulation agree to within 8% of experiment, while the thrust is 12% lower and anode efficiency 20% lower than the measured values. This is a result of the fact that the calibration was based on matching experimental discharge current and ion velocity profile while we did not consider thrust. Comparable levels of reduced thrust and anode efficiency compared to experiment have been reported in previous Hall2De simulations with empirically-based anomalous collision frequency profiles [32].
With that said, as we discussed in section 2.2, despite under-predicting aspects of global performance, empirical profiles calibrated in this way to match ion velocity and current are often treated as surrogate measurements of the true anomalous collision frequency in these devices. Our goal is thus to determine how well closure models calibrated against such profiles perform when implemented directly into a Hall thruster code. We address this comparison in the following section.

Closure models
We consider in this section the results of the closure models introduced in section 2.3 when self-consistently implemented in Hall2De. As we discussed in this previous section, we calibrated these closures by evaluating them on the timeaveraged plasma parameter outputs of the reference simulation (figures 2(b) and 6(b)) and adjusting the model coefficients until we obtained quantitative agreement with the empirical anomalous collision frequency profile ( figure 2(a)). The underlying assumption in this case is that a model that matches the empirical profile will yield similar predictions.
In table 3, we compare the results of the closure models described in section 2.3 to the reference simulation. Note here that in table 3 we have computed the IVE with respect to the reference simulation rather than the experimental data. As can be seen, when we self-consistently implement these models, neither agree with the reference simulation on all five metrics described in section 3.3. The FP model (FP, equation (3)) shows the closest agreement with the reference discharge current and thrust-the thrust is 1 mN lower and discharge current 4.5 A higher than those of the reference simulation. On the other hand, the simulated anode efficiency is 11.5% lower than the reference. In contrast, the DD model (DD, equation (4)) predicts discharge current twice that of reference simulation and thrust 70 mN higher than the reference simulation, with anode efficiency 9% lower. Both simulations have high IVEs-0.307 for the FP model and 0.344 for the DD modelindicating poor agreement with the reference ion velocity curve. Lastly, in both cases the IAEs have nearly doubled compared to the non-self-consistently ('initial') implemented values. This suggests that the anomalous collision frequency profiles have significantly diverged from the original empirical profile. In summary, the major implication of the collected results in table 3 is that, despite the fact the closure models agree with the empirical profile, their self-consistent implementation does not yield results that match the empirical profile's outputs.
We can explain the global metrics reported in table 3 by examining the time-averaged anomalous collision frequency profiles of the two models. To this end, we show in figure 6(a) the resulting anomalous collision frequency profiles varying as a function of position after we have implemented them selfconsistently in the code. For comparison, we also show the empirical profile. We emphasize here that these are the selfconsistent results for anomalous collision frequency and differ from the spatial profiles shown in figure 3. In this earlier case, we simply trained the models against the outputs of the reference simulation with the empirically determined collision frequency profile. As we already had concluded globally from the IAEs, the two models have diverged from the empirical profile. The FP model still qualitatively matches the empirical profile in some respects-namely, the profile exhibits a minimum collision frequency with the right order of magnitudebut the location of this minimum is shifted downstream by half of a channel-length from the original value. The minimum collision frequency predicted by the FP model is also about 50% higher than the minimum of the empirical profile. This leads to enhanced electron transport in this region over the reference simulation, explaining the increased discharge current and reduced anode efficiency reported in table 3. In contrast, the DD profile shows an anomalous collision frequency profile which reaches a maximum, rather than a minimum, downstream from the exit plane with a collision frequency almost ten times higher than the empirical profile in this region. As in the FP model, this has the effect of increasing the current and reducing the efficiency. We can in part explain the differences in IVE between the two models by comparing the spatially resolved ion velocity curves to the anomalous collision frequency profiles in figure 6(a). In figure 6(b), we show the predicted ion velocity curves from the two models next to that of the reference simulation. Commensurate with its lower IVE, the ion velocity of the FP model matches the empirical result more closely than that of the DD model. Specifically, the FP model captures the steep slope of the empirical ion velocity curve while the DD model exhibits gradual ion acceleration. This results from the characteristics of the minimum in anomalous collision frequency exhibited by the FP model ( figure 6(a)). The point of minimum collision frequency corresponds to the point of maximum electron impedance, which in turn promotes a steep electric field. This field leads to rapid ion acceleration at the location coincident with this minimum. However, as can be seen, this profile is shifted downstream from the empirical profile, thus leading to the delayed ion acceleration profile. The relatively gradual acceleration of the DD model can similarly be understood by noting the high collision frequency near the thruster exit plane, which produces low electron impedance, a relaxed electric field, and thus slow ion acceleration.
In summary, we have shown that models trained on timeaveraged data to match empirical profiles diverge from simulation results obtained using the original empirical profile, in terms of both global performance metrics and spatially-resolved plasma properties. The differences in these quantities can be understood in terms of the shapes of the anomalous collision frequency profiles produced by the models.

Discussion
In this section, we discuss possible explanations for our results, outline the implications of our findings for developing models of anomalous transport, and provide recommendations for future Hall thruster model calibration efforts.

Divergence between self-consistent models and empirical profiles
We consider in the following three possible explanations for the divergence between modeling results when we calibrate the closures on an empirical profile versus self-consistent implementation. These include non-linearity of the governing equations, use of time-averaged data for calibration, and nonuniqueness of the empirical profile. As the models were compared to a calibrated reference simulation instead of directly to experimental data, we do not expect that differences between the Hall2De physics model and the true physics governing the thruster operation are able to explain our results.

Non-linearity of the governing equations.
As the dynamics of Hall thruster plasmas are highly non-linear, small initial variations in initial conditions can compound into larger differences over time. As a result, the small differences between our calibrated closure models and the empirical profile (figure 3) may have led to increasingly divergent behavior, resulting in the self-consistent profiles that do not match the initial empirically-informed ones. Moreover, in the absence of any differences between the closure models and the empirical profile (i.e. a perfect match), we suspect plasma oscillations inherent to Hall thrusters may still be sufficient to cause the model to diverge from the empirical profile.

Use of time-averaged data for calibration.
Another contribution to the diverging profiles may stem from the use of time-averaged plasma properties to calibrate our closure models. We can illustrate the problem with a simple example. As Hall thruster discharges are oscillatory, we can consider two plasma quantities, A(t) and B(t), which oscillate at a frequency ω, where B has a constant phase offset from A of ϕ radians. The time average operator is If A(t) and B(t) are sinusoidal functions with zero mean, we find after integration that where ϕ denotes the phase between A and B. The implication of this result that the product of the averages can be fundamentally different than the average of the products. Therefore, when we make calculations using time-averaged oscillatory quantities, the result may differ from what we would find if we had performed our computation at all times and then averaged the result. When evaluating models of anomalous transport of time-averaged data, there are thus likely to be averaging artifacts if the model involves products of oscillating quantities. We can apply this reasoning to illustrate the effect of time averaging on one of our closure models, the DD result (equation (4)). For this analysis, we neglect the ion sound speed c s for simplicity to arrive at We furthermore restrict the problem to one dimension along the axial coordinate z and restrict our analysis to locations immediately downstream of the ion stagnation point, so that ⃗ u i ≈ u i > 0 and ⃗ E ≈ E > 0. We can then collect factors which do not vary in time (e, |B|, and m e ) into a constant K.
After these simplifications, we have We can then assume both the ion velocity and electric field oscillate sinusoidally with a frequency ω and that the electric field oscillates out of phase with the ion velocity with a phase difference of 90 • . This is consistent with previous measurements of the breathing mode [37]. This yields where variables in brackets (⟨..⟩) denote mean quantities and those with tildes (..) denote amplitudes. This expression is also a periodic quantity with period 2π/ω. The time-averaged anomalous collision frequency is thus given by This yields a final result that scales as where ϵ is a positive number which is in general a function of space, E, and u i . This simple illustration demonstrates that we cannot replace the arguments of a closure model with timeaveraged values and assume the result will map to the time averaged collision frequency. This may in part explain why we can achieve quantitative agreement with the empirical profiles (that are implicitly time-averaged) while not matching the results with a self-consistent implementation. In practice, if this is a contributing factor, calibration of closure models would require time-dependent measurements.

Non-uniqueness of empirical profiles.
As a final note, we consider the implications of the fact that the empirical profiles we have used for calibration are non-unique. Mikellides and Lopez-Ortega have demonstrated [22] that stationary empirically-calibrated anomalous collision frequency models can vary significantly in certain aspects with only minor impacts on the simulation output. In particular, changing the anomalous collision frequency far upstream or far downstream of the thruster exit-plane has a smaller impact on the shape of the potential drop in the channel than changing it in the vicinity of the exit plane. Additionally, altering how the anomalous collision frequency changes away from the thruster centerline changes little about the observable device behavior beyond small changes in the plume divergence angle. Practically, this insensitivity is fortuitous as it allows more latitude in calibrating simulations to match data. However, it also poses significant challenges if we wish to use calibrated empirical profiles to inform our understanding of how anomalous transport actually varies in Hall thrusters. The non-uniqueness suggests there may be no effective 'ground truth' for calibrating or validating a closure model in certain regions of the thruster. To this point, given how insensitive the simulation is to the anomalous collision frequency in the region upstream of the exit plane, there may be multiple empirical profiles that give the same results as our reference simulation that we could use for regressing our closure models. These may ultimately have yielded improved predictive capabilities. One possible way to resolve this ambiguity would be to increase the amount of data used for determining the empirical profiles, e.g. additional velocity measurements off centerline and measuring other plasma properties. Alternatively, direct experimental measurements of the anomalous collision frequency would resolve this non-uniqueness [24].

DD generation of empirical profiles
While we have shown the limitations of the self-consistent implementation of closure models calibrated against empirical profiles, there is still a practical application for using these expressions to increase the efficiency of calibrating empirical models against new datasets. Given experimental measurements of the time-averaged plasma properties, we could use a given closure model to estimate what the shape of the empirical anomalous collision frequency profile should be. Indeed, in our previous work [17], we showed that closure models generated by regressing a dataset of empirical profiles were able to predict the empirical profiles of thrusters and operating conditions not in the training dataset. Reinforcing this point, the fit between the DD model and the H9 empirical profile displayed in figure 3 was obtained with no alteration to the model coefficients in [17], despite the fact that the DD model was not trained on data from this thruster. Using DD methods to speed the rate of empirical model calibration is an active area of research [38].

Implications for future closure model calibration and discovery
Self-consistent closure models for electron transport are critical for enabling predictive simulations. As the community continues to propose new models based on FP or DD methods, it is necessary to develop methods for calibrating and validating them. Our work has shown that using 'measurements' of the steady-state electron collision frequency inferred from calibrating against experimental data may not be sufficient for evaluating new models. With this in mind, as we look to future efforts, there are a number of possible strategies that may be employed for generating data sufficient for evaluating and calibrating new closure models. The most useful data for calibrating transport models would be direct experimental measurements of the anomalous collision frequency. Such measurements have been made in the past [24,36,39,40], though they are often subject to simplifying assumptions about the plasma state that introduce uncertainty in the values. Furthermore, these measurements often do not extend upstream of the exit plane, which leaves the problem of non-uniqueness stemming from the value of the transport in the near-anode region unresolved. With that said, as new experimental methods for characterizing the electron properties in situ become available, directly measuring the time and spatially-resolved transport may become more tractable [41].
In complement to developing methods to measure directly the anomalous collision frequency, it may be possible to infer the time-resolved electron dynamics by combining timeresolved data with real-time state estimation techniques [42]. A with the standard, empirical calibration approach we outlined in section 2.2, this method still relies on inferring the collision frequency indirectly by running a model iteratively and comparing the outputs to data. However, this approach has the advantage of the resulting empirical profile being timedependent. In light of the discussion in the previous section, this would provide a much higher-fidelity dataset for the calibration of closure models.
As an alternative approach, we could in principle calibrate and evaluate closure models by self-consistently implementing them in Hall thruster codes and iterating until the results match experiment. This has the advantage that it is not necessary to directly measure the collision frequency for comparison. The drawback is that this approach may be more time-intensive than evaluating the model against empirically calibrated or direct measurements. Reducedorder and multi-fidelity modeling may accelerate this process [43,44].
As a last comment, we remark that while our goal remains to find self-consistent closures for predictive-fluid modeling, empirically-inferred transport profiles remain critical to thruster development and qualification. They have been leveraged, for example, for discovering new designs to that increase thruster lifetime [21] and for evaluating plasma conditions in regions of the discharge which are difficult to measure [23]. With this in mind, we have shown that closure models calibrated against empirical collision frequency datasets may play an important role in accelerating the calibration procedure for new thrusters and operating conditions.

Conclusion
In this paper, we investigated the relationship between static, empirically-inferred profiles of Hall thruster anomalous collision frequency and self-consistent closure models calibrated on these profiles. We implemented two closure modelsone FP and one DD-into a state-of-the-art two-dimensional fluid Hall thruster code and compared the results to a validated reference simulation. We ultimately found that closure models tuned to match empirical profiles when using timeaveraged simulation data diverged from the empirical results when implemented self-consistently into simulations. These results may be attributed to a number of factors such as the non-linearity in the Hall thruster fluid model, averaging artifacts introduced by time-averaging, and the non-uniqueness of empirically-inferred collision frequency profiles. This inability to use steady-state profiles poses a challenge for the development and validation of new closure models for the anomalous transport.

Data availability statement
Some of the data is export-controlled. The data cannot be made publicly available upon publication due to legal restrictions preventing unrestricted public distribution. The data that support the findings of this study are available upon reasonable request from the authors.  in these performance metrics against the grid resolution. We find that the error decreases with increasing resolution. At the penultimate resolution of 3438 cells, all quantities differ by less than 5% from those at the finest studied resolution of 3925 cells. We note that the convergence is not monotonicincreasing the number of cells does not always decrease the error. This is likely due to the manual mesh generation procedure, which makes it difficult to ensure that the resolution smoothly increases across the whole domain as the number of cells is increased. Despite these caveats, our chosen resolution of 3925 cells appears to be sufficiently converged.