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.

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 allowed14,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.

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.

FIG. 1.

The global model is formulated in the zr plane with symmetry about θ. Level curves of ψ [Eq. (5)] are shown in gray with the current loop at (0,1). The plasma density (pink) scales with ψα for some constant α. The control volume V* is bounded by the ψ=ψ* surface, S*. Axially flowing neutrals have freestream density n and velocity u.

FIG. 1.

The global model is formulated in the zr plane with symmetry about θ. Level curves of ψ [Eq. (5)] are shown in gray with the current loop at (0,1). The plasma density (pink) scales with ψα for some constant α. The control volume V* is bounded by the ψ=ψ* surface, S*. Axially flowing neutrals have freestream density n and velocity u.

Close modal

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 (2n), i.e., those generated in i-sn charge exchange (CX) reactions. Conservation of mass and energy for each population j is described by

njt+·(njuj)=sṅjs,
(1)
εjt+·(εjuj)=sε̇js,
(2)

where nj is the number density, εj is the energy density, uj is the bulk fluid velocity, and /t 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 ẋjs is a rate of change of variable x in species j due to process s. The energy density is defined following Meier and Shumlak17 as

εj=12mjnjuj2+pj/(γ1),
(3)

where mj is the particle mass, pj is the species pressure, and γ 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 rc in cylindrical coordinates (r,θ,z). This is described by Jackson18 as

B=1rêθ×ψ,
(4)
ψ=2B0rc2rπ(2κ2)K(κ2)2E(κ2)κ2(rc2+r2)+z2,
(5)
κ2=4rrc(r+rc)2+z2,
(6)

where êθ 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. B0 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:

  1. 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 εsn is, by definition, fixed, and (3) assuming quasineutrality and singly charged ions (ni = ne). The latter argument equating ion and electron continuity implies that plasma diffusion is ambipolar.

  2. We assume a dense (ne>1016 m−3), collisional plasma with electron Hall parameter Ωe1. 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).

  3. 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 

  4. By neglecting parallel currents20,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 

FIG. 2.

Diagrammatic representation of the global model. Ionization of neutral species j (iz,j), charge exchange with neutral species j (cx,j), Coulomb collisions (coll), artificial plasma injection (inj), and diffusion (dif) are the processes contributing to mass and energy transfer between the populations sn, i, e, and 2n. Some of these processes involve both mass and energy transport (shaded arrows) while others only transfer energy (unshaded arrows).

FIG. 2.

Diagrammatic representation of the global model. Ionization of neutral species j (iz,j), charge exchange with neutral species j (cx,j), Coulomb collisions (coll), artificial plasma injection (inj), and diffusion (dif) are the processes contributing to mass and energy transfer between the populations sn, i, e, and 2n. Some of these processes involve both mass and energy transport (shaded arrows) while others only transfer energy (unshaded arrows).

Close modal
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

ṅiiz,j=njneRiz,
(7)

denoting the rate coefficient of a reaction s as Rs. Neutral energy is absorbed by the ion population while electrons incur a loss of ϵion, which is the effective ionization cost including losses due to scattering and excitation,

ε̇iiz,j=32Tjṅiiz,j,
(8)
ε̇eiz,j=ṅiiz,j(32Me/jTjϵion),
(9)

where Mk/jmk/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,

ṅ2ncx,sn=ninsnRcx.
(10)

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

ε̇icx=Qi,2ncx,2nQ2n,icx+12miu2ṅ2ncx,sn,
(11)
ε̇2ncx=Qi,2ncx,2n+Q2n,icx.
(12)

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 Braginskii22 who expresses the transfer rate as

ε̇icoll=ε̇ecoll=neτie(TeTi),
(13)

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 Γinj at fixed temperatures, Ti,inj and Te,inj, such that the particle and energy rates are

ṄiinjΓinjmi,
(14)
Ėjinj=52Tj,injṄiinj,
(15)

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., εjpj/(γ1). Considering only perpendicular diffusion, we have23 

njuj=Djnj,
(16)
εjuj=Tjγ1Djnj,
(17)

where is the field line-normal gradient operator and Tj 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 Dj 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, DB, and collisional coefficient, Dc, are given by Lieberman and Lichtenberg23 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 Lichtenberg23 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 ) are primarily used for normalization. The following scheme is employed:

(ẑ,r̂)=(z,r)/rct̂=tu/rcn̂j=nj/nûj=uj/uε̂j=εj/εT̂j=Tj/TB̂=B/B0ψ̂=ψ/(B0rc2).
(18)

