A simulation study is presented for a collisional, low-temperature dipole plasma encountered by rarefied, hypersonic neutral gas flow. A global model is developed and averaged over a toroidal control volume to describe the mass and energy exchange between the neutral stream and plasma. Simulations over a large range of magnet and freestream parameters reveal three distinct physical regimes that have significant bearing on the magnitude of plasma/flow interaction. The transitions between these regimes exhibit characteristics of resistive critical ionization, whereby the relative kinetic energy between plasma and neutral gas collisionally heats electrons, driving rapid and complete ionization of the gas. Two regime transitions are observed here with sudden exponential increases in plasma density occurring at velocity thresholds that depend on several energy loss mechanisms. The higher-velocity transition is a classical presentation of critical ionization where flow neutrals are ionized directly by plasma electrons. The other is a unique case in which charge exchange between ions and flow neutrals supplies both the particles and energy required to initiate critical ionization. This transition is distinct from any critical ionization effect reported in literature and indicates the existence of a lower critical velocity governed by collisional and diffusive effects as opposed to ionization energy losses only. Drag force on the magnetic field is considered by examining absorption and deflection of ionized flow by the dipole. The critical ionization thresholds increase the force on the magnet by up to two orders of magnitude compared to aerodynamic drag on an equivalently sized flow impediment.

## I. INTRODUCTION

Dipole plasmas interact with particle flows in many settings such as planetary magnetospheres,^{1,2} fusion reactors,^{3} or magnetic sails.^{4} One technology hinging on this interaction is plasma aerocapture, which uses an artificially seeded dipole plasma to initiate ionization, capture, and deflection of *in situ* atmospheric gas to decelerate a spacecraft.^{5–7} As this device flies through the atmosphere, the flow necessarily impinges on the full spatial extent of the plasma. Therefore, a global accounting of the mass, momentum, and energy transferred from the flow is warranted to understand the processes governing this coupled interaction.

Studies of plasma–neutral flow interaction often center on the critical ionization velocity effect (CIV). Critical ionization occurs when plasma and neutrals counterstream in a magnetic field and has been investigated in theory, experiments, and numerical simulations.^{8–12} This phenomenon, whereby background neutrals are rapidly ionized once the plasma flow reaches a threshold velocity, is ubiquitous from space to laboratory plasmas. No single comprehensive theory of critical ionization has yet been established that explains CIV observations across vastly different experimental and simulation conditions. The common mechanism among all observations is an energy feedback loop: ionization deposits the neutral kinetic energy with the ions, which transfer heat to the electrons to drive further ionization. Typically, the ion–electron energy transfer is assumed to come from collective processes.^{13} Some studies, however, have shown differences in the threshold velocity and rate of ionization when collisional heating processes are allowed^{14,15} as in a high density plasma where wave heating processes may be damped.^{12} Furthermore, neither theory nor experiment has addressed CIV in the context of a magnetic dipole geometry. Dipole plasmas are finite and inhomogeneous, whereas most analyses of critical ionization are one-dimensional and assume the plasma is spatially infinite.^{9} An investigation of the effects of finite geometry, diffusion, and collisional processes on CIV is warranted to characterize their role in the dipole-flow interaction.

In Part I,^{16} we investigated the interaction of a dipole plasma with neutral flow from a theoretical standpoint, finding that there are two critical velocity thresholds dependent on the time scales of flow energy capture, electron energy confinement, and ion–electron energy transfer. These time scales were treated as independent parameters although in reality their physics are not only coupled but highly nonlinear. In this paper, we present an analytical model of the flow–plasma interaction globally with a detailed accounting of specified mass and energy transfer processes. This is a control volume (CV) analysis coupling properties of the geometry, magnetic field, plasma, and flow to mass and energy exchange. The goal is to identify general physical scaling that will demonstrate the theoretical findings of Part I, contextualize dipole plasmas in the CIV literature, and inform the feasibility and design of plasma aerocapture techologies.

## II. GLOBAL PLASMA MODEL

The model consists of conservation equations that describe the exchange of mass and energy among the plasma and flow. The equations are normalized, averaged over a control volume, and solved numerically for the time evolution of total populations and energies. The volume-averaging approach eases computational load and therefore enables a broad parameter study from which we can identify physical scaling. Figure 1 provides an overview and reference for the problem formulation that is applied throughout this paper.

### A. Three-fluid model

For a control volume analysis, we start from basic conservation equations which allow us to account for all possible internal sources/sinks as well as flows across the boundary surface. Ions (subscript *i*), electrons (*e*), and neutrals are treated as separate species. We distinguish between stream neutrals (*sn*), which are yet-untouched flow particles, and secondary neutrals (2*n*), i.e., those generated in *i*-*sn* charge exchange (CX) reactions. Conservation of mass and energy for each population *j* is described by

where *n _{j}* is the number density, $\epsilon j$ is the energy density, $uj$ is the bulk fluid velocity, and $\u2202/\u2202t$ is the partial derivative with respect to time

*t*. The right hand sides of Eqs. (1) and (2) represent summation over all contributing sources and losses, where $x\u0307js$ is a rate of change of variable

*x*in species

*j*due to process

*s*. The energy density is defined following Meier and Shumlak

^{17}as

where *m _{j}* is the particle mass,

*p*is the species pressure, and

_{j}*γ*is the ratio of specific heats. With these general equations established, the problem is now constrained to enable formulation of an analytic model.

#### 1. Model constraints

We model the magnetic field as that of a loop of current centered at the origin with radius *r _{c}* in cylindrical coordinates $(r,\theta ,z)$. This is described by Jackson

^{18}as

where $e\u0302\theta $ is the azimuthal unit vector, *ψ* is the magnetic scalar flux function, *K* and *E* are the complete elliptic integrals of the first and second kind, respectively, and *κ* is a substitute variable. *B*_{0} is the magnetic field strength at the center of the loop. This magnetic field topology is shown in Fig. 1.

In order to facilitate a global approach over a detailed computational model, we reduce the problem complexity with a series of assumptions. Although these provide a simplistic approach, they enable us to attain broad physical insight from modeling before adding complexity in later analyses. They are as follows:

The system of equations is reduced by (1) considering only a single atomic constituent for the stream, ions, and secondary neutrals, (2) recognizing that energy conservation of stream neutrals is irrelevant since $\epsilon sn$ is, by definition, fixed, and (3) assuming quasineutrality and singly charged ions (

*n*=_{i}*n*). The latter argument equating ion and electron continuity implies that plasma diffusion is ambipolar._{e}We assume a dense ($ne>1016$ m

^{−3}), collisional plasma with electron Hall parameter $\Omega e\u226b1$. Electrons are thus magnetized and diffuse across field lines due to heavy particle collisions. By assuming diffusion is much slower than collisional time scales, energy distributions are Maxwellian and species act as ideal gases with isotropic temperature ($pj=njTj$).We take the magnetic pressure to dominate that of the plasma (i.e., low beta) and the stream dynamic pressure. For low-temperature plasma (

*T*< 100 eV), this implies a lower bound on applied field strength for a given plasma density, generally ranging from 10–100 mT in this study. This assumption allows us to neglect magnetic fields induced by either plasma motion or diamagnetic currents.^{19}By neglecting parallel currents