The stream is hypersonic with freestream velocity uêz and freestream density n (Fig. 1). Thus, the freestream energy density is primarily kinetic, i.e., εmsnnu2/2. We define an effective temperature of the incoming flow, T, such that ε=3nT/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.

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 is confined. This boundary has been previously derived6 (with more detail in Part I16) by an analysis of ion trajectories in the dipole field. It is a toroidal volume enclosed by a surface of constant ψ̂,

ψ̂*=2ρL,
(19)

where ψ̂* is the normalized value of the magnetic flux function along the boundary contour and

ρL=miueB0rc,
(20)

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.,

nj=nj,r(ψ̂ψ̂r)α,
(21)

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 ni, ne, 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 (njnk) do not result in higher order terms. While these assumptions sacrifice physical fidelity, they simplify the model such that spatial dependencies are functions of ψ 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̂j, and the normalized temperatures, T̂j, of each tracked species as a function of time, t̂. We derive these by normalizing and integrating the conservation equation variables, nj and εj,

N̂j=V*n̂jdV̂,
(22)
Êj=V*ε̂jdV̂=2/3γ1N̂jT̂j.
(23)

Here, Êj is the normalized total energy of species j in the CV. Combining Eqs. (21) and (22) gives the useful relation

n̂j,r=N̂j/Iψα,
(24)

where Iψα is a magnetic field integral defined in  Appendix B Eq. (B1). Relation (24) allows us to describe the spatially resolved density profile, n̂j, from the CV-averaged solution, N̂j, through Eq. (21). Recall that V* (defined by ψ̂>ψ̂*) is a function of ρL only and ψ̂r is fixed.

Interactions between tracked species, of the form ṅjsnjnkRs, depend spatially on the species densities which both follow Eq. (21). Such interactions include 2n ionization and i-2n charge exchange. After volume averaging, these terms have the form Ṅ̂jsn̂j,rn̂k,rIψ2α/τ̂s, where Iψ2α [defined in  Appendix B Eq. (B2)] is a magnetic field integral resulting from taking the product of tracked species densities. τ̂s is a dimensionless parameter describing the characteristic time scale of reaction s compared to its transit time across the magnet, defined as

τ̂s=u/rcRsn.
(25)

Interactions of tracked species with the neutral stream, having the form ṅjsnjnsnRs, are more complicated because nsn does not follow Eq. (21). The resultant volume-averaged interaction has the form Ṅ̂jsn̂e,rIsn/τ̂s, where the so-called “capture integral” Isn involves integrating n̂sn. Its derivation is found in  Appendix B, including an analytic solution to the stream neutral continuity equation yielding n̂sn.

Finally, tracked species diffusion, with the form ṅjdif·(Djnj), depends on the different diffusion models described in Sec. II A 3. Volume averaging results in a common scaling for each, Ṅ̂jdifD̂jn̂j,rID. Here, D̂j is a dimensionless, spatially invariant form of Dj specific to each diffusion model. Similarly, ID are magnetic field integrals that reflect the spatial scalings associated with DB, Dc, D2n, and density gradients. The different scaling with B in each Dj yields three unique integrals IB, Ic, and I2n that are all functions of the magnetic field only [ Appendix B Eqs. (B6)–(B8)].

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̂i,N̂2n,T̂i,T̂e, and T̂2n. 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̂j,r is replaced with N̂j/Iψα 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 ni = ne, 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,

dN̂idt̂=[D̂BT̂eIB+D̂c(T̂e+T̂i)Ic]N̂iIψαdif+Ṅ̂injinj+N̂eN̂2nτ̂izIψ2α(Iψα)2iz,2n+N̂eτ̂izIsnIψαiz,sn,
(26)
dN̂2ndt̂=N̂iτ̂cx,snIsnIψαcx,snD̂2nT̂2n1/2N̂2nI2ndifN̂eN̂2nτ̂izIψ2α(Iψα)2iz,2n.
(27)

Diffusion is a loss for both species, with D̂B,D̂c, and D̂2n 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,

Ṅ̂inj=Ṅiinjnurc2,
(28)

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, Isn, 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 2ns that are reionized. Both of these stream interaction processes are functions of Isn and supply ions from outside the CV with no corresponding loss term from another model equation.

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,

2/3γ1d(N̂iT̂i)dt̂=2/3γ1[D̂BT̂eIB+D̂c(T̂e+T̂i)Ic]N̂iT̂iIψαdif+(1τ̂cx,snQ̂sn,i)N̂iIsnIψαcx,sn+N̂eτ̂izIsnIψαiz,sn+P̂i,injinj+32(Q̂i,2nQ̂2n,i)N̂iN̂2nIψ2α(Iψα)2cx,2n+N̂eN̂2nT̂2nτ̂izIψ2α(Iψα)2iz,2n+2Me/iN̂e2τ̂ieT̂eT̂iT̂e3/2Iψ2α(Iψα)2coll,
(29)
2/3γ1d(N̂eT̂e)dt̂=(Me/snϵ̂ion)N̂eτ̂izIsnIψαiz,sn2Me/iN̂e2τ̂ieT̂eT̂iT̂e3/2Iψ2α(Iψα)2coll+(Me/2nT̂2nϵ̂ion)N̂eN̂2nτ̂izIψ2α(Iψα)2iz,2n+P̂e,injinj2/3γ1[D̂BT̂eIB+D̂c(T̂e+T̂i)Ic]N̂eT̂eIψαdif,
(30)
2/3γ1d(N̂2nT̂2n)dt̂=32(Q̂2n,iQ̂i,2n)N̂iN̂2nIψ2α(Iψα)2cx,2n+N̂iQ̂sn,iIsnIψαcx,snN̂eN̂2nτ̂izT̂2nIψ2α(Iψα)2iz,2n2/3γ1D̂2nT̂2n3/2N̂2nI2ndif.
(31)

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

P̂j,inj=53T̂j,injṄ̂inj,
(32)

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̂j,k and τ̂ie defined in  Appendix A Eqs. (A12) and (A2), respectively]. Although ionization terms are similar to their continuity counterparts for i and 2n, the electrons incur an additional loss of

ϵ̂ion=ϵionmsnu2/2,
(33)

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

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 solver27 is used and outputs the time evolution of N̂i,N̂2n,T̂i,T̂e, and T̂2n. 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,u, rc, B0, and Γ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 ψ̂r=ψ̂* 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 njB 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.

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=1018 m−3, u=12 km/s, B0=0.25 T, rc = 1 m, and Γinj=5 mg/s. At these conditions, t̂=12 corresponds to t =1 ms. The initial conditions are N̂i(0)=N̂2n(0)=1 and T̂i(0)=T̂2n(0)=T̂e(0)=0.1. The populations N̂i and N̂2n 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̂240) although this time decreases if the initial conditions of N̂i and N̂2n are increased. The steady-state populations and temperatures, however, are only sensitive to the scaling parameters (n, 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̂1). This justifies our earlier assumption that density distributions can be described by their equilibrium profiles.

FIG. 3.

Time-resolved numerical solutions to the model outputs: (a) populations N̂i and N̂2n and (b) temperatures T̂i,T̂e, and T̂2n. Simulation conditions are n=1018 m−3, u=12 km/s, B0=0.25 T, rc = 1 m, and Γinj=5 mg/s.

FIG. 3.

Time-resolved numerical solutions to the model outputs: (a) populations N̂i and N̂2n and (b) temperatures T̂i,T̂e, and T̂2n. Simulation conditions are n=1018 m−3, u=12 km/s, B0=0.25 T, rc = 1 m, and Γinj=5 mg/s.

Close modal

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.

FIG. 4.

Spatial distribution of plasma and stream densities as the system evolves to steady state. Three times are sampled: (a) right before the sudden density increase depicted in Fig. 3; (b) right after the sudden density increase; (c) steady state. The control volume boundary is shown as a solid line. Simulation conditions are n=1018 m−3, u=12 km/s, B0=0.25 T, rc = 1 m, and Γinj=5 mg/s.

FIG. 4.

Spatial distribution of plasma and stream densities as the system evolves to steady state. Three times are sampled: (a) right before the sudden density increase depicted in Fig. 3; (b) right after the sudden density increase; (c) steady state. The control volume boundary is shown as a solid line. Simulation conditions are n=1018 m−3, u=12 km/s, B0=0.25 T, rc = 1 m, and Γinj=5 mg/s.

Close modal

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̂<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̂=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.

FIG. 5.