^{20,21}and toroidal drifts, we may assume the perpendicular diffusion is the only bulk fluid motion besides the freestream flow. The freestream is taken to be in free molecular flow and hypersonic so that its energy density is essentially kinetic.

A discussion of the validity and consequences of these assumptions is found in Sec. IV. We now develop the model equations with specified mass and energy transfer processes according to the applied constraints.

#### 2. Sources and losses

For this collisional, partially ionized plasma, the primary sources and sinks for mass and energy are ionization, charge exchange, Coulomb interaction, and plasma injection. Figure 2 serves as a visual guide to these mass and energy transfer processes among the model species. The treatment of each is detailed here and is largely supported by the plasma-neutral fluid model developed by Meier and Shumlak.^{17}

##### a. Ionization (iz)

Both stream neutrals and secondary neutrals can be ionized by electron impact. In general, ionization of neutral species *j* (superscript $iz,j$) occurs at a rate

denoting the rate coefficient of a reaction *s* as *R _{s}*. Neutral energy is absorbed by the ion population while electrons incur a loss of $\u03f5ion$, which is the effective ionization cost including losses due to scattering and excitation,

where $Mk/j\u2261mk/mj$ is the mass ratio of species *k* and *j*.

##### b. Charge Exchange (cx)

Ions undergo charge exchange (CX) with both neutral species, but only CX with the stream produces a net mass change by creating secondary neutrals from the ion population,

In general, each particle will maintain the energy it had before the reaction,

Here, $Qj,ks$ is a heat transfer from species *k* to *j* associated with reaction *s* [see Appendix A Eq. (A12)]. The first two terms in each of Eqs. (11) and (12) represent the exchange of thermal energy between ions and secondary neutrals during CX reactions. The third term in Eq. (11) representing power absorbed by the ions during stream CX is therefore an energy source for the system.

##### c. Collisions (coll)

Energy transfer in Coulomb collisions between ions and electrons is described by Braginskii^{22} who expresses the transfer rate as

where *τ _{ie}* is Braginskii's electron collision time [see Appendix A Eq. (A1)]. Note that this is the only energy channel between ions and electrons; we explicitly neglect collective processes, which are usually invoked to heat electrons in CIV models. Still, the model can be adapted to consider turbulent heating in later work by simply changing this thermalization term.

##### d. Injection (inj)

Finally, we account for the possibility of an artificial plasma source inside the control volume. In practice, this may represent a mechanism to seed ionization or modulate the plasma–flow interaction. It is modeled as a source of total particles *N* and total energy *E* that are distributed about the populations rapidly with respect to the time scales of interest. Ions and electrons are injected with a mass flow rate $\Gamma inj$ at fixed temperatures, $Ti,inj$ and $Te,inj$, such that the particle and energy rates are

where *j* represents ions or electrons and perfect gas enthalpy, $5nT/2$, is assumed.^{23}

#### 3. Diffusion

The flux divergence terms in Eqs. (1) and (2) account for any fluid motion across the control volume boundary. Recall that we assume diffusion is ambipolar, perpendicular to the field, and governed by electrons. Recall also that diffusion is the only bulk fluid motion and therefore the only boundary flux. The energy lost is thermal only, i.e., $\epsilon j\u2248pj/(\gamma \u22121)$. Considering only perpendicular diffusion, we have^{23}

where $\u2207\u22a5$ is the field line-normal gradient operator and *T _{j}* is the species temperature in energy units. Note that energy flux is scalar with particle flux due to the ideal gas assumption, $pj=njTj$.

In order to derive diffusion terms for the three populations, we must determine the proper coefficient *D _{j}* for each. With electrons governing the ambipolar plasma diffusion rate, it is appropriate to consider anomalous Bohm effects and electron collisions with both ions and secondary neutrals. Therefore, the effective diffusion coefficient for electrons (and, by symmetry, ions) is $De=Di=DB+Dc$. The Bohm diffusion coefficient,

*D*, and collisional coefficient,

_{B}*D*, are given by Lieberman and Lichtenberg

_{c}^{23}and Goldston and Rutherford,

^{24}respectively. Detailed descriptions are provided in Appendix A Eqs. (A3)–(A8). Secondary neutral diffusion arises from momentum-transfer collisions with ions and gas-kinetic collisions with other neutrals. We use a general form for gas diffusion given by Lieberman and Lichtenberg

^{23}and adapt it for these two collision types to obtain the secondary neutral diffusion coefficient, $D2n$, which is detailed in Appendix A Eqs. (A19)–(A20).

#### 4. Normalization

The three-fluid model will be simulated for a wide range of flow and plasma conditions. Therefore, we normalize the model so that self-similar physical scaling laws are readily identified. Because they are fixed inputs, freestream parameters (subscript $\u221e$) are primarily used for normalization. The following scheme is employed:

The stream is hypersonic with freestream velocity $u\u221ee\u0302z$ and freestream density $n\u221e$ (Fig. 1). Thus, the freestream energy density is primarily kinetic, i.e., $\epsilon \u221e\u2248msnn\u221eu\u221e2/2$. We define an effective temperature of the incoming flow, $T\u221e$, such that $\epsilon \u221e=3n\u221eT\u221e/2$. The normalization scheme in (18) yields a number of dimensionless parameters that govern the system physics as we will see in Secs. II B and II C.

### B. Control Volume Averaging

To describe the global system, the normalized equations are integrated over the control volume (CV). The CV representing the dipole plasma describes the spatial extent to which newly formed ions are trapped by the magnetic field. In other words, the CV is enclosed by a boundary inside which an ion with initial axial velocity $u\u221e$ is confined. This boundary has been previously derived^{6} (with more detail in Part I^{16}) by an analysis of ion trajectories in the dipole field. It is a toroidal volume enclosed by a surface of constant $\psi \u0302$,

where $\psi \u0302*$ is the normalized value of the magnetic flux function along the boundary contour and

is the characteristic ion Larmor radius normalized by magnet coil radius. Figure 1 shows the azimuthal cross section of the CV boundary surface, $S*$, which encloses the volume of the toroid, $V*$. Equations (19) and (20) demonstrate simple scaling between model parameters and the CV: as *ρ _{L}* increases, the CV shrinks.

The volume-averaging process is simplified considerably by assuming a particular form for the density distribution. Dipole plasmas are challenging to model with high physical fidelity. Power-law density scaling is commonly employed,^{20,25,26} i.e.,

where $nj,r$ is the invariant density on an arbitrary flux contour *ψ _{r}* and the parameter

*α*depends on a variety of factors,

^{25}but is usually between ∼1 and 10. Equation (21) implies that density is constant along magnetic field lines and, therefore, along our CV boundary. We assume that

*n*,

_{i}*n*, and $n2n$ are described by this profile, which is represented in Fig. 1. Any perturbations to this profile are assumed to be damped on time scales much shorter than the characteristic system evolution, so interspecies reactions ($\u221dnjnk$) do not result in higher order terms. While these assumptions sacrifice physical fidelity, they simplify the model such that spatial dependencies are functions of

_{e}*ψ*and

*B*only. This, combined with the fixed CV size, facilitates global modeling by reducing integrals over the CV to fixed values for a given set of input conditions. We leverage this to volume-average the model equations and drastically reduce computational load. A discussion of the validity of Eq. (21) and ways to build complexity in future analyses is provided in Sec. IV.

The dimensionless model equations after CV-averaging describe the rate of change of the total population and energy present in the volume. The model variables are the total normalized populations in the CV, $N\u0302j$, and the normalized temperatures, $T\u0302j$, of each tracked species as a function of time, $t\u0302$. We derive these by normalizing and integrating the conservation equation variables, *n _{j}* and $\epsilon j$,

Here, $E\u0302j$ is the normalized total energy of species *j* in the CV. Combining Eqs. (21) and (22) gives the useful relation

where $I\psi \alpha $ is a magnetic field integral defined in Appendix B Eq. (B1). Relation (24) allows us to describe the spatially resolved density profile, $n\u0302j$, from the CV-averaged solution, $N\u0302j$, through Eq. (21). Recall that $V*$ (defined by $\psi \u0302>\psi \u0302*$) is a function of *ρ _{L}* only and $\psi \u0302r$ is fixed.

Interactions between tracked species, of the form $n\u0307js\u221dnjnkRs$, depend spatially on the species densities which both follow Eq. (21). Such interactions include 2*n* ionization and *i*-2*n* charge exchange. After volume averaging, these terms have the form $N\u0307\u0302js\u221dn\u0302j,rn\u0302k,rI\psi 2\alpha /\tau \u0302s$, where $I\psi 2\alpha $ [defined in Appendix B Eq. (B2)] is a magnetic field integral resulting from taking the product of tracked species densities. $\tau \u0302s$ is a dimensionless parameter describing the characteristic time scale of reaction *s* compared to its transit time across the magnet, defined as

Interactions of tracked species with the neutral stream, having the form $n\u0307js\u221dnjnsnRs$, are more complicated because *n _{sn}* does not follow Eq. (21). The resultant volume-averaged interaction has the form $N\u0307\u0302js\u221dn\u0302e,rIsn/\tau \u0302s$, where the so-called “capture integral”

*I*involves integrating $n\u0302sn$. Its derivation is found in Appendix B, including an analytic solution to the stream neutral continuity equation yielding $n\u0302sn$.

_{sn}Finally, tracked species diffusion, with the form $n\u0307jdif\u221d\u2207\xb7(Dj\u2207\u22a5nj)$, depends on the different diffusion models described in Sec. II A 3. Volume averaging results in a common scaling for each, $N\u0307\u0302jdif\u221dD\u0302jn\u0302j,rID$. Here, $D\u0302j$ is a dimensionless, spatially invariant form of *D _{j}* specific to each diffusion model. Similarly,

*I*are magnetic field integrals that reflect the spatial scalings associated with

_{D}*D*,

_{B}*D*, $D2n$, and density gradients. The different scaling with

_{c}*B*in each

*D*yields three unique integrals

_{j}*I*,

_{B}*I*, and $I2n$ that are all functions of the magnetic field only [ Appendix B Eqs. (B6)–(B8)].

_{c}### C. Global model equations

The final global model equations are the result of applying the volume-averaging in Sec. II B to the conservation equations developed in Sec. II A. This system consists of five nonlinear ordinary differential equations and is solved for the time evolution of five unknowns: $N\u0302i,\u2009N\u03022n,\u2009T\u0302i,\u2009T\u0302e$, and $T\u03022n$. The equations describe ion/electron continuity (26), secondary neutral continuity (27), ion energy (29), electron energy (30), and secondary neutral energy (31) which are coupled via mass and energy transport terms. Note that the appearance of $n\u0302j,r$ is replaced with $N\u0302j/I\psi \alpha $ using Eq. (24) to highlight the model dependence on the five dimensionless variables and the CV integrals. Underbrackets are used as an aid to the reader indicating the physical process represented by each term.

#### 1. Continuity equations

Recall from Sec. II A that *n _{i}* =

*n*, so ions and electrons are described by the same continuity equation (presented here for ions). Recall also that the solution to the stream neutral continuity equation is already known [Eqs. (B4)–(B5)]. Thus, global model continuity consists only of that for ions and secondary neutrals,

_{e}Diffusion is a loss for both species, with $D\u0302B,\u2009D\u0302c$, and $D\u03022n$ being dimensionless, spatially invariant forms of the diffusion coefficients. These are defined in Appendix A Eqs. (A7), (A8), and (A20), respectively. The particle injection rate,

represents the volume-averaged plasma source injection as a fraction of the characteristic stream neutral flux.

Equations (26) and (27) reveal how the capture integral, *I _{sn}*, acts to couple the flow and plasma. There is a direct source of plasma from ionization of stream neutrals as well as an indirect source when CX with the stream creates 2

*n*s that are reionized. Both of these stream interaction processes are functions of

*I*and supply ions from outside the CV with no corresponding loss term from another model equation.

_{sn}#### 2. Energy equations

Conservation of energy is considerably more complex. Despite eliminating the need for *sn* energy balance in Sec. II A 1, tracking system energy requires three distinct, nonlinear equations describing ions, electrons, and secondary neutrals,

Notice that many energy terms result simply from the transfer of internal energy associated with the continuity terms. The dimensionless injected power,

is scalar with dimensionless particle injection, following Eqs. (14) and (15). Several more terms in Eqs. (26) and (27) result from energy transfer processes without an associated net particle transfer. For instance, CX reactions and Coulomb collisions conserve particle populations but are mechanisms for exchanging thermal energy between species [with $Q\u0302j,k$ and $\tau \u0302ie$ defined in Appendix A Eqs. (A12) and (A2), respectively]. Although ionization terms are similar to their continuity counterparts for *i* and 2*n*, the electrons incur an additional loss of

which is the effective ionization cost [ Appendix A Eq. (A16)] normalized by stream kinetic energy.

## III. RESULTS

The interaction is simulated for argon flow and plasma by solving Eqs. (26)–(27) and (29)–(31) over a time scale long enough to reach equilibrium. Details of the argon model are provided in Appendix A–C. The LSODA differential solver^{27} is used and outputs the time evolution of $N\u0302i,\u2009N\u03022n,\u2009T\u0302i,\u2009T\u0302e$, and $T\u03022n$. We first use these to examine the nature of the system as it evolves toward a steady state in a particular simulation. Then, we run an array of simulations across broad ranges of inputs $n\u221e,\u2009u\u221e$, *r _{c}*,