Temporal solution of the contributing terms in the control volume equations. Plots correspond to the dimensionless model equations as follows: (a) Eq. (26), (b) Eq. (27), (c) Eq. (29), (d) Eq. (31), and (e) Eq. (30). The simulation conditions are n=1018 m−3, u=12 km/s, B0=0.25 T, rc = 1 m, and Γinj=5 mg/s.

FIG. 5.

Temporal solution of the contributing terms in the control volume equations. Plots correspond to the dimensionless model equations as follows: (a) Eq. (26), (b) Eq. (27), (c) Eq. (29), (d) Eq. (31), and (e) Eq. (30). The simulation conditions are n=1018 m−3, u=12 km/s, B0=0.25 T, rc = 1 m, and Γinj=5 mg/s.

Close modal

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̂* (i.e., on ψ̂*) and average temperatures T̂ on the freestream velocity u. In contrast to the three tracked species, n̂sn is not invariant along field lines, so we compute n̂sn* 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's8 critical velocity,

ucr=2Uizmsn,
(34)

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̂sn*1) and few ions are present. The model species are thermalized with the stream (T̂1) by CX interaction. As the velocity increases past u0.4ucr, there is a sudden exponential increase in plasma and secondary neutral density, while n̂sn* 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̂2n* roughly 10% of n̂i*. Increasing velocity past u1.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.

FIG. 6.

(a) Reference densities, n̂*, and (b) average temperatures, T̂, of all system species. n̂i,e* and n̂2n* are invariant on the reference surface S* while n̂sn* is computed at the downstream-most point of the reference surface. The fixed conditions are n=1018 m−3, B0=0.25 T, rc = 1 m, and Γinj=5 mg/s.

FIG. 6.

(a) Reference densities, n̂*, and (b) average temperatures, T̂, of all system species. n̂i,e* and n̂2n* are invariant on the reference surface S* while n̂sn* is computed at the downstream-most point of the reference surface. The fixed conditions are n=1018 m−3, B0=0.25 T, rc = 1 m, and Γinj=5 mg/s.

Close modal

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. 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.

FIG. 7.

Steady state of significant contributing terms in the control volume equations as a function of stream velocity. Plots correspond to the dimensionless model equations as follows: (a) Eq. (26), (b) Eq. (27), (c) Eq. (29), (d) Eq. (31), and (e) Eq. (30). Black dotted lines distinguish the three physical regimes I (left), CX (middle), and CIV (right). The fixed conditions are n=1018 m−3, B0=0.25 T, rc = 1 m, and Γinj=5 mg/s.

FIG. 7.

Steady state of significant contributing terms in the control volume equations as a function of stream velocity. Plots correspond to the dimensionless model equations as follows: (a) Eq. (26), (b) Eq. (27), (c) Eq. (29), (d) Eq. (31), and (e) Eq. (30). Black dotted lines distinguish the three physical regimes I (left), CX (middle), and CIV (right). The fixed conditions are n=1018 m−3, B0=0.25 T, rc = 1 m, and Γinj=5 mg/s.

Close modal

Injection (I) regime: At low velocities, the only energy deposited into the plasma comes from the power injection term P̂inj [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 2n 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 2ns, is transferred to the electron population via Coulomb interaction. Electrons subsequently have enough energy to ionize some of the 2n 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 2ns 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 ϵ̂ion (i.e., ionization and radiation) and 2n 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̂i* now decreases with u due to a combination of increased Bohm diffusion as T̂e 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̂e,r/τ̂cap is the critical term that couples the plasma properties to the flow through Isn [Eq. (B3)]. Whereas earlier this parameter was defined for an arbitrary contour ψr, here we define it on the CV boundary, ψr=ψ*, to represent the lower bound of reaction rate in the CV. We examine the behavior of ζs* in the three regimes where

ζs*=Rsne*u/rc=n̂e*/τ̂s,
(35)

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, ζs*1 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 ζtot*=ζcx,sn*+ζ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 ζtot* in the I regime. The CX regime is characterized by a sharp and monotonic increase in ζcx,sn*, which drives ζtot* high enough to yield significant interaction between the plasma and flow. ζiz* remains quite low until a drastic jump denoting the CIV regime. Here, stream ionization is so strong it drives ζ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.

FIG. 8.

The reaction rates, ζs*, reveal three distinct physical regimes. Plasma-flow interaction strengthens with increasing injection power (Ėeinj=60 W–1.2 kW) in the I regime. The CX regime is dominated by charge exchange reactions while the CIV regime exhibits strong ionization of the stream. The reference surface is ψr=ψ* and the fixed conditions are n=1018 m−3, B0=0.25 T, rc = 1 m, and Γinj=5 mg/s (corresponding to the 0.3 kW curve).

FIG. 8.

The reaction rates, ζs*, reveal three distinct physical regimes. Plasma-flow interaction strengthens with increasing injection power (Ėeinj=60 W–1.2 kW) in the I regime. The CX regime is dominated by charge exchange reactions while the CIV regime exhibits strong ionization of the stream. The reference surface is ψr=ψ* and the fixed conditions are n=1018 m−3, B0=0.25 T, rc = 1 m, and Γinj=5 mg/s (corresponding to the 0.3 kW curve).

Close modal

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 (τ̂e), particle capture time scale (τ̂cap), ion–electron energy transfer time scale (τ̂ie), and the ratio of the effective ionization energy cost to the kinetic energy of stream neutrals (ϵ̂ion). 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̂i as a function of u for a range of magnetic field strengths. For four cases of B025 mT, the plasma undergoes both I-CX and CX-CIV mode transitions as u is increased. Interestingly, the I-CX threshold velocity decreases as B0 increases, exhibiting values both less than and greater than ucr. Unlike the I-CX mode transition, the threshold value of u for the CX-CIV mode transition is always greater than ucr and insensitive to B0. 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 B0 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 studies9 may be understood from energy confinement considerations without invoking collective heating; this is discussed in Sec. IV.

FIG. 9.

(a) Ion population, N̂i, as a function of neutral flow velocity, u, for magnetic field strengths from B0=0.01 to 1 T. The fixed conditions are n=1018 m−3, rc = 1 m, and Γinj=3 mg/s. Data points indicate: (○) I regime; (●) CX regime; and (●) CIV regime. Dashed vertical lines denote ucr [Eq. (34)] and the CIV threshold velocity, u* [Eq. (37)]. Results from panel (a) are plotted on the flow regime map in panel (b). Solid black lines mark the flow regime boundaries derived in Part I.16 

FIG. 9.

(a) Ion population, N̂i, as a function of neutral flow velocity, u, for magnetic field strengths from B0=0.01 to 1 T. The fixed conditions are n=1018 m−3, rc = 1 m, and Γinj=3 mg/s. Data points indicate: (○) I regime; (●) CX regime; and (●) CIV regime. Dashed vertical lines denote ucr [Eq. (34)] and the CIV threshold velocity, u* [Eq. (37)]. Results from panel (a) are plotted on the flow regime map in panel (b). Solid black lines mark the flow regime boundaries derived in Part I.16 

Close modal

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 τ̂e/τ̂cap,τ̂ie/τ̂cap, and ϵ̂ion. The regime mapping in the limit τ̂ieτ̂cap (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 B0, increasing u within the I regime decreases both ϵ̂ion (leftward shift) and T̂i (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̂i and increased τ̂e. The slight rightward shift immediately after the transition can be explained by increased collisional excitation losses associated with a decrease in T̂eAppendix A Eqs. (A16)–(A18)]. Further increasing u, 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. 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.

The CX-CIV threshold velocity can be derived analytically using the energy constraint for the CIV mode. Detailed in Part I16 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 τ̂ieτ̂cap,T̂eT̂i, in which case the dimensional CIV mode energy balance can be written,

12msnu2=ϵion(Te)+3Teϵtot(Te),
(36)

where we have defined ϵtot as the sum of energy lost to ionization and thermal diffusion of ions and electrons. Because ϵion(Te) decreases monotonically with TeAppendix A Eq. (A16)], ϵtot(Te) has a definite minimum. Therefore, there exists a minimum value of u to satisfy Eq. (36), below which the equality does not hold and the CIV regime cannot exist. This threshold can be written symbolically as

u*=2min{ϵtot(Te)}msn.
(37)

For the argon model used here, we find the minimum occurs around Te6 eV and has a value ϵtot51 eV. Using Eq. (37), the threshold velocity is calculated to be u*15.6 km/s. This value, marked in Fig. 9(a), accurately predicts the CX-CIV transition observed from the global model results.

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 exhibiting I-CX mode transitions. Pcap and ṁcap 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 B0 is the variable parameter instead of u. 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, τ̂e>τ̂e*, and massive intake of flow energy is possible. The shift of the mode transition to lower B0 at higher n 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.

FIG. 10.

Rate of energy (Pcap) and flow mass (ṁcap) absorption by the plasma. The fixed conditions are u=12 km/s, rc = 1 m, Γinj=5 mg/s, and Pinj=300 W. The model solution is intractable at n=1016 m−3 for B0<25 mT, implying the dipole is unsustainable below this field strength for the low density case.

FIG. 10.

Rate of energy (Pcap) and flow mass (ṁcap) absorption by the plasma. The fixed conditions are u=12 km/s, rc = 1 m, Γinj=5 mg/s, and Pinj=300 W. The model solution is intractable at n=1016 m−3 for B0<25 mT, implying the dipole is unsustainable below this field strength for the low density case.

Close modal

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 I16 to derive the total drag FD on the dipole as a function of ζtot* and ρL. This force is plotted in Fig. 11(a) for three cases of n. We observe a sharp jump in FD when the plasma transitions from the I to CX regime, reaching hundreds of N at moderate field strengths for n=1018 m−3. To contextualize these results, we compare the drag on a dipole plasma to that of a rigid body by introducing the normalized drag,

F̂D=2FDmsnnu2πrc2.
(38)
FIG. 11.

(a) Drag force, FD, and (b) normalized drag force, F̂D, on the magnet are modulated by the applied field strength, B0. The fixed conditions are u=12 km/s, rc = 1 m, and Γinj=5 mg/s. The model solution is intractable at n=1016 m−3 for B0<25 mT, implying the dipole is unsustainable below this field strength for the low density case.

FIG. 11.

(a) Drag force, FD, and (b) normalized drag force, F̂D, on the magnet are modulated by the applied field strength, B0. The fixed conditions are u=12 km/s, rc = 1 m, and Γinj=5 mg/s. The model solution is intractable at n=1016 m−3 for B0<25 mT, implying the dipole is unsustainable below this field strength for the low density case.

Close modal

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̂D<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.

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 Varma29 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 B0 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 Goertz11 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>1022 m−3) are well outside the free molecular neutral flow regime assumed by our model (n1019 m−3). Venkataramani and Mattoo33 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äs34 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<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 known20 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̂2n 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 aerocapture7 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.

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.

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

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

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