*B*

_{0}, and $\Gamma inj$, extracting the equilibrium solutions to examine physical scaling with these parameters.

In the following, the reference surface for computations such as Eq. (24) is $\psi \u0302r=\psi \u0302*$ unless otherwise noted. We use *α* = 4 for the density profiles given by Eq. (21). While this approach is a simplification, there are several justifications. First, this value of *α* commonly appears in the literature for the theoretical stationary state of a dipole plasma.^{21,25,28} Second, a constant *α* enables computation of the CV integrals independent of the numerical evolution of the system, a stated motivation for our global model approach. Finally, we observe that the results to follow are not qualitatively sensitive to the specific density profile employed. Previous work using $nj\u221dB$ exhibited the same general character of the plasma–flow interaction.^{6,7} Thus, the present implementation provides a lens to obtain physical insight into the interaction such that appropriate complexity can be developed later.

### A. Temporal behavior

To gain intuition on how the coupled mass and energy transfer processes in Sec. II A drive the system toward equilibrium, we first examine the initial time evolution of the interaction during a single simulation. Figure 3 shows the output of a simulation with $n\u221e=1018$ m^{−3}, $u\u221e=12$ km/s, $B0=0.25$ T, *r _{c}* = 1 m, and $\Gamma inj=5$ mg/s. At these conditions, $t\u0302=12$ corresponds to

*t*=

*1 ms. The initial conditions are $N\u0302i(0)=N\u03022n(0)=1$ and $T\u0302i(0)=T\u03022n(0)=T\u0302e(0)=0.1$. The populations $N\u0302i$ and $N\u03022n$ both build up gradually before suddenly increasing exponentially to their equilibria. This sudden increase is accompanied by a sharp drop in all species temperatures relative to the stream energy. The system effectively reaches equilibrium within about 20 ms ($t\u0302\u223c240$) although this time decreases if the initial conditions of $N\u0302i$ and $N\u03022n$ are increased. The steady-state populations and temperatures, however, are only sensitive to the scaling parameters ($n\u221e$, etc.) and not the initial conditions. Note that the time scale of evolution in Fig. 3 is much slower than the characteristic particle transit time ($t\u0302\u226b1$). This justifies our earlier assumption that density distributions can be described by their equilibrium profiles.*

Figure 4 shows what this density evolution looks like spatially. At first, the control volume contains very little plasma, and the flow is only affected by the extreme magnetic fields at the ring current [Fig. 4(a)]. After the sharp jump [Fig. 4(b)], the density has increased substantially and the flow is impeded in the region around the magnet. The plasma continues to absorb flow particles, eventually reaching a steady state [Fig. 4(c)] where its density is orders of magnitude higher than the freestream throughout the CV and it produces a wake several times the size of the magnet. Note at steady state that the neutral flow is entirely absorbed over the radial extent of the control volume.

Figure 5 reveals the underlying processes behind this evolution to steady state. Charge exchange between ions and stream neutrals begins depositing secondary neutrals in the control volume, which are lost entirely to diffusion for $t\u0302<100$ [Fig. 5(b)]. Meanwhile, ions and electrons build up slowly from ionization of stream neutrals during this time [Fig. 5(a)]. As the ion population grows, charge exchange becomes more frequent, producing more secondary neutrals. This results in a runaway effect where the growing reservoir of neutrals becomes an ionization source that causes an apparent critical ionization effect around $t\u0302=100$. Here the plasma population exponentially increases [Fig. 5(a)], driven by ionization of secondary neutrals produced from charge exchange with the stream. In steady state, the energetic electrons sustain this secondary neutral ionization [Fig. 5(e)]. Their energy is derived from collisions with ions that are heated by charge exchange with the stream [Fig. 5(c)]. In this manner, charge exchange provides both the neutral reservoir and the energy source for a critical ionization effect to take place in the control volume.

### B. Steady state regimes

Together, Figs. 3–5 provide comprehensive intuition about the nonlinear interdependence of the model equations and the coupled interaction between the plasma and flow. We now seek to identify physical regimes and scaling of the interactions by investigating the steady-state solutions to an array of simulations, each with different magnet and flow conditions. Figure 6 shows the dependence of the boundary density $n\u0302*$ (i.e., on $\psi \u0302*$) and average temperatures $T\u0302$ on the freestream velocity $u\u221e$. In contrast to the three tracked species, $n\u0302sn$ is not invariant along field lines, so we compute $n\u0302sn*$ at the furthest downstream point on $S*$ as a characteristic measure of *sn* density. Note in Fig. 6 that velocity is normalized by Alfvén's^{8} critical velocity,

where $Uiz$ is the first ionization energy ($ucr=8.7$ km/s for argon). Like in the temporal solutions, a critical ionization effect is apparent. At low velocities, very little of the flow is converted to plasma ($n\u0302sn*\u223c1$) and few ions are present. The model species are thermalized with the stream ($T\u0302\u223c1$) by CX interaction. As the velocity increases past $u\u221e\u223c0.4ucr$, there is a sudden exponential increase in plasma and secondary neutral density, while $n\u0302sn*$ begins gradually sinking to zero. This is accompanied by a sharp drop in the species temperatures. While plasma dominates in this velocity range, there is a significant fraction of secondary neutrals present, with $n\u03022n*$ roughly 10% of $n\u0302i*$. Increasing velocity past $u\u221e\u223c1.8ucr$, a third distinct regime is apparent where the plasma becomes fully ionized. Secondary neutral density effectively reduces to zero though the plasma maintains its high density. The plasma density spikes for a brief window where stream CX remains a significant plasma source (evidenced later by Fig. 8), and then gradually declines with increasing velocity due to increased diffusion and a shrinking CV [Eq. (20)]. Species temperatures approach the order of the stream energy and continue to rise with increasing velocity.

As with the temporal solutions of Figs. 3–5, the steady-state regimes observed in Fig. 6 are governed by the physical processes that connect the model equations. In Fig. 7, we show the steady-state contributions of mass and energy transfer terms that drive the equilibrium solutions as a function of the input $u\u221e$. In other words, all of the curves shown correspond to the steady-state values of right hand-side terms in Eqs. (26), (27), and (29)–(31). Black dotted lines indicate the transitions between the three regimes that we now characterize according to the main source of plasma sustainment.

**Injection (I) regime:** At low velocities, the only energy deposited into the plasma comes from the power injection term $P\u0302inj$ [Fig. 7(e)]. All injected mass and energy diffuses from the control volume. There is substantial charge exchange with the stream [Fig. 7(b)] but it makes no net contribution because the ionization time is much longer than the 2*n* diffusion time. Interestingly, even though stream charge exchange supplies new ions, the kinetic energy of flow neutrals is too low to heat the ion population. In fact, CX is mostly an energy sink in this regime [Fig. 7(c)]. Because the injected power and mass dominate over any plasma/flow interaction, we define this as the “injection regime.” Note that plasma sustainment in the I regime requires an artificial source, without which the plasma would extinguish.