Braginskii22 defines the ion–electron collision time as

τie=3.44×1011lnΛe3/2Te3/2ne,
(A1)

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

τ̂ie=3.44×1011lnΛe3/2T3/2unrc.
(A2)
B. Plasma diffusion

The Bohm diffusion coefficient DB is given by Lieberman and Lichtenberg23 as

DB=Te16eB,
(A3)

where e is the elementary charge. For electron diffusion due to collisions with heavy particles, the diffusion coefficient Dc is given by Goldston and Rutherford24 as

Dc=me(νei+νen)e2B2(Te+Ti),
(A4)

where νei and νen are the electron–ion and electron–neutral collision frequencies, respectively. The electron–ion collision frequency is given by Goldston and Rutherford,24 Eq. (11.22), as a function of ni and Te,

νei=2e4lnΛ12π3/2ε02meniTe3/2,
(A5)

where ε0 is the vacuum permittivity. The electron–neutral collision frequency is modeled by Itikawa37 as

νen=83(2πTeme)1/2σenn2n,
(A6)

where σen is the average electron–neutral momentum transfer cross section. Energy-dependent values of σen are sourced from Harstad.38 The normalized forms of these diffusion coefficients appearing in Eqs. (26), (29), and (30) are

D̂B=T16eB0rcu,
(A7)
D̂c=meTne2B02rcuIψα(ν¯eiN̂i+ν¯enN̂2n),
(A8)

where ν¯ej=νej/nj is defined to remove the spatial dependence of Dc. Note that D̂c depends on model outputs N̂i and N̂2n 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 σv, they obtain

ṅicxninnσcxVcx,
(A9)

where nn 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

Vcx,sn=4πvT,i2+u2,
(A10)
Vcx,2n=4π(vT,i2+vT,2n2),
(A11)

for thermal velocities vT,j. The CX cross section, σcx(Vcx) is sourced from Phelps39 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 Shumlak17 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

Q̂j,k=v̂T,kσcxnrc4πv̂T,j+649πv̂T,k+v̂j,k2,
(A12)
v̂T,j=2Tj/mju,
(A13)
v̂j,k=ujuku,
(A14)

where σcx is that for jk charge exchange.

2. Ionization

The ionization reaction rate is given analytically by Table 3.3 of Lieberman and Lichtenberg23 as

Riz=2.34×1014Te0.59exp(17.44/Te),
(A15)