**Charge exchange (CX) regime:** In the mid-velocity regime, a large spike is observed in the stream CX contribution to the ion energy equation [Fig. 7(c)], heating the ions with flow kinetic energy. This energy, in addition to that from CX with 2*n*s, is transferred to the electron population via Coulomb interaction. Electrons subsequently have enough energy to ionize some of the 2*n* reservoir resulting from the original stream CX reactions [Fig. 7(e)]. This forms a critical ionization feedback loop: ion–stream CX produces abundant secondary neutrals and simultaneously heats the electron population; electron ionization of 2*n*s augments the ion population, further increasing CX with the stream. This process is distinct from most models of critical ionization in that it is driven by charge exchange, so we refer to this as the “CX regime.” Note that negligible energy is lost in plasma diffusion, instead leaving the system entirely through $\u03f5\u0302ion$ (i.e., ionization and radiation) and 2*n* diffusion [Figs. 7(e) and 7(d), respectively]. This effect eases the energy burden on electrons to sustain critical ionization, ultimately enabling the low threshold velocity observed, consistent with the theoretical results of Part I.^{16}

**Critical ionization velocity (CIV) regime:** The second regime transition at high velocities occurs when the plasma starts ionizing the stream directly. Ionization supplies the bulk of ion energy by converting stream kinetic energy to the plasma [Fig. 7(c)]. This again serves to heat the electrons through Coulomb interaction, enabling them to continue ionizing stream particles directly. This feedback loop of energy transfer is the classic manifestation of the critical ionization phenomenon, so we term this the “CIV regime.” Because the electron temperature is higher in this regime, failed ionizations (scattering and excitations) are less likely so ionization losses are reduced. Therefore, diffusion becomes an equally important energy loss mechanism for electrons. Note also that $n\u0302i*$ now decreases with $u\u221e$ due to a combination of increased Bohm diffusion as $T\u0302e$ rises and shrinking of the ion-trapping region according to Eqs. (19)–(20).

Figures 6 and 7 clearly demonstrate that plasma sustainment is driven by injection, charge exchange, and ionization in three distinct regimes. We are primarily interested in the extent to which coupling between the flow and plasma drives this behavior. Recall from Sec. II B that $n\u0302e,r/\tau \u0302cap$ is the critical term that couples the plasma properties to the flow through *I _{sn}* [Eq. (B3)]. Whereas earlier this parameter was defined for an arbitrary contour

*ψ*, here we define it on the CV boundary, $\psi r=\psi *$, to represent the lower bound of reaction rate in the CV. We examine the behavior of $\zeta s*$ in the three regimes where

_{r}is the rate of reaction *s* (either CX or ionization) on stream neutrals at the CV boundary compared with the neutral rate of transit across the magnet. Thus, $\zeta s*\u22731$ implies significant flow–plasma interaction (and therefore mass and energy absorption) within the control volume.

Since the stream only undergoes CX or ionization by the plasma, let us define the total reaction parameter $\zeta tot*=\zeta cx,sn*+\zeta iz*$. These three rates are shown as a function of freestream velocity in Fig. 8. In the injection regime, overall reaction with the plasma is quite low, though CX dominates over ionization. This explains why Figs. 7(a) and 7(b) show substantial charge exchange with the stream but no significant accumulation of plasma in the injection regime. Figure 8 also demonstrates that increasing injected power directly impacts $\zeta tot*$ in the I regime. The CX regime is characterized by a sharp and monotonic increase in $\zeta cx,sn*$, which drives $\zeta tot*$ high enough to yield significant interaction between the plasma and flow. $\zeta iz*$ remains quite low until a drastic jump denoting the CIV regime. Here, stream ionization is so strong it drives $\zeta tot*$ to its maximum even while the CX rate declines. Figures 7 and 8 offer an understanding of the governing physics within these three regimes; we now examine the mechanisms by which the system transitions between these states.

### C. Regime transitions

The existence and nature of the flow interaction regimes described in Sec. III B strongly depend on the properties of the neutral flow and magnet. In Part I,^{16} we derived threshold criteria for the I (injection), CX, and CIV regimes as a function of the normalized electron diffusion timescale ($\tau \u0302e$), particle capture time scale ($\tau \u0302cap$), ion–electron energy transfer time scale ($\tau \u0302ie$), and the ratio of the effective ionization energy cost to the kinetic energy of stream neutrals ($\u03f5\u0302ion$). Here, we examine these thresholds from the perspective of the global modeling results to understand the transition between regimes now that particular mass and energy transfer processes are considered.

In Fig. 9(a), we show how the I, CX, and CIV regimes manifest in $N\u0302i$ as a function of $u\u221e$ for a range of magnetic field strengths. For four cases of $B0\u226525$ mT, the plasma undergoes both I-CX and CX-CIV mode transitions as $u\u221e$ is increased. Interestingly, the I-CX threshold velocity decreases as *B*_{0} increases, exhibiting values both less than and greater than $ucr$. Unlike the I-CX mode transition, the threshold value of $u\u221e$ for the CX-CIV mode transition is always greater than $ucr$ and insensitive to *B*_{0}. It is also clear from Fig. 9(a) that neither of the CX or CIV flow-sustained regimes exists below a critical magnetic field strength, irrespective of the flow velocity. Notably, a modest reduction in *B*_{0} from 25 mT to 10 mT can have profound consequences on the steady-state nature of the flow/plasma interaction for the conditions considered in Fig. 9(a). This result implies that the lower limit on magnetic field strength observed in CIV studies^{9} may be understood from energy confinement considerations without invoking collective heating; this is discussed in Sec. IV.

In Part I,^{16} we derived a general set of physical requirements for each flow regime and presented a map of these requirements as a function of the quantities $\tau \u0302e/\tau \u0302cap,\u2009\tau \u0302ie/\tau \u0302cap$, and $\u03f5\u0302ion$. The regime mapping in the limit $\tau \u0302ie\u226a\tau \u0302cap$ (which is valid across all simulation cases) is reproduced here in Fig. 9(b), onto which the data presented in Fig. 9(a) are projected. Figure 9 demonstrates clear agreement between the theoretical framework describing the steady-state regimes and the results of empirical data-supported global modeling. For large *B*_{0}, increasing $u\u221e$ within the I regime decreases both $\u03f5\u0302ion$ (leftward shift) and $T\u0302i$ (upward shift). Eventually the I-CX mode transition boundary is reached where the energy deposited into the plasma via stream charge exchange is large enough to outweigh ionization losses. The curves shift upward following the I-CX mode transition due to a combination of decreased $T\u0302i$ and increased $\tau \u0302e$. The slight rightward shift immediately after the transition can be explained by increased collisional excitation losses associated with a decrease in $T\u0302e$ [ Appendix A Eqs. (A16)–(A18)]. Further increasing $u\u221e$, the flow eventually transitions into the CIV regime where the energy deposited into the plasma via stream neutral ionization outweighs losses due to both ionization and thermal diffusion. The regime requirements mapping also explains why the $B0=0.01$ T case in Fig. 9(a) does not transition into a flow-sustained regime with increasing $u\u221e$. Specifically, since the diffusion coefficients scale inversely with *B*, the energy lost to particle diffusion is greater than the energy captured from the neutral flow for all values of $u\u221e$.

The CX-CIV threshold velocity can be derived analytically using the energy constraint for the CIV mode. Detailed in Part I^{16} and summarized here, the total energy balance for the system is observed to consist of the kinetic energy gained from ionizing a stream particle and the subsequent loss due to ionization cost and ultimate diffusion of the resulting charged particles. In the limit $\tau \u0302ie\u226a\tau \u0302cap,\u2009T\u0302e\u2248T\u0302i$, in which case the dimensional CIV mode energy balance can be written,

where we have defined $\u03f5tot$ as the sum of energy lost to ionization and thermal diffusion of ions and electrons. Because $\u03f5ion(Te)$ decreases monotonically with *T _{e}* [ Appendix A Eq. (A16)], $\u03f5tot(Te)$ has a definite minimum. Therefore, there exists a minimum value of $u\u221e$ to satisfy Eq. (36), below which the equality does not hold and the CIV regime cannot exist. This threshold can be written symbolically as

For the argon model used here, we find the minimum occurs around $Te\u22486$ eV and has a value $\u03f5tot\u224851$ eV. Using Eq. (37), the threshold velocity is calculated to be $u\u221e*\u224815.6$ km/s. This value, marked in Fig. 9(a), accurately predicts the CX-CIV transition observed from the global model results.

### D. Plasma/flow coupling

The global model approach thus far provides a detailed look at the net exchange of mass and energy between the flow and plasma. This total exchange is specifically important to plasma aerocapture technologies which propose to leverage the atmospheric flow for plasma sustainment. However, the control volume by definition only represents particles trapped by the dipole field. We seek to distinguish these so-called “absorbed” particles from those that may be ionized outside the control volume, which are not accounted for in the model but may still influence the coupling between the plasma and stream.

We examine flow absorption by summing the stream ionization and stream CX contributions to mass/energy in the global model equations. Figure 10 shows this for three different cases of $n\u221e$ exhibiting I-CX mode transitions. $Pcap$ and $m\u0307cap$ represent the rates of energy and mass absorption, respectively, from the stream by the plasma. The transition to the CX regime is apparent as an order-of-magnitude enhancement in flow absorption. Interestingly, the mode transition is observed here where *B*_{0} is the variable parameter instead of $u\u221e$. This is consistent with the scaling predicted in Part I;^{16} the characteristic electron diffusion time increases with magnetic field strength, eventually reaching a threshold where the diffusion time exceeds a critical value, $\tau \u0302e>\tau \u0302e*$, and massive intake of flow energy is possible. The shift of the mode transition to lower *B*_{0} at higher $n\u221e$ is also expected, since faster characteristic diffusion can be tolerated with the increase in captured energy. Consistent with our understanding of the I regime, we see that the flow mass absorbed before the mode transition is on the order of the injected mass, while after the transition the flow sustains the plasma with several orders more mass and many kW of power. This has implications for the potential of a dipole plasma to effectively utilize *in situ* atmospheric flow to fuel a plasma aerocapture maneuver.

Although absorption by the plasma is substantial in the CX and CIV regimes, we notice in Fig. 4 that much of the ionizing interaction occurs outside the control volume, where ions are not trapped. This produces significant deflected flow as these new ions are turned away by the field. A force is imparted by this deflection resulting from the interaction of the applied magnetic field with the diamagnetic current produced by newly formed ions.^{19} Recall that the magnetic field induced by this current does not significantly alter the applied field provided that the dipole magnetic pressure dominates freestream dynamic pressure.

A single-particle modeling approach was presented in Part I^{16} to derive the total drag *F _{D}* on the dipole as a function of $\zeta tot*$ and

*ρ*. This force is plotted in Fig. 11(a) for three cases of $n\u221e$. We observe a sharp jump in

_{L}*F*when the plasma transitions from the I to CX regime, reaching hundreds of N at moderate field strengths for $n\u221e=1018$ m

_{D}^{−3}. To contextualize these results, we compare the drag on a dipole plasma to that of a rigid body by introducing the normalized drag,

This is a measure of momentum transfer to the plasma compared to the force expected on a generic magnet-sized body under the same flow conditions. Figure 11(b) shows that before the CX regime, the force imparted on the plasma is weak ($F\u0302D<1$). Transitioning to the CX regime yields significant enhancement to the momentum transferred from flow to magnet. The drag can be tens of times larger than that on a rigid body immersed in this flow. It is notable that the general shape of Fig. 11(b) is consistent with that predicted theoretically in Part I.

## IV. DISCUSSION OF RESULTS

Our analysis reveals that critical ionization occurs as a result of global interaction between a neutral flow and dipole plasma, leading to an order-of-magnitude increase in flow absorption by the plasma. This analysis differs from most CIV models which typically use one spatial dimension and are limited to describing the plasma/flow interaction locally (or assuming infinite extent). The CIV phenomenon is usually described for a plasma flowing relative to a stationary neutral gas and magnetic field (a notable exception is the model of McKenzie and Varma^{29} although they explicitly omit consideration of the mechanism for energy transfer to the electrons). Prevailing theory states that the counterstreaming plasma induces turbulent heating of the electrons that drives ionization of the gas.^{13} In this study, we instead describe a confined plasma, considering only collisional energy transfer between species. Furthermore, the plasma is stationary with respect to the magnetic field rather than impinging as a beam. Our observation of critical ionization even after changing the classical assumptions adds to the wealth of evidence in literature for the ubiquity of the CIV phenomenon.

This model supports the effects anticipated by the theory of Part I.^{16} We observe that changing the velocity, density, and magnetic field strength has significant impact on the mode of interaction as mapped in Fig. 9(b). These results are a direct consequence of including diffusion as a primary loss mechanism. For example, the I-CX threshold velocity is reduced as *B*_{0} is increased because electron energy is confined longer and therefore less flow heating of electrons is required to sustain the plasma. The effect of diffusion on critical ionization in a closed field geometry is profound, but it remains uncharacterized in the literature. One exception is an early analysis of homopolar device experiments by Lehnert,^{30} but this uses parallel diffusive loss to argue that thermal electron-impact ionization could not account for observed CIV effects. As we have shown, this argument does not necessarily hold when plasma is confined along the field lines and instead diffuses perpendicularly. Our study is the first comprehensive investigation of critical ionization physics for such a closed-field geometry where perpendicular diffusion dominates the particle and energy balance.