for Te in eV. The effective ionization cost, ϵion, is also given in eV units by Lieberman and Lichtenberg in Eq. (3.5.8) as

ϵion=(RizUiz+RexUex+3Me/iRelTe)/Riz,
(A16)

for Te 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

Rel=2.336×1014Te1.609×exp[0.0618ln(Te)20.1171ln(Te)3],
(A17)
Rex=2.48×1014Te0.33exp(12.78/Te),
(A18)

for Te in eV.

3. Neutral diffusion

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

D2n=(niσQm+n2nσgk)1(πT2nm2n)1/2,
(A19)

where σQm and σgk are the cross sections of ion–neutral momentum transfer and neutral–neutral collisions, respectively. The Ar–Ar+ momentum transfer cross section, σQm, is sourced from Phelps.39 We use the highest tabulated value, 1.57×1018 m2, 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, σgk=5×1019 m2, is given by Lieberman and Lichtenberg23 in Table 9.3, reproduced from Smirnov.40 The normalized diffusion coefficient appearing in Eqs. (27) and (31) is

D̂2n=(4π3)1/2(N̂iλ̂Qm+N̂2nλ̂gk)1.
(A20)

Note again that D̂2n is dependent on model outputs N̂i and N̂2n and thus varies with time.

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 Isn, 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̂j and reference density n̂j,r in relation (24) is defined as

Iψα=V*(ψψr)αdV̂,
(B1)

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

Iψ2α=V*(ψψr)2αdV̂.
(B2)

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

Isn=V*(ψ̂ψ̂r)αn̂sndV̂.
(B3)

Computing this integral requires knowledge of n̂sn. We derive n̂sn 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

·(nsnuêz)=Rtotnensn,
(B4)

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

n̂sn=exp[n̂e,rτ̂capẑ(ψψr)αdẑ].
(B5)

The stream capture time scale, τ̂cap, is defined as in Eq. (25) for reaction rate Rtot. Therefore, n̂e,r/τ̂cap represents the characteristic rate at which stream neutrals interact with the plasma. Isn must be computed at each time step of the numerical solution, unlike the other CV integrals. A functional fit for Isn(ρL,n̂e,r/τ̂cap) is provided in Part I16 that facilitates rapid computation with reasonable accuracy.

Volume averaging of diffusion terms results in magnetic field integrals that reflect the spatial scalings associated with DB, Dc, and D2n. The respective integrals are

IB=S*4r̂ψ̂(ψ̂ψ̂r)αdŜ,
(B6)
Ic=S*4r̂ψ̂B̂(ψ̂ψ̂r)2αdŜ,
(B7)
I2n=S*r̂B̂ψ̂dŜ.
(B8)

Note these are surface integrals over S* owing to Gauss' theorem. Once again, IB, Ic, and I2n are all functions of the magnetic field only. It is the different scaling with B in each Dj that yields three unique integrals.