Möbius *et al.*^{15} previously discussed the effects of various energy transfer processes on the interaction of a plasma beam with a neutral gas. Notably, a threshold velocity was derived that implied critical ionization could initiate below $ucr$ provided that the CX reaction rate dominates ionization. They concluded that CX can thus be a “triggering process” for critical ionization because it adds energy to the plasma beam without an associated depletion of tail electron energy. McNeil *et al.*^{14} similarly identified the importance of CX to enabling critical ionization when electron energy is scarce. The CX regime presented here is distinct from these observations, however, because it is not an enhancement to CIV; it exhibits characteristically different physics. Specifically, in the CX regime, the neutral reservoir feeding critical ionization is not the flow but rather the secondary neutrals arising from the ion population. Charge exchange thus acts as the primary source of both neutrals and the energy to ionize them, a unique distinction from present theory.

No well-characterized experimental observation of this CX regime transition exists in literature. This is likely due to open field geometries that do not confine electron energy well. In traditional CIV experiments, heating of the electrons to ∼100 eV is observed,^{31} which is necessary to achieve ionization times faster than the rapid loss time scale. At these energies, flow ionization is so efficient that the CX regime would collapse and the observed transition would be from I to CIV, as indicated by Fig. 9(b) of Part I.^{16} However, some anomalous experimental results are potentially explained by collisional effects. For example, Himmel *et al.*^{32} observed three regimes in their rotating plasma investigation of critical ionization, with threshold velocities below $ucr$ measured at very high background pressures. They made no explicit connection to collisional heating, but Machida and Goertz^{11} briefly theorized that the observation could be explained by an electron energy balance between resistive heating and effective ionization losses, which is consistent with our results. However, we note that the neutral densities in these experiments ($n\u221e>1022$ m^{−3}) are well outside the free molecular neutral flow regime assumed by our model ($n\u221e\u22721019$ m^{−3}). Venkataramani and Mattoo^{33} also found that plasma deceleration through a neutral gas cloud was far more rapid than could be explained by ionization alone. The authors postulated that resonant charge exchange could explain the enhanced momentum loss at higher neutral densities. This seems to agree with Machida and Goertz' simulation study,^{11} which predicted that resistive processes should dominate collective ones at high neutral densities. Indeed, in our study, we observe starker critical ionization thresholds to the CX regime as density is increased (see Fig. 10, for example). The species, velocities, densities, and field strengths in Venkataramani's experiment are remarkably similar to the parameters used in our study. Their conclusions may serve as a useful backdrop for future laboratory investigations of the CX regime transition, but ultimately novel experimentation is needed to clarify and validate the physical regimes identified here.

Brenning and Axnäs^{34} have discussed some of the “unsolved problems” in CIV theory. One in particular is the effect of magnetic field strength on the interaction. They indicated there is evidence the flow must be subalfvénic ($u\u221e<vA$) in order for turbulent heating not to be damped, imposing a lower bound on the magnetic field strength depending on neutral density. We have shown a lower bound on magnetic field arising even when only collisional heating is allowed. This suggests the link between critical ionization and field strength may be nonspecific to turbulent heating while still being consistent with experimental observations. Similarly, they described two distinct regimes within the collective-heating models that depend on neutral density. We show here this distinction exists with collisional heating as well, suggesting an underlying mechanism agnostic to the energy transfer channel.

Future work is required to expand the physical fidelity and complexity of this model. The density profile in Eq. (21), which is supported by theory,^{21,35} was primarily adopted to simplify the model formulation. However, it may be modified to more closely represent realistic dipole plasmas. Pressure profiles are known^{20} to depend on plasma beta and, in the laboratory, on anisotropy, finite coil effects, and plasma currents,^{25} which are all neglected in this study. Even the presence of a mechanical support, which will be required in an experimental validation of this study, modifies the density scaling with *ψ* considerably.^{25} An additional simplification is the assumption that secondary neutrals share the same density profile as the plasma. In a real system, they are more likely to exist at the interaction front since they result from CX with the stream. Figure 4 demonstrates that the plasma density is a reasonable geometric proxy for the interaction region and therefore secondary neutral production. Future higher fidelity simulations must examine the steady state $n\u03022n$ accounting for spatially accurate neutral production and depletion. We reiterate that the general results found here are insensitive to a change in the assumed density profile, as evidenced by similar previous results using a different density scheme.^{6,7}

Finally, we note the importance of these results for application to plasma aerocapture technology. Figure 10 indicates that in order to sustain a critical ionization discharge in the dipole, tens of kW of plasma heating is required. Artificially supporting this discharge from a spacecraft-based power source would negate the mass savings that motivate the use of aerocapture. Therefore, these results provide bounds on the magnet design and atmospheric flight trajectory for a mission using plasma aerocapture. The demands of aerocapture can be extreme; drag forces on the order of 100 kN are needed to decelerate massive spacecraft by several km/s in a few minutes. A systems analysis of plasma aerocapture^{7} shows that at the outer planets, without excessive onboard fuel and power, this may only be achieved by utilizing the flow-sustained regimes (CX and CIV). Additionally, the variation of drag force with magnetic field strength shown in Fig. 11 implies the capability for continuously variable drag modulation, an aerocapture flight control scheme that dramatically improves accuracy and reduces risk.^{36} This would be a significant technological leap for aerocapture concepts as modern proposed devices have no mechanism for this type of drag modulation. Further development is required to represent atmospheric species of interest for aerocapture missions, such as $CO2,H2$, or $N2$, involving considerably more conservation equations of greater intricacy.

## V. CONCLUSION

This study details an analytical model of the global interaction between a rarefied, hypersonic neutral flow and magnetic dipole plasma. By simulating over a large parameter range, the global model reveals three distinct physical regimes governed by injection (I), charge exchange (CX), and ionization (CIV). These regimes, discussed in Sec. III B, are equivalent to those derived theoretically in Part I.^{16} While the theory treated various time scales and reaction rates as independent parameters, this analysis demonstrates the existence of the two critical ionization thresholds (I-CX and CX-CIV) even when accounting for nonlinearly coupled physics. Notably, this global model has not invoked plasma instabilities or other collective electron heating processes; instead, we demonstrate that critical ionization is not only possible with but enhanced by charge exchange and Coulomb collisions. As anticipated in Part I, the threshold velocity for the I-CX regime transition can be considerably lower than $ucr$ while the CX-CIV transition is necessarily higher. These thresholds induce a significant enhancement of the global interaction between the plasma and stream. Substantial flow mass and energy are deposited to naturally sustain the plasma discharge while large forces are coupled to the magnetic field.

The use of a global model in the investigation of a hypersonic flow interacting with magnetized plasma offers certain advantages. Numerical studies of confined plasmas are usually computationally intensive and thus must be constrained to very specific parameters for their application. In the model derived here, geometric assumptions allow us to numerically solve a large system in a few seconds. This enables its use for a broad array of simulation conditions as seen in this study. The species interactions are input from experimental data or functional approximations thereof, allowing for the addition of far more species and reactions so long as the limiting assumptions are not violated. These features of the model make it applicable to many types of investigations, such as validation of experimentally observed plasma/neutral interactions, simulation of space plasmas, or performance scaling of plasma aerocapture and propulsion concepts.