1.
E. T.
Karlson
,
Phys. Fluids
6
,
708
(
1963
).
2.
O.
Buneman
,
T.
Neubert
, and
K.-I.
Nishikawa
,
IEEE Trans. Plasma Sci.
20
,
810
(
1992
).
3.
A.
Hasegawa
,
L.
Chen
, and
M. E.
Mauel
,
Nucl. Fusion
30
,
2405
(
1990
).
4.
R. M.
Winglee
,
J.
Slough
,
T.
Ziemba
, and
A.
Goodson
,
J. Geophys. Res.: Space Phys.
105
,
21067
, (
2000
).
5.
J.
Slough
,
D.
Kirtley
, and
A.
Pancotti
, in
32nd International Electric Propulsion Conference
(Wiesbaden, Germany,
2011
).
6.
C. L.
Kelly
and
J. M.
Little
, in
IEEE Aerospace Conference
(Big Sky, MT,
2019
).
7.
C. L.
Kelly
and
J. M.
Little
, in
36th International Electric Propulsion Conference
(Vienna, Austria,
2019
).
8.
H.
Alfvén
,
On the Origin of the Solar System
edited by
N. F.
Mott
and
E. C.
Bullard
(
Oxford University Press
,
London
,
1954
).
9.
N.
Brenning
,
Space Sci. Rev.
59
,
209
(
1992
).
11.
S.
Machida
and
C. K.
Goertz
,
J. Geophys. Res.
91
,
11965
, (
1986
).
12.
S.
Machida
and
C. K.
Goertz
,
J. Geophys. Res.
93
,
11495
, (
1988
).
13.
M. A.
Raadu
,
Astrophys. Space Sci.
55
,
125
(
1978
).
14.
W. J.
McNeil
,
S. T.
Lai
, and
E.
Murad
,
J. Geophys. Res.
95
,
10345
, (
1990
).
15.
E.
Möbius
,
K.
Papadopoulos
, and
A.
Piel
,
Planet. Space Sci.
35
,
345
(
1987
).
16.
J.
Little
and
C. L.
Kelly
,
Phys. Plasmas
27
(
11
),
113512
(
2020
).
17.
E. T.
Meier
and
U.
Shumlak
,
Phys. Plasmas
19
,
072508
(
2012
).
18.
J. D.
Jackson
,
Classical Electrodynamics
(
John Wiley & Sons, Inc
.,
New York
,
1962
).
19.
K.
Takahashi
and
A.
Ando
,
Phys. Rev. Lett.
118
,
225002
(
2017
).
20.
S. I.
Krasheninnikov
,
P. J.
Catto
, and
R. D.
Hazeltine
,
Phys. Rev. Lett.
82
,
2689
(
1999
).
21.
J.
Kesner
,
D. T.
Garnier
, and
M. E.
Mauel
,
Phys. Plasmas
18
,
050703
(
2011
).
22.
S. I.
Braginskii
, in Reviews of Plasma Physics (Consultants Bureau, New York, 1965), Vol. 1, p. 205.
23.
M. A.
Lieberman
and
A. J.
Lichtenberg
,
Principles of Plasma Discharges and Materials Processing
, 2nd ed. (
John Wiley & Sons, Inc.
,
Hoboken, NJ
,
2005
).
24.
R. J.
Goldston
and
P. H.
Rutherford
,
Introduction to Plasma Physics
, Vol.
1
(
Institute of Physics Publishing
,
Bristol/Philadelphia
,
1995
).
25.
M. S.
Davis
,
M. E.
Mauel
,
D. T.
Garnier
, and
J.
Kesner
,
Plasma Phys. Controlled Fusion
56
,
095021
(
2014
).
26.
M.
Vellante
,
M.
Piersanti
, and
E.
Pietropaolo
,
J. Geophys. Res.: Space Phys.
119
,
2623
, (
2014
).
27.
A. C.
Hindmarsh
and
L. R.
Petzold
, “
LSODA, ordinary differential equation solver for stiff or non-stiff system
,” NEA, 2005, see https://www.osti.gov/etdeweb/biblio/21352532
28.
D. T.
Garnier
,
A. C.
Boxer
,
J. L.
Ellsworth
,
J.
Kesner
, and
M. E.
Mauel
,
Nucl. Fusion
49
,
055023
(
2009
).
29.
J. F.
McKenzie
and
R. K.
Varma
,
J. Plasma Phys.
25
,
491
(
1981
).
30.
B.
Lehnert
,
Phys. Fluids
9
,
774
(
1966
).
31.
L.
Danielsson
,
Phys. Fluids
13
,
2288
(
1970
).
32.
G.
Himmel
,
E.
Möbius
, and
A.
Piel
, “
Zeitschrift fur Naturforschung—Section A
,”
J. Phys. Sci.
31
,
934
(
1976
).
33.
N.
Venkataramani
and
S. K.
Mattoo
,
Phys. Lett. A
79
,
393
(
1980
).
34.
N.
Brenning
and
I.
Axnäs
,
Astrophys. Space Sci.
144
,
15
(
1988
).
35.
S.
Kobayashi
,
B. N.
Rogers
, and
W.
Dorland
,
Phys. Rev. Lett.
103
,
055003
(
2009
).
36.
Z. R.
Putnam
and
R. D.
Braun
,
J. Spacecr. Rockets
51
,
139
(
2014
).
37.
Y.
Itikawa
,
Phys. Fluids
16
,
831
(
1973
).
38.
K. G.
Harstad
, “
Transport equations for gases and plasmas obtained by the 13-moment method—A summary
,” Report No. NASA-CR-97148, NASA Jet Propulsion Lab, Pasadena, CA,
1968
.
39.
A. V.
Phelps
,
J. Phys. Chem. Ref. Data
20
,
557
(
1991
).
40.
B. M.
Smirnov
,
Introduction to Plasma Physics
(
Mir
,
Moscow
,
1977
).