## ACKNOWLEDGMENTS

This work was supported by a NASA Space Technology Research Fellowship.

## DATA AVAILABILITY

The data that support the findings of this study are available from the corresponding author upon reasonable request.

### APPENDIX A: PLASMA MODEL

Here we detail the implementation of the plasma model in this study. All units are SI (with temperatures in J) unless otherwise specified.

##### A. Coulomb interaction

Braginskii^{22} defines the ion–electron collision time as

after converting from cgs units to SI, where $ln\u2009\Lambda $ is the Coulomb logarithm. The normalized form of this appearing in Eqs. (29) and (30) is

##### B. Plasma diffusion

The Bohm diffusion coefficient *D _{B}* is given by Lieberman and Lichtenberg

^{23}as

where *e* is the elementary charge. For electron diffusion due to collisions with heavy particles, the diffusion coefficient *D _{c}* is given by Goldston and Rutherford

^{24}as

where *ν _{ei}* and

*ν*are the electron–ion and electron–neutral collision frequencies, respectively. The electron–ion collision frequency is given by Goldston and Rutherford,

_{en}^{24}Eq. (11.22), as a function of

*n*and

_{i}*T*,

_{e}where $\epsilon 0$ is the vacuum permittivity. The electron–neutral collision frequency is modeled by Itikawa^{37} as

where *σ _{en}* is the average electron–neutral momentum transfer cross section. Energy-dependent values of

*σ*are sourced from Harstad.

_{en}^{38}The normalized forms of these diffusion coefficients appearing in Eqs. (26), (29), and (30) are

where $\nu \xafej=\u27e8\nu ej\u27e9/nj$ is defined to remove the spatial dependence of *D _{c}*. Note that $D\u0302c$ depends on model outputs $N\u0302i$ and $N\u03022n$ and thus changes with time.

##### C. Argon Model

Empirical data on reaction rates and cross sections for argon are compiled from various sources. We identify these here.

###### 1. Charge exchange

The charge exchange reaction rate is modeled after Meier and Shumlak.^{17} After integrating over the distribution function to find $\u27e8\sigma v\u27e9$, they obtain

where *n _{n}* is the neutral species density and $Vcx$ is an effective collision velocity. In our model there are two neutral species. The two effective velocities are

for thermal velocities $vT,j$. The CX cross section, $\sigma cx(Vcx)$ is sourced from Phelps^{39} using the effective velocity to calculate the collision energy. The heat transfer associated with CX between species *j* and *k* is described by Meier and Shumlak^{17} with the parameter $Qj,k$ that depends on thermal and fluid velocities. After normalization, the heat transfer term appearing in Eqs. (29) and (31) is

where *σ _{cx}* is that for

*j*–

*k*charge exchange.

###### 2. Ionization

The ionization reaction rate is given analytically by Table 3.3 of Lieberman and Lichtenberg^{23} as

for *T _{e}* in eV. The effective ionization cost, $\u03f5ion$, is also given in eV units by Lieberman and Lichtenberg in Eq. (3.5.8) as

for *T _{e}* in eV. Here, $Uex=12.14$ eV is the first excitation energy of argon and $Uiz=15.76$ eV is the first ionization energy. $Rel$ and $Rex$ are the elastic scattering and first excitation reaction rates, respectively, and are also given in Table 3.3 as

for *T _{e}* in eV.

###### 3. Neutral diffusion

Lieberman and Lichtenberg^{23} give a general form for gas diffusion, which we adapt for two collision types, 2*n*–2*n* gas-kinetic and *i*-2*n* momentum transfer. This yields the secondary neutral diffusion coefficient,

where $\sigma Qm$ and *σ _{gk}* are the cross sections of ion–neutral momentum transfer and neutral–neutral collisions, respectively. The Ar–Ar

^{+}momentum transfer cross section, $\sigma Qm$, is sourced from Phelps.

^{39}We use the highest tabulated value, $1.57\xd710\u221218$ m

^{2}, to represent a worst-case diffusion of secondary neutrals, noting that this value changes by less than a factor of two over the energy range of interest for this study. The gas-kinetic cross section, $\sigma gk=5\xd710\u221219$ m

^{2}, is given by Lieberman and Lichtenberg

^{23}in Table 9.3, reproduced from Smirnov.

^{40}The normalized diffusion coefficient appearing in Eqs. (27) and (31) is

Note again that $D\u03022n$ is dependent on model outputs $N\u0302i$ and $N\u03022n$ and thus varies with time.

### APPENDIX B: CONTROL VOLUME INTEGRALS

The volume-averaging process results in several integrals appearing in the global model equations (26)–(27) and (29)–(31). We note that with the exception of *I _{sn}*, these integrals are only functions of the magnetic field profile and the CV size. Both these parameters are static for a given simulation, so these integrals need only be computed once per full numerical solution. They also lend themselves well to functional fitting. These characteristics enable drastic reduction of the computational intensity associated with solving the model.

The integral that relates the species population $N\u0302j$ and reference density $n\u0302j,r$ in relation (24) is defined as

reflecting the density profile in Eq. (21). Similarly, when volume-averaging interaction terms between tracked species, the density profile is squared and integrated, yielding

Volume-averaging interactions of tracked species with the neutral stream results in the capture integral,

Computing this integral requires knowledge of $n\u0302sn$. We derive $n\u0302sn$ by recognizing that the only stream–plasma interactions are ionization and CX and assuming the interaction is steady-state. Thus, the stream continuity equation has the form

where $Rtot=Riz+Rcx,sn$ is the total stream reaction rate. Equation (B4) is normalized and solved to yield

The stream capture time scale, $\tau \u0302cap$, is defined as in Eq. (25) for reaction rate $Rtot$. Therefore, $n\u0302e,r/\tau \u0302cap$ represents the characteristic rate at which stream neutrals interact with the plasma. *I _{sn}* must be computed at each time step of the numerical solution, unlike the other CV integrals. A functional fit for $Isn(\rho L,n\u0302e,r/\tau \u0302cap)$ is provided in Part I

^{16}that facilitates rapid computation with reasonable accuracy.

Volume averaging of diffusion terms results in magnetic field integrals that reflect the spatial scalings associated with *D _{B}*,

*D*, and $D2n$. The respective integrals are

_{c}Note these are surface integrals over $S*$ owing to Gauss' theorem. Once again, *I _{B}*,

*I*, and $I2n$ are all functions of the magnetic field only. It is the different scaling with

_{c}*B*in each

*D*that yields three unique integrals.

_{j}## References

*Reviews of Plasma Physics*(Consultants Bureau, New York, 1965), Vol. 1, p. 205.