# Rapid Mass Parameter Estimation of
Binary Black Hole Coalescences Using Deep Learning

Alistair McLeod<EMAIL_ADDRESS>
<EMAIL_ADDRESS>Daniel Jacobs  Chayan Chatterjee  Linqing Wen
<EMAIL_ADDRESS>Fiona Panther  Department of Physics, OzGrav-UWA, The
University of Western Australia, 35 Stirling Hwy, Crawley, Western Australia
6009, Australia

###### Abstract

Deep learning can be used to drastically decrease the processing time of
parameter estimation for coalescing binaries of compact objects including
black holes and neutron stars detected in gravitational waves (GWs). As a
first step, we present two neural network models trained to rapidly estimate
the posterior distributions of the chirp mass and mass ratio of a detected
binary black hole system from the GW strain data of LIGO Hanford and
Livingston Observatories. Using these parameters the component masses can be
predicted, which has implications for the prediction of the likelihood that a
merger contains a neutron star. The results are compared to the ‘gold
standard’ of parameter estimation of gravitational waves used by the LIGO-
Virgo Collaboration (LVC), LALInference. Our models predict posterior
distributions consistent with that from LALInference while using orders of
magnitude less processing time once the models are trained. The median
predictions are within the 90% credible intervals of LALInference for all
predicted parameters when tested on real binary black hole events detected
during the LVC’s first and second observing runs. We argue that deep learning
has strong potential for low-latency high-accuracy parameter estimation
suitable for real-time GW search pipelines.

††preprint: APS/123-QED

##  I Introduction

Gravitational waves are ripples in spacetime, first predicted by Einstein [1],
and first successfully detected by LIGO [2] in 2015. The majority of events
detected by the LIGO-Virgo Collaboration (LVC) [3, 4] are binary black hole
(BBH) mergers [5, 6, 7], but binary neutron star (BNS) mergers [8] and neutron
star - black hole mergers (NSBH) [9] have also been detected. Most significant
gravitational wave events were detected first in online real-time searches
which use either matched filtering or burst search methods for signal
filtering [10, 11, 12, 13, 14]. Detecting these events in real-time is of
importance as it enables rapid follow-up observations of the event by
electromagnetic telescopes. To this end, rapid parameter estimation is needed
to assess the likelihood that the merger contains a neutron star that is
considered more likely to be associated with electromagnetic radiation.

These compact binary systems are described by a 15 dimensional parameter
space, consisting of intrinsic attributes of the two merging compact objects
(component mass, spin etc), and extrinsic parameters such as the orientation
of the merging system to a detector. Once a gravitational wave is detected in
detector data by a real-time search pipeline, the pipeline performs parameter
estimation on the detected event [10]. Parameter estimation is currently
performed as a Bayesian analysis with Markov Chain Monte Carlo and nested
sampling of the parameter posteriors [15, 16, 17, 18]. While this is extremely
accurate, it is also extremely computationally intensive, and can take hours
or days for a full parameter estimation of a binary black hole system to be
calculated [19]. This is insufficient for BNS mergers, as prompt
electromagnetic radiation could be emitted within seconds of the merger [20].
Rapid parameter estimation of gravitational wave signals is needed to help
capture these prompt electromagnetic emissions.

Several machine learning systems have been developed for BBH parameter
estimation [21, 22, 23]. Neural networks of varying architectures have been
shown to be able to accurately approximate Bayesian posteriors for source
parameters as a Gaussian mixture [24], and even perform full source parameter
inference [25, 26]. Shen  _et al._ [27] showcases a Bayesian neural network
that could produce Bayesian posteriors by repeated sampling with uncertainty
of network weights. This method produces accurate posteriors significantly
faster than traditional methods, but due to the repeated sampling method,
slower than simpler neural network architectures could. Gabbard  _et al._ [28]
use conditional variational autoencoders to accurately predict BBH source
parameters, however they have not tested their network on signals embedded in
real noise yet. Dax  _et al._ [26] applied normalising flows to accurately
infer the posteriors of all 15 source parameters, and can produce posteriors
that almost perfectly match those produced by LALInference. However, they use
frequency domain gravitational wave signals. Frequency domain signals need to
accumulate data, which could introduce additional latency. Additionally, in
order to achieve this level of accuracy, their network needs to be trained for
several weeks, and when predicting must sample the posterior several thousand
times which takes several seconds. This may not be fast enough for online
parameter inference of GW mergers which will likely be the next target of ML-
based inference, to enable rapid EM followups.

In this work, we present a simple system of neural networks which can rapidly
predict the chirp mass and mass ratio of a BBH signal in detector noise. BBH
signals were used instead of BNS signals, due to the greater length of BNS
signals leading to impractically large input layers. These networks are
designed for speed and simplicity, meaning they can be rapidly re-trained and
are are well suited for integration into online search pipelines [10]. Our
method produces a mean and standard deviation prediction for each parameter
(thus assuming a Gaussian distribution), which can be used to avoid the
computational cost of sampling to produce a posterior, with no loss in
prediction power over Bayesian prediction methods. By estimating both the
chirp mass and mass ratio, it is possible to estimate the two component
masses, $m_{1}$ and $m_{2}$. We apply this system to the 10 BBH events from
the first two gravitational wave observing runs, and report on the network’s
performance.

##  II Network architecture

The two neural network prediction systems presented here for chirp mass and
mass ratio respectively use the same general architecture, with slightly
different parameters and training as outlined in section III. Both prediction
systems consist of a denoising autoencoder, whose output is passed to a dense
multi-branched convolutional neural network which predicts the relevant
parameter.

The first component of both systems is a denoising autoencoder [29], which
cleans the GW time series data from LIGO Hanford and Livingston detectors to
remove as much noise as possible. Autoencoders are neural networks which aim
to reconstruct their own inputs, and are designed to learn a dimension-
efficient representation of the input data. An autoencoder is very effective
at denoising signals, recovering the ‘clean’ structured waveform while most
contaminating noise is removed through the compression of the data. Through
trial and error, we found that at high signal-to-noise ratios (SNRs), SNR
$>30$ , the dense estimator network can predict the parameters reasonably well
on its own. However, as the SNR drops into the level normally seen in real
signals, the noise becomes overwhelmingly high and the network struggles to
identify the actual waveform among the noise by itself. The denoising
autoencoder assists here by handling the noise removal, outputting a cleaner
signal which is fed to the estimator network, greatly improving the
estimator’s performance. Additional information on the architecture, training
and performance of the autoencoder can be found at Chatterjee  _et al._ [29].

We then apply a second multi-branched dense neural network to predict the
parameter of interest from the output of the autoencoder (Fig. 1). The network
used is a convolutional neural network (CNN), created using the Keras library
[30] with Tensorflow [31] backend. The network receives denoised GW signals as
input, and produces a mean and standard deviation prediction for each signal.
First, the detector signals are passed to two separate input branches, which
consist of multiple 1D convolution layers and residual blocks. Residual blocks
are a type of neural network structure in which the output of a layer ‘skips’
several layers, and is then added to a later output. The advantage of this
type of layer structure is that the intervening layers learn the residual of
the output and input, which has been shown to improve deep network performance
[32]. These branches are then concatenated into a single branch, which
consists of several dense layers. We found that the chirp mass estimator
performed best with 256 neurons in the dense layer, while the mass ratio
estimator performed best with 64 neurons per layer. Dropout is included in the
dense layers to reduce overfitting and to improve the network’s performance on
unseen signals [33]. After the dense layers, an IndependentNormal Tensorflow
layer is used to produce a mean and standard deviation prediction as the
network’s output. The networks were trained using the ADAM optimiser [34], and
the negative log likelihood loss function was used. Training typically took 1
hour to complete when done using an NVIDIA P100 GPU.

Figure 1: Flow chart representing our CNN design. The input data are the
Hanford and Livingston strains that have been denoised by the autoencoder.
Additional information on autoencoder architecture can be found in Chatterjee
_et al._ [29]. Repeated sections of the network have been replaced with a
‘x3’, for simplicity.

###  II.1 Training Set

Both parameter estimator networks were trained on 500,000 BBH gravitational
wave events. These samples were generated using the SEOBNRv4 waveform
approximant [35], and were injected into Gaussian noise coloured by the
detectors’ power spectral density (PSD) from the second observing run, O2.
These signals were then passed through the autoencoder, with the output of the
autoencoder being the training set for the parameter estimator networks. These
signals are each a quarter of a second in length and sampled at 2048 Hz, which
ensures the waveform covers all detectable BBH signals in our mass range. Each
network’s training set was uniformly distributed in the relevant parameter
(i.e. chirp mass and mass ratio) to avoid the networks becoming biased by the
training set. While this is not a realistic astrophysical distribution for
real gravitational wave events (high mass ratios are strongly favoured over
low mass ratios in detected BBH events [6], for example), it ensures the
networks are generalised enough to make accurate predictions on future events.
Both networks were trained with parameters from Tab. 1. Note that forcing a
uniform distribution in either the chirp mass or mass ratio leads to a non-
uniform distribution in both $m_{1}$ and $m_{2}$. The coalescence time was not
varied in the samples, as it was deemed not relevant for mass estimation.

Parameter | Range
---|---
Mc (M⊙) | [6, 70]
q | [0.1, 1]
m1,m2 (M⊙) | [7, 100]
distance (Mpc) | [40, 3000]
right ascension | [0, 2$\pi$]
declination | [0, $\pi$]
SNR | [10, 30]
spin | [-1, 0.99]
inclination | [0, $\pi$]
Table 1:  Parameters used to create the training sets. All non-mass parameters
were distributed uniformly, and the chirp mass ($M_{c}$) and mass ratio (q)
were distributed uniformly for their respective networks. The target SNR was
achieved by calculating a scaling factor from the optimal matched filtering
SNR, then scaling the signal amplitude by this factor [36].

The chirp mass and mass ratio of a gravitational wave event are defined in
Equations 1,2. These parameters determine the frequency evolution of the
compact binary coalescence (CBC), and are therefore easier to measure from the
gravitational wave signal than the individual masses, $m_{1}$ and $m_{2}$
(where $m_{1}>m_{2}$) [37]. The component masses can still be recovered by
solving these chirp mass and mass ratio equations for $m_{1}$ and $m_{2}$.

$M_{c}=\frac{(m_{1}m_{2})^{\frac{3}{5}}}{(m_{1}+m_{2})^{\frac{1}{5}}}$ (1)

$q=\frac{m_{2}}{m_{1}}$ (2)

##  III Prediction Results

The parameter estimator systems were tested on a set of 25,000 simulated
signals with SNRs between 10 and 30, which were not part of the training set.
The signals were injected into Gaussian noise coloured by the O2 PSD. The
network’s ability to produce fast and accurate predictions and useful
probability distributions were tested. The networks were able to produce a
prediction for each sample in $\sim$ 1.4 milliseconds when tested with a
NVIDIA P100 GPU and an Intel Gold 6140 CPU, a prediction speed which is faster
than that shown by Shen  _et al._ [27], Chua and Vallisneri [24] and Green
_et al._ [38]. The networks were shown to produce useful confidence intervals
using a probability-probability plot (Fig. 2). This plot was constructed by
taking the smallest confidence interval produced by the networks that
contained the true value for each sample, and then binning the data and
producing a cumulative distribution. This shows that the networks are able to
produce useful confidence intervals as every confidence level contained a
number of true parameters consistent with the expected count, within 5$\%$.

Figure 2: p-p plot for the 25,000 predictions made on the Gaussian test set
for both parameters. Each curve shows the cumulative distribution of the
confidence levels containing the true parameter.

The chirp mass estimator had a lower average error than the mass ratio
estimator (Fig. 3), which was expected as the mass ratio of a signal is not as
strongly constrained as its chirp mass [37]. The chirp mass estimator showed
no significant bias in overestimating/underestimating, with a 51/49 split. The
mass ratio estimator showed a slight bias toward underestimating, however,
with a 40/60 split. This can be seen in the error histograms for both
networks, as the chirp mass predictions appear symmetrically distributed,
while the mass ratio predictions are skewed towards underestimating.

Figure 3: Error histograms for the two networks. Note that the average error
on the mass ratio network is higher, and skewed towards underestimating the
mass ratio.

As expected, the accuracy of both networks improved with increasing SNR (Fig.
4). This shows that the networks perform better on higher SNR signals,
presumably due to the autoencoder producing more accurate waveforms. Since the
SNR of a signal is dependent on the chirp mass (due to the chirp mass
affecting the signal strength), this improvement may partially be due to the
chirp mass of the signal as seen in Fig 5.

Figure 4: Binned median parameter estimator prediction error versus network
SNR on test set.

###  III.1 Estimator network results

The chirp mass system achieved a mean error of $10.53\%$ on the Gaussian test
set, with the mean error decreasing with increasing chirp mass (see Fig. 5).
The greater prediction error seen at low chirp masses is partially due to the
autoencoder, which produces signal representations with poorer overlaps at
lower chirp masses [29]. Even at low SNRs the network was able to produce
predictions with relatively low average error of 16% compared to the mean of
10.53%.

The network was then tested on the real data from the 10 O1 and O2 BBH events,
that had been passed through the denoising autoencoder. The chirp mass
estimator produced predictions with a $6.20\%$ average error between the
predicted median and LAL simulations. Every 90% confidence interval produced
by the network contained the true value produced by LAL simulations.

The second estimator network was trained to predict the mass ratio q. The
architecture of the network remained largely the same, but with one-quarter
the amount of neurons in the dense layers than were used in the chirp mass
network.

The mass ratio system achieved a mean error of $24.26\%$ on the Gaussian test
set (Fig. 6), There was a strong negative correlation with SNR, with
prediction accuracy dropping off rapidly below an SNR of 20. We also find that
with increasing SNR the network error decreases, as with the chirp mass
network (4).

The network was then tested on the real data from the 10 O1 and O2 BBH events,
that had been passed through the denoising autoencoder. The mass ratio
estimator produced predictions with a $8.50\%$ average error between the
predicted median and LAL simulations. This error is significantly lower than
on the Gaussian test set, due to the 10 real events having predicted mass
ratios of $0.56<q<0.86$, which corresponds to the lowest error region of the
network’s predictions (Fig. 6).

Figure 5: Binned median chirp mass error for the Gaussian test set. Figure 6:
Binned median mass ratio error for the Gaussian test set.

###  III.2 Component Masses

By using the mass ratio and chirp mass predicting networks together it is
possible to derive the component masses of a CBC event. This is done by making
a prediction for both the mass ratio and chirp mass as described above, and
then solving equations 1 and 2 for the component masses. Repeated sampling is
therefore necessary to create a posterior for the component masses, should it
be required. The predicted posteriors for the component masses were then
compared to the LALInference predictions, with the percent error in the
primary mass being $6.73\%$ and the percent error in the secondary mass being
$7.87\%$. The median predictions for both masses fell within the LALInference
credible intervals for every event, and conversely the predicted 90%
confidence intervals captured the LALInference median predictions for every
event. This shows that the networks have learned an accurate representation of
the signals that can compete with Bayesian inference.

A full analysis of the two component mass results are shown in tables 2 and 3.
It can be seen in Fig. 7 that the predicted credible intervals and median
predictions closely match that of the LALInference simulations, with the other
event prediction figures available in the appendix.

Figure 7: Plot of the parameter estimation for GW150914. The solid red vertical lines are the median estimates produced by LAL simulation, and the gray areas are the 90% credible intervals from the LAL simulations. The vertical dashed lines are the 90% credible intervals for our predictions. Every plot has a histogram of the network’s predictions for that parameter with the component mass predictions being calculated using the chirp mass and mass ratio predictions. These histograms were produced by sampling the chirp mass and mass ratio distributions 10,000 times, as this is necessary to produce the component mass distributions Note that the mass ratio posteriors are truncated at 1 (due to greater mass ratios being unphysical) and 0.1 (due to very small mass ratios leading to extreme component mass predictions). Additional plots of the other O1 and O2 events can be found in the appendix. Table 2:  The results of predicting the component mass $m_{1}$ on the first 10 BBH events detected by LIGO. The predictions and credible intervals are in units of solar masses ($M_{\odot}$) | LAL $m_{1}$ median | LAL 90% |  $m_{1}$ deep | Deep learning |  | Deep learning | Deep learning
---|---|---|---|---|---|---|---
Signal | prediction | CI | learning median | 90% CI | Error % | median within | cred. interval
|  |  |  |  |  | LAL cred. bounds? | captures LAL prediction?
GW150914 | 35.6 | [32.6, 40.4] | 36.91 | [30.11, 47.77] | 3.68 | yes | yes
GW151012 | 23.3 | [17.8, 37.3] | 21.32 | [16.76, 31.77] | 8.51 | yes | yes
GW151226 | 13.7 | [11.5, 22.5] | 15.95 | [11.46, 28.13] | 16.44 | yes | yes
GW170104 | 31.0 | [25.4, 38.2] | 27.44 | [22.29, 36.43] | 11.48 | yes | yes
GW170608 | 10.9 | [9.2, 16.2] | 11.34 | [7.68, 18.04] | 4.03 | yes | yes
GW170729 | 50.6 | [40.4, 67.2] | 50.37 | [38.16, 70.56] | 0.45 | yes | yes
GW170809 | 35.2 | [29.2, 43.5] | 32.36 | [26.77, 41.92] | 8.08 | yes | yes
GW170814 | 30.7 | [27.7, 36.4] | 31.75 | [26.46, 40.59] | 3.43 | yes | yes
GW170818 | 35.5 | [30.8, 43.0] | 33.23 | [26.89, 44.75] | 6.38 | yes | yes
GW170823 | 39.6 | [33.0, 49.6] | 37.71 | [30.31, 52.24] | 4.77 | yes | yes
|  |  |  |  | 6.73 | 10/10 | 10/10
Table 3:  The results of predicting the component mass $m_{2}$ on the first 10 BBH events detected by LIGO. The predictions and credible intervals are in units of solar masses ($M_{\odot}$) | LAL $m_{2}$ median | LAL 90% |  $m_{2}$ deep | Deep learning |  | Deep learning | Deep learning
---|---|---|---|---|---|---|---
Signal | prediction | CI | learning median | 90% CI | Error % | median within | cred. interval
|  |  |  |  |  | LAL cred. bounds? | captures LAL prediction?
GW150914 | 30.6 | [26.2, 33.6] | 26.68 | [21.06, 32.29] | 12.82 | yes | yes
GW151012 | 13.6 | [8.8, 17.7] | 13.57 | [9.47, 17.04] | 0.19 | yes | yes
GW151226 | 7.7 | [5.1, 9.9] | 8.55 | [5.32, 11.52] | 11.04 | yes | yes
GW170104 | 20.1 | [15.6, 25.0] | 19.15 | [14.73, 23.30] | 4.70 | yes | yes
GW170608 | 7.6 | [5.5, 8.9] | 6.69 | [4.18, 9.33] | 12.01 | yes | yes
GW170729 | 34.3 | [24.2, 43.4] | 34.54 | [24.67, 44.57] | 0.70 | yes | yes
GW170809 | 23.8 | [18.7, 29.0] | 23.22 | [18.19, 27.84] | 2.46 | yes | yes
GW170814 | 25.3 | [21.2, 28.2] | 23.45 | [18.66, 27.97] | 7.32 | yes | yes
GW170818 | 26.8 | [21.6, 31.1] | 22.93 | [17.50, 28.07] | 14.43 | yes | yes
GW170823 | 29.4 | [22.3, 35.7] | 25.58 | [19.15, 31.37] | 12.99 | yes | yes
|  |  |  |  | 7.87 | 10/10 | 10/10

##  IV Discussion

While our networks do not reach the accuracy achieved by other neural network
BBH parameter estimators [26], their utility lies in their speed, in both
training and prediction. Our estimator networks are able to produce Gaussian
posterior predictions for the chirp mass and mass ratio of a CBC within 1.4
milliseconds, and can be trained in an hour on 500,000 GW samples. The
component mass errors for these predictions relative to LALInference are
within 10% when tested on real events, and these features combined make our
networks ideal to implement in an online search pipeline such as SPIIR [10].
This would be especially useful for early warning [39], where the rapid
classification of GW signals is necessary. Early warning in gravitational wave
search pipelines would be of interest for transient astronomers searching for
electromagnetic counterparts [40].

The simulated signals used in this work were coloured with the PSD of the
detectors during the O2 science run. To improve our results we can train the
network on signals coloured by O3 detector noise, and to also test on the BBH
signals observed during O3. The three input variant including input from the
Virgo detector should also be tested on the dataset coloured with the O3 PSD.
It is possible that the inclusion of Virgo may improve the model’s performance
and bring it even closer to LALInference’s level of accuracy.

In this work the neural networks were only trained and tested on BBH systems,
however, the LSC has also detected CBCs of systems containing neutron stars
(e.g. GW170817). These signals typically exist in the LSC detector sensitivity
band for as long as 500 seconds, and thus are much more difficult for a neural
network to analyse. However, the neural networks may be able to predict
sufficiently small chirp masses on systems containing neutron stars to infer
the existence of a neutron star in a signal with minimal modification.

The denoising autoencoder model produces a single point estimate of a pure GW
waveform from noisy strain data. However, there is no unique partition of the
data into signal and noise. The CNN is therefore unable to learn the inherent
stochasticity in the data distribution, which may explain why our predicted
mass posteriors do not exactly match those of LALInference, which models the
data distribution as stationary and Gaussian. To improve this, the models
could be joined together as a single model, with an additional ‘uncertainty’
parameter being passed to the parameter estimator along with the extracted
waveform. This may lead to posteriors which more closely match the
LALInference posteriors, as the parameter estimation network would be better
equipped to learn the appropriate uncertainty in the relevant parameter.

##  V Conclusion

This work presented a deep learning system composed of a denoising autoencoder
and an estimator that predicts the posterior distribution of the chirp mass of
a binary black hole system given the gravitational wave signal emitted during
coalescence. The system was trained and tested on binary black hole
coalescence gravitational wave signals immersed in Gaussian noise coloured by
the PSD of the LVC interferometers during the O2 science run. The system has
an average error of $10.53\%$ for chirp mass and $24.26\%$ for the mass ratio
across the test set of signals of SNRs between 10 and 30, but is several
orders of magnitudes faster than LALInference, taking $\sim$ 1 millisecond to
make a prediction on each event. This work shows the promise of deep learning
applied to gravitational wave parameter estimation by approaching the accuracy
of LALInference despite being much lower latency. This shows that deep
learning can take on the role of a rapid first estimation for parameters while
the costly LALInference estimation can be performed offline for a more
accurate estimate.

This work could be extended to estimating the chirp mass and mass ratio of
neutron star binary systems with minor changes to the input dimensions. The
system design presented here should also be trained and tested using real
detector noise from the O3 science run. If the system also performs well on
real detector noise or can be extended to do so, it will be ready to be
implemented in the LVC parameter estimation pipeline as a rapid initial
parameter estimator.

##  VI Acknowledgements

AM and DJ utilized the OzSTAR national facility at Swinburne University of
Technology. The OzSTAR program receives funding in part from the Astronomy
National Collaborative Research Infrastructure Strategy (NCRIS) allocation
provided by the Australian Government. AM and DJ acknowledge support from
Australian Research Council Centre of Excellence for Gravitational Wave
Discovery (OzGrav).

We thank Leïla Haegel and Stephen Green for the helpful comments on the paper.

This research has made use of data, software and/or web tools obtained from
the Gravitational Wave Open Science Center (https://www.gw-openscience.org/),
a service of LIGO Laboratory, the LIGO Scientific Collaboration and the Virgo
Collaboration. LIGO Laboratory and Advanced LIGO are funded by the United
States National Science Foundation (NSF) as well as the Science and Technology
Facilities Council (STFC) of the United Kingdom, the Max-Planck-Society (MPS),
and the State of Niedersachsen/Germany for support of the construction of
Advanced LIGO and construction and operation of the GEO600 detector. Virgo is
funded, through the European Gravitational Observatory (EGO), by the French
Centre National de Recherche Scientifique (CNRS), the Italian Istituto
Nazionale di Fisica Nucleare (INFN) and the Dutch Nikhef, with contributions
by institutions from Belgium, Germany, Greece, Hungary, Ireland, Japan,
Monaco, Poland, Portugal, Spain. This material is based upon work supported by
NSF’s LIGO Laboratory which is a major facility fully funded by the National
Science Foundation.

CC, LW and FP acknowledge funding and support from Australian Research Council
Centre of Excellence for Gravitational Wave Discovery (OzGrav) under grant
CE170100004.

##  VII Appendices

Figure 8: Plot of the parameter estimation for GW151012. Figure 9: Plot of the
parameter estimation for GW151226. Figure 10: Plot of the parameter estimation
for GW170104. Figure 11: Plot of the parameter estimation for GW170608. Figure
12: Plot of the parameter estimation for GW170729. Figure 13: Plot of the
parameter estimation for GW170809. Figure 14: Plot of the parameter estimation
for GW170814. Figure 15: Plot of the parameter estimation for GW170818. Figure
16: Plot of the parameter estimation for GW170823.

## References

  * Einstein [1916] A. Einstein, Die Grundlage der allgemeinen Relativitätstheorie, Annalen der Physik 354, 769 (1916).
  * Abbott  _et al._ [2016] B. Abbott _et al._ (LIGO Scientific, Virgo), Observation of Gravitational Waves from a Binary Black Hole Merger, Phys. Rev. Lett. 116, 061102 (2016), arXiv:1602.03837 [gr-qc] .
  * LIGO Scientific Collaboration  _et al._ [2015] LIGO Scientific Collaboration, J. Aasi, _et al._ , Advanced LIGO, Classical and Quantum Gravity 32, 074001 (2015), arXiv:1411.4547 [gr-qc] .
  * Acernese  _et al._ [2015] F. Acernese _et al._ , Advanced Virgo: a second-generation interferometric gravitational wave detector, Classical and Quantum Gravity 32, 024001 (2015), arXiv:1408.3978 [gr-qc] .
  * Abbott  _et al._ [2019] B. Abbott _et al._ (LIGO Scientific, Virgo), GWTC-1: A Gravitational-Wave Transient Catalog of Compact Binary Mergers Observed by LIGO and Virgo during the First and Second Observing Runs, Phys. Rev. X 9, 031040 (2019), arXiv:1811.12907 [astro-ph.HE] .
  * Abbott  _et al._ [2021a] R. Abbott _et al._ (LIGO Scientific, Virgo), GWTC-2: Compact Binary Coalescences Observed by LIGO and Virgo During the First Half of the Third Observing Run, Phys. Rev. X 11, 021053 (2021a), arXiv:2010.14527 [gr-qc] .
  * The LIGO Scientific Collaboration  _et al._ [2021] The LIGO Scientific Collaboration, the Virgo Collaboration, the KAGRA Collaboration, R. Abbott, _et al._ , GWTC-3: Compact Binary Coalescences Observed by LIGO and Virgo During the Second Part of the Third Observing Run, arXiv e-prints , arXiv:2111.03606 (2021), arXiv:2111.03606 [gr-qc] .
  * Abbott  _et al._ [2017] B. Abbott _et al._ (LIGO Scientific, Virgo), GW170817: Observation of Gravitational Waves from a Binary Neutron Star Inspiral, Phys. Rev. Lett. 119, 161101 (2017), arXiv:1710.05832 [gr-qc] .
  * Abbott  _et al._ [2021b] R. Abbott _et al._ (LIGO Scientific, VIRGO, KAGRA), Observation of gravitational waves from two neutron star-black hole coalescences, Astrophys. J. Lett. 915, L5 (2021b), arXiv:2106.15163 [astro-ph.HE] .
  * Chu  _et al._ [2020] Q. Chu _et al._ , The SPIIR online coherent pipeline to search for gravitational waves from compact binary coalescences, arXiv e-prints , arXiv:2011.06787 (2020), arXiv:2011.06787 [gr-qc] .
  * Messick  _et al._ [2017] C. Messick _et al._ , Analysis framework for the prompt discovery of compact binary mergers in gravitational-wave data, Phys. Rev. D 95, 042001 (2017), arXiv:1604.04324 [astro-ph.IM] .
  * Nitz  _et al._ [2017] A. H. Nitz, T. Dent, T. Dal Canton, S. Fairhurst, and D. A. Brown, Detecting Binary Compact-object Mergers with Gravitational Waves: Understanding and Improving the Sensitivity of the PyCBC Search, Astrophys. J. 849, 118 (2017), arXiv:1705.01513 [gr-qc] .
  * Adams  _et al._ [2016] T. Adams, D. Buskulic, V. Germain, G. M. Guidi, F. Marion, M. Montani, B. Mours, F. Piergiovanni, and G. Wang, Low-latency analysis pipeline for compact binary coalescences in the advanced gravitational wave detector era, Classical and Quantum Gravity 33, 175012 (2016), arXiv:1512.02864 [gr-qc] .
  * Klimenko  _et al._ [2016] S. Klimenko _et al._ , Method for detection and reconstruction of gravitational wave transients with networks of advanced detectors, Phys. Rev. D 93, 042004 (2016), arXiv:1511.05999 [gr-qc] .
  * Thrane and Talbot [2019] E. Thrane and C. Talbot, An introduction to Bayesian inference in gravitational-wave astronomy: Parameter estimation, model selection, and hierarchical models, Publications of the Astron. Soc. of Australia 36, e010 (2019), arXiv:1809.02293 [astro-ph.IM] .
  * Veitch  _et al._ [2015] J. Veitch _et al._ , Parameter estimation for compact binaries with ground-based gravitational-wave observations using the lalinference software library, Phys. Rev. D 91, 042003 (2015).
  * Ashton  _et al._ [2019] G. Ashton _et al._ , BILBY: A User-friendly Bayesian Inference Library for Gravitational-wave Astronomy, Astrophys. J., Supp. 241, 27 (2019), arXiv:1811.02042 [astro-ph.IM] .
  * Lange  _et al._ [2018] J. Lange, R. O’Shaughnessy, and M. Rizzo, Rapid and accurate parameter inference for coalescing, precessing compact binaries, arXiv e-prints , arXiv:1805.10457 (2018), arXiv:1805.10457 [gr-qc] .
  * Berry  _et al._ [2015] C. P. L. Berry _et al._ , Parameter estimation for binary neutron-star coalescences with realistic noise during the advanced ligo era, Astrophys. J. 804, 114 (2015).
  * Goldstein  _et al._ [2017] A. Goldstein _et al._ , An Ordinary Short Gamma-Ray Burst with Extraordinary Implications: Fermi-GBM Detection of GRB 170817A, Astrophys. J., Lett. 848, L14 (2017), arXiv:1710.05446 [astro-ph.HE] .
  * Cuoco  _et al._ [2020] E. Cuoco _et al._ , Enhancing Gravitational-Wave Science with Machine Learning, arXiv e-prints , arXiv:2005.03745 (2020), arXiv:2005.03745 [astro-ph.HE] .
  * Chatterjee  _et al._ [2019] C. Chatterjee, L. Wen, K. Vinsen, M. Kovalam, and A. Datta, Using deep learning to localize gravitational wave sources, Phys. Rev. D 100, 103025 (2019), arXiv:1909.06367 [astro-ph.IM] .
  * Haegel and Husa [2020] L. Haegel and S. Husa, Predicting the properties of black-hole merger remnants with deep neural networks, Classical and Quantum Gravity 37, 135005 (2020), arXiv:1911.01496 [gr-qc] .
  * Chua and Vallisneri [2020] A. J. Chua and M. Vallisneri, Learning bayesian posteriors with neural networks for gravitational-wave inference, Phys. rev. lett. 124, 041102 (2020).
  * Green and Gair [2020] S. R. Green and J. Gair, Complete parameter inference for GW150914 using deep learning, arXiv e-prints , arXiv:2008.03312 (2020), arXiv:2008.03312 [astro-ph.IM] .
  * Dax  _et al._ [2021] M. Dax, S. R. Green, J. Gair, J. H. Macke, A. Buonanno, and B. Schölkopf, Real-time gravitational-wave science with neural posterior estimation, arXiv e-prints , arXiv:2106.12594 (2021), arXiv:2106.12594 [gr-qc] .
  * Shen  _et al._ [2019] H. Shen, E. A. Huerta, E. O’Shea, P. Kumar, and Z. Zhao, Statistically-informed deep learning for gravitational wave parameter estimation, arXiv e-prints , arXiv:1903.01998 (2019), arXiv:1903.01998 [gr-qc] .
  * Gabbard  _et al._ [2019] H. Gabbard, C. Messenger, I. S. Heng, F. Tonolini, and R. Murray-Smith, Bayesian parameter estimation using conditional variational autoencoders for gravitational-wave astronomy, arXiv e-prints , arXiv:1909.06296 (2019), arXiv:1909.06296 [astro-ph.IM] .
  * Chatterjee  _et al._ [2021] C. Chatterjee, L. Wen, F. Diakogiannis, and K. Vinsen, Extraction of binary black hole gravitational wave signals from detector data using deep learning, Phys. Rev. D 104, 064046 (2021), arXiv:2105.03073 [gr-qc] .
  * Chollet [2015] F. Chollet, keras, https://github.com/fchollet/keras (2015).
  * Abadi  _et al._ [2015] M. Abadi _et al._ , TensorFlow: Large-scale machine learning on heterogeneous systems (2015), software available from tensorflow.org.
  * He  _et al._ [2015] K. He, X. Zhang, S. Ren, and J. Sun, Deep residual learning for image recognition (2015), arXiv:1512.03385 [cs.CV] .
  * Hinton  _et al._ [2012] G. E. Hinton, N. Srivastava, A. Krizhevsky, I. Sutskever, and R. R. Salakhutdinov, Improving neural networks by preventing co-adaptation of feature detectors (2012), arXiv:1207.0580 [cs.NE] .
  * Kingma and Ba [2017] D. P. Kingma and J. Ba, Adam: A method for stochastic optimization (2017), arXiv:1412.6980 [cs.LG] .
  * Bohé  _et al._ [2017] A. Bohé _et al._ , Improved effective-one-body model of spinning, nonprecessing binary black holes for the era of gravitational-wave astrophysics with advanced detectors, Phys. Rev. D 95, 044028 (2017), arXiv:1611.03703 [gr-qc] .
  * Gebhard  _et al._ [2019] T. D. Gebhard, N. Kilbertus, I. Harry, and B. Schölkopf, Convolutional neural networks: A magic bullet for gravitational-wave detection?, Phys. Rev. D 100, 063015 (2019), arXiv:1904.08693 [astro-ph.IM] .
  * Cutler and Flanagan [1994] C. Cutler and É. E. Flanagan, Gravitational waves from merging compact binaries: How accurately can one extract the binary’s parameters from the inspiral waveform\?, Phys. Rev. D 49, 2658 (1994), arXiv:gr-qc/9402014 [gr-qc] .
  * Green  _et al._ [2020] S. R. Green, C. Simpson, and J. Gair, Gravitational-wave parameter estimation with autoregressive neural network flows, Phys. Rev. D 102, 104057 (2020), arXiv:2002.07656 [astro-ph.IM] .
  * Magee  _et al._ [2021] R. Magee _et al._ , First Demonstration of Early Warning Gravitational-wave Alerts, Astrophys. J., Lett. 910, L21 (2021), arXiv:2102.04555 [astro-ph.HE] .
  * James  _et al._ [2019] C. W. James, G. E. Anderson, L. Wen, J. Bosveld, Q. Chu, M. Kovalam, T. J. Slaven-Blair, and A. Williams, Using negative-latency gravitational wave alerts to detect prompt radio bursts from binary neutron star mergers with the Murchison Widefield Array, Monthly Notices of the RAS 489, L75 (2019), arXiv:1908.08688 [astro-ph.HE] .


# Partial-Ionization Deconfinement Effect in Magnetized Plasma

M. E. Mlodik<EMAIL_ADDRESS>E. J. Kolmes  I. E. Ochs  T. Rubin  N. J.
Fisch  Department of Astrophysical Sciences, Princeton University, Princeton,
New Jersey, USA, 08543 and
Princeton Plasma Physics Laboratory, Princeton, New Jersey, USA, 08540

###### Abstract

In partially ionized plasma, where ions can be in different ionization states,
each charge state can be described as a different fluid for the purpose of
multi-ion collisional transport. In the case of two charge states, transport
pushes plasma toward equilibrium which is found to be a combination of local
charge state equilibrium and generalized pinch relations between ion fluids
representing different charge states. Combined, these conditions lead to a
dramatic deconfinement of ions. This deconfinement happens on the timescale
similar but not identical to the multi-ion cross-field transport timescale, as
opposed to electron-ion transport timescale in fully ionized plasma.
Deconfinement occurs because local charge state equilibration enforces the
disparity in diamagnetic drift velocities of ion fluid components, which in
turn leads to the cross-field transport due to ion-ion friction.

##  I Introduction.

A new deconfinement mechanism is identified in partially ionized and
magnetized plasma. Partially ionized plasma is a plasma where some atomic
nuclei retain some of their bound electrons. In partially ionized plasma where
ions can be in different ionization states, these states can be modeled as
fluids of different species for the purpose of describing multi-ion
collisional transport. However, these fluids can transform into one another
via ionization and recombination which is not the case in fully ionized
plasma. As such, transport in partially ionized plasma is different from
transport in fully ionized plasma. In fact, we show here that ionization and
recombination combined with multi-ion cross-field transport lead to a dramatic
increase of ion deconfinement. We also uncover the new physical mechanism
behind it.

Significant progress in understanding of transport in partially ionized
plasmas has been made in the parameter regime of the mix of single-ionized
ions and neutrals. In particular, relevant transport properties were derived
in Ref. Helander  _et al._ (1994); Catto (1994); Hazeltine  _et al._ (1992).
These results became a basis to study tokamak scrape-off layer, such as in
Ref.Feng  _et al._ (2021); Valovi  _et al._ (2004); Zhang and Krasheninnikov
(2020); Casson  _et al._ (2014); McDermott  _et al._ (2021); Angioni and
Helander (2014); Angioni (2021); Angioni  _et al._ (2011); Viezzer  _et al._
(2012). More recently, the same case was also studied in Ref.Meier and Shumlak
(2012).

The case of multiple charge states is less explored in the literature than is
the case of singly ionized plasma and neutrals. More often than not, transport
in such cases is found by using ad-hoc diffusion coefficients which match
experimental observations, such as in Ref.Shurygin (2006); Lee  _et al._
(2022); Odstrcil  _et al._ (2017). In addition to multiple-charge-state
effects playing a role in hot fusion-grade plasmas, particularly in high-Z
impurity transport, as indicated above, these effects may also be expected in
low temperature plasma devices, particularly, mass separation devices
employing partially-ionized magnetized plasmaLehnert (1971); Krishnan (1983);
Geva  _et al._ (1984); Bittencourt and Ludwig (1987); Dolgolenko and Yu. A.
Muromkin (2017); Ochs  _et al._ (2017); Gueroult  _et al._ (2018a); Gueroult
and Fisch (2014); Gueroult  _et al._ (2015, 2018b, 2019); Zweben  _et al._
(2018); Rax and Gueroult (2016). Therefore, there is a need to provide a
first-principle explanation of transport in partially ionized plasma where
ions can be in multiple charge states. This paper identifies the key new
mechanisms at play and shows how large and unanticipated new effects may
occur.

The difficulty of ions with varying charge states to remain in local
thermodynamic equilibrium while respecting the momentum conservation in
collisions between magnetized ions results in two kinds of mechanisms. First,
there is a larger net transport or deconfinement of all ions. Second, the
relative local densities of different charge states may differ from those in
thermal equilibrium, affecting inferences of plasma parameters such as
electron temperature Gregorian  _et al._ (2005a, b). Thus, deviations from
local charge state equilibrium may affect a variety of properties of plasma,
such as the time-dependent radial distributions of the ion velocities Foord
_et al._ (1994), the magnetic field Mikitchuk  _et al._ (2019); Rosenzweig
_et al._ (2020); Davara  _et al._ (1998), the charge-state composition
Gregorian  _et al._ (2003), the electron temperature Gregorian  _et al._
(2005a), the ion temperature Alumot  _et al._ (2019), and the electron
densities Gregorian  _et al._ (2005b); Doron  _et al._ (2016) in Z-pinch
devices. The electron density can be determined from the absolute intensities
of spectral lines Gregorian  _et al._ (2003) and from the ionization times in
the plasma Gregorian  _et al._ (2005b).

If ions are restricted to only two charge states, we prove here that, in fact,
ion densities are in local charge state equilibrium, signifying no impact of
multi-ion transport on spectroscopic inferences. However, the deviation of
relative densities from local plasma equilibrium can be anticipated to appear
in the plasma once ions can be in three or more charge states. As such, it
could be an important extension of the results presented here to estimate the
size of the impact of multi-ion transport in partially ionized plasma on the
densities of constituents. The case of three or more charge states, however,
is out of the scope of this paper; here, we are content to identify the
possibly huge deconfinement effect that is already present in the two ion
charge state case.

The paper is organized as follows. In Sec. II, equilibrium in the two-charge-
state case is considered. In Sec. III, cross-field transport timescale is
derived. Sec. IV summarizes our results and discusses potential applications.

##  II Equilibrium Conditions on Ion Density Profiles.

In the case of ions being in only two charge states the density profiles of
both charge state fluids can be found exactly. The equilibrium density
profiles balance the tension between two types of processes that can change
the number of ions in a given charge state at a given spot: non-local and
local. The non-local process, which is the cross-field transport, is ambipolar
to the leading order, i.e. up to $\mathcal{O}(\sqrt{m_{e}/m_{i}})$ or
$\mathcal{O}(\rho_{i}/L)^{2}$, meaning that if ions can be in charge states
$Z_{a}$ and $Z_{b}$ ($Z_{a}<Z_{b}$) and $\Gamma_{s}=n_{s}\mathbf{u}_{s}$ is
the cross-field particle flux of ions in charge state $s$, then
$Z_{a}\Gamma_{a}+Z_{b}\Gamma_{b}=0$. This comes from the following
observation. In magnetized plasma, force $\mathbf{F}_{s}$ acting on species
$s$ leads to an $\mathbf{F}_{s}\times\mathbf{B}$ drift, particle flux being

$\displaystyle\Gamma_{s}=\frac{n_{s}\mathbf{F}_{s}\times\mathbf{B}}{Z_{s}eB^{2}}.$
(1)

If $\mathbf{F}_{ab}$ is the rate of momentum transfer from fluid $b$ to fluid
$a$, then $\mathbf{F}_{ba}=-\mathbf{F}_{ab}$. Consequently, particle flux due
to interaction between fluids $s$ and $s^{\prime}$, which is

$\displaystyle\Gamma_{ss^{\prime}}=\frac{n_{s}\mathbf{F}_{ss^{\prime}}\times\mathbf{B}}{Z_{s}eB^{2}},$
(2)

obeys the following relation:
$Z_{s}\Gamma_{ss^{\prime}}+Z_{s^{\prime}}\Gamma_{s^{\prime}s}=0$. When plasma
is out of equilibrium, ion-ion friction force is much larger than electron-ion
friction force and viscous forces. Therefore, the corresponding particle flux
is also larger and to the leading order $\Gamma_{a}$ is equal to $\Gamma_{ab}$
while $\Gamma_{b}$ is equal to $\Gamma_{ba}$. It also implies that when plasma
is pushed out of equilibrium due to a change in external forces acting on
plasma, there is timescale separation between ion-ion transport timescale and
slower electron-ion transport timescale, as described in Ref.Mlodik  _et al._
(2020).

Figure 1: Multi-ion cross-field transport in two-charge-state plasma between
locations $x_{1}$ and $x_{2}$. Change in density can happen due to ionization
or recombination (black lines), transport of ions in charge state $Z_{a}$
(blue) and $Z_{b}$ (red). In equilibrium, all the lines should represent the
same rate of density change. However, multi-ion cross-field transport obeys
ambipolarity condition $Z_{a}\Gamma_{a}+Z_{b}\Gamma_{b}=0$. As such,
$\Gamma_{a}=\Gamma_{b}=0$ and $\nu_{i}n_{a}-\nu_{r}n_{b}=0$ in the
equilibrium.

The local processes are ionization and recombination, and they conserve the
number of ions, i.e. $n_{a}+n_{b}$, locally. If $s$ is the net rate of change
of particle density in charge state $b$, then the continuity equation for
fluids $a$ and $b$ can be written in equilibrium as

$\displaystyle\nabla\cdot\Gamma_{a}=-s;$ (3)
$\displaystyle\nabla\cdot\Gamma_{b}=s.$ (4)

They can be combined to get

$\displaystyle\nabla\cdot\left(\Gamma_{a}+\Gamma_{b}\right)=0.$ (5)

Together with the ambipolarity condition $Z_{a}\Gamma_{a}+Z_{b}\Gamma_{b}=0$,
this leads to $\Gamma_{a}=const$. If the boundary condition is no ion flux
through the boundary, then particle flux vanishes everywhere, i.e.
$\Gamma_{a}=\Gamma_{b}=0$. Therefore, both non-local and local processes are
in the state of dynamic equilibrium when plasma itself is in equilibrium as
visualized in Fig. 1. In particular, in equilibrium the rate of ionization and
recombination are the same. Also in equilibrium there is no cross-field
transport due to ion-ion friction. As such, density profiles of charge state
fluids satisfy two constraints. The first is generalized pinch relations which
are satisfied in the plasma with no net cross-field transport due to ion-ion
friction (see Ref.Taylor (1961); Spitzer (1952); Braginskii (1965); Mlodik
_et al._ (2020, 2021); Kolmes  _et al._ (2020, 2021) for information on
generalized pinch relations and Ref.Kolmes  _et al._ (2018) for the
derivation). In particular, in the absence of a temperature gradient, they
take the form

$\displaystyle\big{(}n_{a}e^{\Phi_{a}/T}\big{)}^{1/Z_{a}}\propto\big{(}n_{b}e^{\Phi_{b}/T}\big{)}^{1/Z_{b}},$
(6)

where $\Phi_{s}$ is the potential energy of an ion in charge state $s$. The
second constraint is the local charge state equilibrium, which is ensured by
the rates of ionization and recombination being equal to each other.
Mathematically, local charge state equilibrium can be found in the following
way in two limits. In a dilute plasma the corona model

$\displaystyle\frac{n_{a}}{n_{b}}=f(T)$ (7)

describes local charge state equilibrium. It can be applied when
$10^{12}t_{I}^{-1}<n_{e}<10^{16}T_{e}^{7/2}~{}cm^{-3}$, where
$t_{I}=\left(\alpha_{r}n_{b}\right)^{-1}$ is the ionization time,
$\alpha_{r}=2.7\times 10^{-13}Z_{b}^{2}T_{e}^{-1/2}cm^{3}/s$, $T_{e}$ is in
the units of $eV$, according to the NRL formularyRichardson (2019). In a dense
plasma, the Saha equilibrium

$\displaystyle\frac{n_{b}n_{e}^{Z_{b}-Z_{a}}}{n_{a}}=f(T)$ (8)

describes local charge state equilibrium. Simultaneously, in this case local
charge state equilibrium is equivalent to local thermodynamic equilibrium.
According to the NRL formularyRichardson (2019), electron density is required
to be $n_{e}\gtrsim 7\times
10^{18}Z_{b}^{7}n^{-17/2}(T_{e}/E_{\infty}^{Z})^{1/2}~{}cm^{-3}$ in this case
if initially an ion in charge state $a$ is in the state $n$ and the ionization
energy of that ion is $E_{\infty}^{Z}$. In the case of intermediate density,
there are no simplifying assumptions such as the prevalence of two-body or
three-body recombination, so codes like FLYCHKChung  _et al._ (2005) need to
be used in order to find the relative abundance of charge states.

These constraints describe an interesting set of density profiles. In
thermodynamic equilibrium, when charge state fluids and electron fluid are all
in equilibrium,

$\displaystyle\frac{\nabla
n_{e}}{n_{e}}=-\frac{\nabla\Phi}{T}\cdot\frac{\langle Z\rangle}{\langle
Z(Z+1)\rangle},$ (9)

where $\langle...\rangle$ denotes ion charge state average. In two-charge-
state plasma, $\langle
X\rangle=\left(X_{a}n_{a}+X_{b}n_{b}\right)/\left(n_{a}+n_{b}\right)$. Eqs.
(6) and (9) can be combined with the appropriate local charge state
equilibrium condition to provide an explicit form of the ion density profiles.

##  III Cross-Field Transport in Two-Charge-State Fluid Model.

Suppose that partially ionized plasma is subjected to a change such as an
application of potential $\Phi$ and the goal is to find how plasma reacts to
that change. In order to isolate the ion deconfinement effect and find out the
timescale on which it happens, consider the case when ion-ion transport is
much faster than electron-ion transport (as for why this ordering of transport
mechanisms is desirable to simplify the model, see Appendix A). It is easiest
to see the combined effects of collisional cross-field transport and
ionization and recombination if plasma possesses the following qualities.
First, assume that plasma is an isothermal slab immersed in uniform magnetic
field with all the gradients being perpendicular to the slab boundaries.
Second, assume that there are no neutral particles in the plasma and ions can
be in one of two charge states $a$ and $b$. Third, assume that strength of the
magnetic field is such that plasma $\beta$ is low ($\beta\ll 1$) and ion Hall
parameter $\Omega_{a}/\nu_{ab}$ is large ($\Omega_{a}/\nu_{ab}\gg 1$). Fourth,
assume that the effects of “charge-exchange-like” collisions, ionization and
recombination on transport coefficients are limited to an effective collision
frequency such that momentum exchange rate between charge state fluids $a$ and
$b$ is $\nu_{ab}m_{a}n_{a}\left(\mathbf{u}_{b}-\mathbf{u}_{a}\right)$ where
$\nu_{ab}=\nu_{ab,elastic}+\nu_{ab,inelastic}+\dot{n_{a}}/n_{a}$. Note that in
some cases it is impossible to separate elastic and inelastic parts in the
collision frequency since ultimately it is a process which obeys the principle
of quantum indistinguishabilityKrstic and Schultz (1999). The fourth
assumption also implies that the ions which have changed their identity have
the fluid velocity of the fluid they were a part of before ionization or
recombination. Note, however, that regardless of the details of momentum
transfer between species $a$ and $b$, the particle flux due to ion-ion
collisions is going to relax to zero in equilibrium because ion-ion transport
is ambipolar as long as momentum transfer takes place between the fluids at
the same point. There could be a change of shape of equilibrium density
profiles by analogy to the effect of thermal force on multi-ion transport.
Overall, taking into account the velocity disparity is not going to change the
conclusions of Sec. II and is not going to change the nature of the
conclusions of Sec. III. Fifth, assume that the change of the external
potential $\Phi$ is happening fast compared to ion-ion collisional transport
timescale but slow compared to faster timescales in the system, and the size
of this change is small ($\Delta\Phi\ll T$).

Once all of these assumptions are satisfied, plasma can be treated as a
collection of three fluids, two of them representing ions in charge states $a$
and $b$, and the third representing electrons. Note that the timescale of a
given cross-field transport mechanism is inversely proportional to the size of
particle flux $\Gamma$ generated by it. Therefore, the fastest collisional
transport timescale is the ion-ion transport timescale, which is to be derived
in this Section. Because the electron-ion friction force and viscous force are
much smaller than the ion-ion friction force when plasma is out of
equilibrium, the corresponding timescales will be much longer than the ion-ion
transport timescale. The change of magnetic field can be neglected due to the
low-$\beta$ assumption. Since both charge state fluids are comprised of ions
of the same mass, temperature equilibration between them can be expected to
happen quickly. Therefore both fluids are assumed to have the same
temperature. Also, linearization of fluid equations around the equilibrium
found in Sec. II is possible.

Despite the fact that the list of necessary assumptions is long, this set of
assumptions describes common situations or close approximations to common
situations. In particular, in a typical tokamak plasma with tungsten
impurities all of the conditions aside from the transport timescale and the
assumption that there are only two charge states are satisfied. Moreover, the
model can accurately describe collisional transport in helium plasma if the
assumptions listed in this Section are satisfied. For example, if parameters
similar to an LAPD discharge aside from ion temperature are taken, i.e.
$n_{e}=10^{12}~{}cm^{-3}$, $T_{i}=7~{}eV$, $B=500~{}G$, then gyroradius of He+
is $\rho_{a}\approx 1.08~{}cm$, gyrofrequency of He+ is $\Omega_{a}\approx
1.2\cdot 10^{6}~{}s^{-1}$, collision frequency of elastic collisions between
He+ and He2+ is $\nu_{ab}\approx 1.8\cdot 10^{4}~{}s^{-1}$, plasma beta for
ions is $2\mu_{0}\left(n_{a}+n_{b}\right)T_{i}/B^{2}\approx 1.2\cdot 10^{-3}$,
charge state abundances, according to corona equilibrium, are around $50\%$.
The net ion charge transport timescale, which is derived later in this
Section, becomes $\tau_{ci}\approx
23\cdot\left(\nu_{ab,elastic}/\nu_{ab}\right)~{}ms$ for a cylindrical plasma
with radius $r=30~{}cm$, which is comparable to the lifetime of an LAPD
discharge, which is on the order of $\sim 10~{}ms$.

The following notation is assumed: ions of mass $m$ can be either in charge
state $Z_{a}$ or $Z_{b}$, density of ions in these states is $n_{a}$ and
$n_{b}$, respectively. All the equations are linearized around the
equilibrium, which is comprised of local charge state equilibrium and
generalized pinch relations.

Figure 2: Transport in fully ionized and partially ionized plasma. In fully
ionized plasma (i) imbalance in flow velocities leads to ion-ion friction and
corresponding cross-field transport. Due to this process ion densities adjust
as to relax the difference of flow velocities. In partially ionized plasma
(ii), once ions in a charge state move across field lines due to ion-ion
friction, they ionize or recombine in order to maintain local charge state
equilibrium. Consequently, densities of ion fluids are no longer independent.
Because the drift velocities are thus constrained, the plasma must move as a
whole.

Assume that heat transport across the system is infinitely fast, such that
plasma is isothermal. Also assume that density perturbation is small. With
these assumptions, the momentum equation for particles of type $s$ is

$\displaystyle
m_{s}\frac{d_{s}\mathbf{u}_{s}}{dt}=q_{s}\mathbf{E}+q_{s}\mathbf{u}_{s}\times\mathbf{B}-T\frac{\nabla
n_{s}}{n_{s}}+\frac{\sum_{s^{\prime}}\mathbf{R}_{ss^{\prime}}}{n_{s}}+\mathbf{F}_{s}.$
(10)

Here $\mathbf{u}_{s}$ is the flow velocity, $d_{s}/dt=\partial/\partial
t+\mathbf{u}_{s}\cdot\nabla$, $m_{s}$ is the mass, $q_{s}$ is the charge, and
$\mathbf{R}_{ss^{\prime}}=\nu_{ss^{\prime}}m_{s}n_{s}(\mathbf{u}_{s^{\prime}}-\mathbf{u}_{s})$
is the friction force (note that thermal friction is zero under the isothermal
plasma assumption), $\mathbf{F}_{s}=-\nabla\Phi_{s}$ is the external force
acting on species $s$. It is important to note that the collision frequency
$\nu_{ss^{\prime}}$ in Eq. (10) includes all processes that participate in the
momentum transfer between fluids $s$ and $s^{\prime}$. Eq. (10) can be
rewritten as

$\displaystyle\mathbf{u}_{s}\times\hat{\mathbf{b}}=-\frac{\mathbf{E}}{B}+\frac{1}{\Omega_{s}}\frac{d_{s}\mathbf{u}_{s}}{dt}+\frac{T\nabla
n_{s}}{m_{s}n_{s}\Omega_{s}}-\frac{\sum_{s^{\prime}}\mathbf{R}_{ss^{\prime}}}{m_{s}n_{s}\Omega_{s}}-\frac{\mathbf{F}_{s}}{m_{s}\Omega_{s}}.$
(11)

The continuity equation for particles in charge state $s$ is

$\displaystyle\frac{\partial n_{s}}{\partial
t}+\nabla\cdot\Gamma_{s}=\left(\frac{\partial n_{s}}{\partial
t}\right)_{i/r}.$ (12)

Here $\left(\partial n_{s}/\partial t\right)_{i/r}$ is the change in density
of particles in charge state $s$ due to ionization and recombination. Eqs.
(11) and (12) are linearized around final global equilibrium which is
described in Sec. II. This final equilibrium includes the effect of an
external force $\mathbf{F}_{s}$ on density profiles. The density of fluid $s$
is split into two components $n_{s}+\widetilde{n}_{s}$ where $n_{s}$ is
density in equilibrium and $\widetilde{n}_{s}$ is the difference between the
density at the given moment of time and the density in the equilibrium. Under
the assumptions made in this section $\widetilde{n}_{s}/n_{s}\ll 1$
everywhere. All the other quantities are split into two parts, equilibrium
value and perturbation, in the same way, and the perturbed part is denoted by
a tilde. Under the assumption that the external potential is small, i.e.
$\Phi/T\ll 1$, density gradients are also small $L\nabla n_{s}/n_{s}\ll 1$ so
the linearized continuity equations are

$\displaystyle\frac{\partial\widetilde{n}_{a}}{\partial
t}=-\nabla\cdot\widetilde{\Gamma}_{ab}+\nu_{r}\widetilde{n}_{b}-\nu_{i}\widetilde{n}_{a},$
(13) $\displaystyle\frac{\partial\widetilde{n}_{b}}{\partial
t}=-\nabla\cdot\widetilde{\Gamma}_{ba}+\nu_{i}\widetilde{n}_{a}-\nu_{r}\widetilde{n}_{b}.$
(14)

Here $\nu_{i}$ is the effective ionization rate and $\nu_{r}$ is the effective
recombination rate, both of which are calculated as linear terms of
$\left(\partial n_{s}/\partial t\right)_{i/r}$ in density perturbation, taken
at equilibrium. $\widetilde{\Gamma}_{ss^{\prime}}$ is particle flux of fluid
$s$ due to momentum transfer with fluid $s^{\prime}$. It is the only component
of particle flux that can have non-zero divergence in low-$\beta$ plasma on
the ion-ion transport timescale. Given that in the equilibrium
$\mathbf{u}_{a}=\mathbf{u}_{b}$, $\widetilde{\Gamma}_{ab}$ is

$\displaystyle\widetilde{\Gamma}_{ab}=\frac{\nu_{ab}}{\Omega_{a}}n_{a}\left(\widetilde{\mathbf{u}}_{b}-\widetilde{\mathbf{u}}_{a}\right)\times\hat{b}.$
(15)

Here $\nu_{ab}$ is momentum transfer frequency between fluids $a$ and $b$.
Note that $\nu_{ab}$ is higher than it would be if $a$ and $b$ were ions of
different species because momentum can be transferred due to ionization or
recombination, as well as due to collisions similar to “charge-exchange” which
are resonant since $a$ and $b$ are ions of same species.
$\widetilde{\Gamma}_{ba}=-Z_{a}\widetilde{\Gamma}_{ab}/Z_{b}$ since ion-ion
transport is ambipolar.

In the projection of Eq. (11) on the direction of gradients the terms
$\propto\mathbf{F}_{s}$ and $\propto\mathbf{E}$ are included in the leading
order, momentum transfer term $\propto\mathbf{R}_{ss^{\prime}}$ is small,
polarization drift term $\propto 1/\Omega_{s}\cdot d_{s}\mathbf{u}_{s}/dt$ is
small. Therefore in the first order the only remaining term is diamagnetic
drift term. The perturbation of the drift velocity up to the first order is

$\displaystyle\widetilde{\mathbf{u}}_{s}=-\frac{\nabla\widetilde{n}_{s}\times\hat{b}}{n_{s}}\frac{T}{m_{s}\Omega_{s}}.$
(16)

Linearized continuity equations become

$\displaystyle\frac{\partial\widetilde{n}_{a}}{\partial
t}=\frac{\nu_{ab}T}{m_{a}\Omega_{a}^{2}}\nabla^{2}\widetilde{n}_{a}-\frac{\nu_{ba}T}{m_{b}\Omega_{a}\Omega_{b}}\nabla^{2}\widetilde{n}_{b}+\nu_{r}\widetilde{n}_{b}-\nu_{i}\widetilde{n}_{a},$
(17) $\displaystyle\frac{\partial\widetilde{n}_{b}}{\partial
t}=\frac{\nu_{ba}T}{m_{b}\Omega_{b}^{2}}\nabla^{2}\widetilde{n}_{b}-\frac{\nu_{ab}T}{m_{a}\Omega_{a}\Omega_{b}}\nabla^{2}\widetilde{n}_{a}+\nu_{i}\widetilde{n}_{a}-\nu_{r}\widetilde{n}_{b}.$
(18)

Eqs. (17) and (18) can be solved by spatial spectral decomposition and then
individually for each eigenmode of Laplace operator $\nabla^{2}$ in the domain
of interest. For any particular eigenmode of Laplace operator which has
eigenvalue of $-k^{2}$ (for example, in slab geometry $N$-th mode has
eigenvalue $k_{N}=\pi N/L$), a substitution $\nabla^{2}\to-k^{2}$ can be done
so the equations for $\widetilde{n}_{a}$ and $\widetilde{n}_{b}$ become

$\displaystyle\frac{\partial\widetilde{n}_{a}}{\partial
t}=-\left(\frac{\nu_{ab}T}{m_{a}\Omega_{a}^{2}}k^{2}+\nu_{i}\right)\widetilde{n}_{a}+\left(\frac{\nu_{ba}T}{m_{b}\Omega_{a}\Omega_{b}}k^{2}+\nu_{r}\right)\widetilde{n}_{b},$
(19) $\displaystyle\frac{\partial\widetilde{n}_{b}}{\partial
t}=-\left(\frac{\nu_{ba}T}{m_{b}\Omega_{b}^{2}}k^{2}+\nu_{r}\right)\widetilde{n}_{b}+\left(\frac{\nu_{ab}T}{m_{a}\Omega_{a}\Omega_{b}}k^{2}+\nu_{i}\right)\widetilde{n}_{a}.$
(20)

Eqs. (19) and (20) form a system of two coupled linear ODEs. As such, they
have two eigenvalues $\lambda$ which correspond to two exponentially changing
solutions ($\propto\exp(\lambda t)$). In particular, if Eqs. (19) and (20) are
rewritten in the form

$\displaystyle\frac{\partial\widetilde{n}_{a}}{\partial
t}=-A_{aa}\widetilde{n}_{a}+A_{ab}\widetilde{n}_{b},$ (21)
$\displaystyle\frac{\partial\widetilde{n}_{b}}{\partial
t}=A_{ba}\widetilde{n}_{a}-A_{bb}\widetilde{n}_{b},$ (22)

then rates of change of density perturbations are

$\displaystyle\lambda_{\pm}=-\frac{A_{aa}+A_{bb}}{2}\pm\sqrt{\left(\frac{A_{aa}+A_{bb}}{2}\right)^{2}+A_{ab}A_{ba}-A_{aa}A_{bb}}.$
(23)

Momentum conservation implies $m_{a}n_{a}\nu_{ab}=m_{b}n_{b}\nu_{ba}$. Charge
state equilibrium implies $\nu_{i}n_{a}=\nu_{r}n_{b}$, although this identity
is exact when $\nu_{i}$ and $\nu_{r}$ are actual ionization and recombination
rate, respectively, while we define them as effective rates. Note that this
identity still stands in many cases when plasma is close to charge state
equilibrium. Given these identities, the following expressions can be
obtained:

$\displaystyle
A_{aa}+A_{bb}=\left(\frac{\nu_{ab}T}{m_{a}\Omega_{a}^{2}}+\frac{\nu_{ba}T}{m_{b}\Omega_{b}^{2}}\right)k^{2}+\nu_{i}+\nu_{r},$
(24)

$\displaystyle
A_{ab}A_{ba}-A_{aa}A_{bb}=-\frac{k^{2}T}{m}\nu_{r}\nu_{ab}\left(\frac{1}{\Omega_{b}}-\frac{1}{\Omega_{a}}\right)^{2}.$
(25)

A useful corollary of these expressions is that $\lambda_{\pm}\leq 0$ in all
cases, and $\lambda_{\pm}=0$ only if $\nu_{i}=\nu_{r}=0$ (no ionization or
recombination). Therefore, the densities of individual plasma components reach
equilibrium and there is no instability. Plasma has two modes that decay to
equilbrium at the rates $|\lambda_{-}|$ and $|\lambda_{+}|$.

In the limiting case $\nu_{i}=\nu_{r}=0$ (when $a$ and $b$ are different
species): $\lambda_{+}=0$,

$\displaystyle\lambda_{-}=-k^{2}\left(\frac{\nu_{ab}T}{m_{a}\Omega_{a}^{2}}+\frac{\nu_{ba}T}{m_{b}\Omega_{b}^{2}}\right),$
(26)

which coincides with Ref.Mlodik  _et al._ (2021). $\lambda_{-}$ has the
meaning of the inverse of multi-ion collisional transport timescale in this
case.

Another limiting case is a sufficiently large system, such that

$\displaystyle\left(\frac{\nu_{ab}T}{m_{a}\Omega_{a}^{2}}+\frac{\nu_{ba}T}{m_{b}\Omega_{b}^{2}}\right)k^{2}\ll\nu_{i}+\nu_{r}.$
(27)

Physically that means that ionization and recombination happen much faster
than collisional cross-field transport. Equivalently, this is the case of the
infinite, homogeneous plasma ($k\to 0$) where transport is absent. Then

$\displaystyle\lambda_{-}\approx-\left(\nu_{i}+\nu_{r}\right),$ (28)
$\displaystyle\lambda_{+}\approx-\frac{k^{2}T}{m}\frac{\nu_{r}}{\nu_{i}+\nu_{r}}\nu_{ab}\left(\frac{1}{\Omega_{b}}-\frac{1}{\Omega_{a}}\right)^{2}.$
(29)

In this case, both eigenvalues correspond to processes with a clear physical
meaning. To see that, note that if Eq. (27) is satisfied and
$\lambda=\lambda_{-}$ is used in Eqs. (19) and (20), then

$\displaystyle\frac{\widetilde{n}_{a}+\widetilde{n}_{b}}{\widetilde{n}_{a}}=\frac{\lambda_{-}+A_{aa}+A_{ab}}{A_{ab}}=\mathcal{O}\left(\frac{k^{2}T}{m_{a}\Omega_{a}^{2}}\frac{\nu_{ab}}{\nu_{i}+\nu_{r}}\right).$
(30)

It can be seen that $\lambda_{-}$ represents local charge state equilibration,
while $\lambda_{+}$ represents global equilibration associated with cross-
field transport. Since $|\lambda_{-}|\gg|\lambda_{+}|$, plasma quickly
approaches local charge state equilibrium everywhere, and then adjusts total
densities in order to reach global equilibrium. Another corollary of this
statement is that spectroscopic inferences from the ratio of densities remain
the same even if multi-ion transport is included as long as ions can be only
in two charge states and plasma obeys condition set in Eq. (27).

Ion charge transport across magnetic field lines is characterized by the
change in total ion charge
$\widetilde{\rho}_{ci}=Z_{a}\widetilde{n}_{a}+Z_{b}\widetilde{n}_{b}$ in the
slower-varying mode. To see the size of ion charge transport, consider
$\lambda=\lambda_{+}$ so

$\displaystyle\frac{Z_{a}\widetilde{n}_{a}+Z_{b}\widetilde{n}_{b}}{\widetilde{n}_{a}}=\frac{Z_{a}\lambda_{+}+Z_{a}A_{aa}+Z_{b}A_{ab}}{A_{ab}}.$
(31)

In the large system where Eq. (27) is satisfied,

$\displaystyle\frac{Z_{a}\widetilde{n}_{a}+Z_{b}\widetilde{n}_{b}}{\widetilde{n}_{a}}\approx
Z_{a}\frac{\nu_{i}}{\nu_{r}}+Z_{b}.$ (32)

This highlights a major difference in the nature of transport in partially
ionized plasma. Ion charge density is no longer conserved in such plasma even
if cross-field transport happens due to ion-ion friction. This ion
deconfinement happens on the timescale

$\displaystyle\tau_{ci}=\frac{\widetilde{\rho}_{ci}}{\partial\widetilde{\rho}_{ci}/\partial
t}=|\lambda_{+}|^{-1}.$ (33)

Alternatively it can be written as

$\displaystyle\tau_{ci}=\frac{e^{2}B^{2}}{k^{2}Tm}\frac{Z_{a}^{2}Z_{b}^{2}}{\left(Z_{b}-Z_{a}\right)^{2}}\frac{n_{a}+n_{b}}{n_{a}\nu_{ab}}\left[1+\mathcal{O}\left(\frac{k^{2}T}{m_{a}\Omega_{a}^{2}}\frac{\nu_{ab}}{\nu_{i}+\nu_{r}}\right)\right].$
(34)

Net ion charge moves across field lines dramatically
($\mathcal{O}(\sqrt{m_{i}/m_{e}})$) faster in partially ionized plasma
compared to fully ionized plasma. To see it, note that a lower bound on the
electron-ion transport timescale in fully ionized plasma is
$\tau_{ie,tr}=\left(\tau_{ae,tr}^{-1}+\tau_{be,tr}^{-1}\right)^{-1}$, where
$\tau_{ae,tr}$ and $\tau_{be,tr}$ can be found, following Ref.Mlodik  _et al._
(2021), as RHS of Eq. (26) with substitutions $b\to e$ and $a\to e$,
respectively. After some simplifications using the assumptions $m_{a}=m_{b}$
and $m_{e}/m_{a}\ll 1$,

$\displaystyle\frac{\tau_{ci}}{\tau_{ie,tr}}=\left(\frac{m_{a}}{2m_{e}}\right)^{1/2}\frac{\nu_{ab}}{\nu_{ab,elastic}}\frac{n_{a}n_{b}\left(Z_{b}-Z_{a}\right)^{2}}{\left(2n_{e}+n_{a}Z_{a}^{2}+n_{b}Z_{b}^{2}\right)\left(n_{a}+n_{b}\right)}.$
(35)

Note that Eq. (35) describes the ratio of net ion charge transport timescales
in partially ionized plasma and fully ionized plasma. Unlike in fully ionized
plasma where multi-ion transport leads to ion stratification while conserving
local ion charge densityKolmes  _et al._ (2018); Mlodik  _et al._ (2020), in
partially ionized plasma multi-ion transport leads to plasma moving across
magnetic field lines as a whole. Moreover, this result relies only on the
existence of an ion density perturbation in partially ionized plasma. As such,
ion deconfinement is going to happen whenever the ion density is out of
equilibrium. For example, if plasma consists of background ion species and an
impurity species which is present in multiple charge states, then, if the
background ion species moves across field lines, the impurity is going to move
across magnetic field lines due to the interplay between ionization and
collisions between ions in different charge states.

##  IV Discussion.

Partial-ionization deconfinement effect in magnetized plasma has been
identified. Partially ionized plasma has both cross-field transport due to
ion-ion friction and ionization and recombination. Combined, these processes
result in plasma as a whole being deconfined on multi-ion transport timescale
as opposed to fully ionized plasma which is deconfined only on electron-ion
transport timescale, $\mathcal{O}(\sqrt{m_{e}/m_{i}})$ slower. Virtually all
plasmas that include high-$Z$ species are not fully ionized but much of the
existing analytic transport theory does not include the possibility of
transition between charge states. This paper describes the physical phenomena
occurring due to these transitions using first-principles.

The main differences in cross-field transport in partially ionized plasma
compared to fully ionized plasma are twofold. First, charge state fluids can
exchange momentum between each other not only due to the usual Coulomb
collisions but also due to charge-exchange collisions and pick-up current.
Second, partially ionized plasma has fewer degrees of freedom to relax the
imbalance in flow velocities of fluids which comprise the plasma due to local
charge equilibrium being enforced by ionization and recombination which
results in extra transport. Perhaps the most important difference between
transport in partially ionized plasma and fully ionized plasma is that
movement of net ion charge in low-$\beta$ plasma becomes possible even on the
ion-ion transport timescale, as opposed to fully ionized plasma where it
happens on the electron-ion timescale. This can result in dramatic
deconfinement whenever there is a change in external force, such as
centrifugal force like in plasma mass filters, acting on plasma. While the
deconfinement effect identified here may be the most dramatic effect, it is
also noteworthy that there is a charge transport effect that results in net
current across field lines, an effect that in and of itself may be quite
significant because, in general, so few processes lead to such currents.

In this paper the simplest model is considered in order to isolate ion
deconfinement effect. However, there are many ways in which it can be
extended. While in two-charge-state case plasma reaches equilibrium when there
is local charge state equilibrium and no ion-ion cross-field transport, that
would not necessarily be the case if ions can be in three or more charge
states. Another extension is adding temperature gradient to the model. Yet
another extension could be a problem of background species and partially
ionized species. This type of analysis could be important for understanding of
high-$Z$ impurity transport.

###### Acknowledgements.

The authors thank R. Gueroult, S. Davidovits, E. Mitra, and J. Waybright for
useful conversations. This work was supported by Cornell NNSA 83228-10966
[Prime No. DOE (NNSA) DE-NA0003764] and by DOE DE-SC0016072.

## Data Availability Statement

Data sharing is not applicable to this article as no new data were created or
analyzed in this study.

##  Appendix A Electron-Ion Transport Considerations.

Electron-ion transport and viscous transport are not included in the model
considered in this paper in order to isolate the deconfinement effect and not
focus on the full dynamics which are complicated due to the following reason.
The number of the pathways to change the density of individual charge states
limits the problems that can be solved analytically. Consider continuity
equation in multispecies plasma:

$\displaystyle\frac{\partial n_{s}}{\partial
t}=-\nabla\cdot\sum_{s^{\prime}}\Gamma_{ss^{\prime}}+s_{s}.$ (36)

Here $\Gamma_{ss^{\prime}}$ is particle flux of species $s$ due to collisions
with species $s^{\prime}$, $s_{s}$ is the source term due to
ionization/recombination, different charge states are treated as separate
species. Suppose that all ions are of the same chemical element and can be in
two charge states $a$ and $b$. If $s$ is the local net rate of ionization and
recombination,

$\displaystyle\frac{\partial n_{a}}{\partial
t}=-\nabla\cdot\left(\Gamma_{aa}+\Gamma_{ab}+\Gamma_{ae}\right)-s.$ (37)
$\displaystyle\frac{\partial n_{b}}{\partial
t}=-\nabla\cdot\left(\Gamma_{ba}+\Gamma_{bb}+\Gamma_{be}\right)+s.$ (38)
$\displaystyle\frac{\partial n_{e}}{\partial
t}=-\nabla\cdot\left(\Gamma_{ea}+\Gamma_{eb}+\Gamma_{ee}\right)+s\left(Z_{b}-Z_{a}\right).$
(39)

Clearly, if viscous transport is significant ($\Gamma_{ss}$ is not
subdominant), then the problem already becomes complicated. However, since
cross-field viscosity is an FLR effect there can be some plasmas where viscous
transport is small compared to the frictional transport, e.g. in Ref.Mlodik
_et al._ (2020). Then

$\displaystyle\frac{\partial n_{a}}{\partial
t}=-\nabla\cdot\left(\Gamma_{ab}+\Gamma_{ae}\right)-s.$ (40)
$\displaystyle\frac{\partial n_{b}}{\partial
t}=-\nabla\cdot\left(\Gamma_{ba}+\Gamma_{be}\right)+s.$ (41)
$\displaystyle\frac{\partial n_{e}}{\partial
t}=-\nabla\cdot\left(\Gamma_{ea}+\Gamma_{eb}\right)+s\left(Z_{a}-Z_{b}\right).$
(42)

However, there is a degree of freedom remaining, since there are three
constraints on four variables ($\Gamma_{ab}$, $\Gamma_{ae}$, $\Gamma_{be}$,
$s$). Therefore, even in the presence of the electron-ion frictional
transport, individual mechanisms are no longer required to be in the detailed
balance (i.e. they don’t have to vanish individually in the equilibrium). In
turn this means that there can be conveyor-belt-type equilibrium: ionization
at the one end of the system, recombination at the other end of the system,
and collisional cross-field transport to move ions between those ends. As
such, if electron-ion transport is not negligible compared to ion-ion
transport then additional physics may arise, thus obscuring the ion
deconfinement effect.

## References

  * Helander  _et al._ (1994) P. Helander, S. I. Krasheninnikov, and P. J. Catto, Phys. Plasmas 1, 3174 (1994).
  * Catto (1994) P. J. Catto, Phys. Plasmas 1, 1936 (1994).
  * Hazeltine  _et al._ (1992) R. D. Hazeltine, M. D. Calvin, P. M. Valanju, and E. R. Solano, Nuclear fusion 32, 3 (1992).
  * Feng  _et al._ (2021) Y. Feng, M. Jakubowski, R. König, M. Krychowiak, M. Otte, F. Reimold, D. Reiter, O. Schmitz, D. Zhang, C. Beidler, C. Biedermann, S. Bozhenkov, K. Brunner, A. Dinklage, P. Drewelow, F. Effenberg, M. Endler, G. Fuchert, Y. Gao, J. Geiger, K. Hammond, P. Helander, C. Killer, J. Knauer, T. Kremeyer, E. Pasch, L. Rudischhauser, G. Schlisio, T. S. Pedersen, U. Wenzel, V. Winters, and W7-X team, Nuclear Fusion 61, 086012 (2021).
  * Valovi  _et al._ (2004) M. Valovi, R. Budny, L. Garzotti, X. Garbet, A. A. Korotkov, J. Rapp, R. Neu, O. Sauter, P. deVries, B. Alper, M. Beurskens, J. Brzozowski, D. McDonald, H. Leggate, C. Giroud, V. Parail, I. Voitsekhovitch, and JET EFDA contributors, Plasma Physics and Controlled Fusion 46, 1877 (2004).
  * Zhang and Krasheninnikov (2020) Y. Zhang and S. I. Krasheninnikov, Plasma Physics and Controlled Fusion 62, 115018 (2020).
  * Casson  _et al._ (2014) F. J. Casson, C. Angioni, E. A. Belli, R. Bilato, P. Mantica, T. Odstrcil, T. Pütterich, M. Valisa, L. Garzotti, C. Giroud, J. Hobirk, C. F. Maggi, J. Mlynar, and M. L. Reinke, Plasma Physics and Controlled Fusion 57, 014031 (2014).
  * McDermott  _et al._ (2021) R. McDermott, C. Angioni, M. Cavedon, A. Kappatou, R. Dux, R. Fischer, P. Manas, and the ASDEX Upgrade Team, Nuclear Fusion 62, 026006 (2021).
  * Angioni and Helander (2014) C. Angioni and P. Helander, Plasma Phys. Control. Fusion 56, 124001 (2014).
  * Angioni (2021) C. Angioni, Plasma Physics and Controlled Fusion 63, 073001 (2021).
  * Angioni  _et al._ (2011) C. Angioni, R. M. McDermott, F. J. Casson, E. Fable, A. Bottino, R. Dux, R. Fischer, Y. Podoba, T. Pütterich, F. Ryter, and E. Viezzer (ASDEX Upgrade Team), Phys. Rev. Lett. 107, 215003 (2011).
  * Viezzer  _et al._ (2012) E. Viezzer, T. Pütterich, R. Dux, and R. M. McDermott, Review of Scientific Instruments 83, 103501 (2012).
  * Meier and Shumlak (2012) E. T. Meier and U. Shumlak, Physics of Plasmas 19, 072508 (2012).
  * Shurygin (2006) V. A. Shurygin, Phys. Plasmas 13, 082506 (2006).
  * Lee  _et al._ (2022) H. Lee, H. Lee, Y. S. Han, J. Song, E. A. Belli, W. Choe, J. Kang, J. Lee, J. Candy, and J. Lee, Physics of Plasmas 29, 022504 (2022).
  * Odstrcil  _et al._ (2017) T. Odstrcil, T. Pütterich, C. Angioni, R. Bilato, A. Gude, M. Odstrcil, ASDEX Upgrade Team, and the EUROfusion MST1 team, Plasma Physics and Controlled Fusion 60, 014003 (2017).
  * Lehnert (1971) B. Lehnert, Nucl. Fusion 11, 485 (1971).
  * Krishnan (1983) M. Krishnan, Phys. Fluids 26, 2676 (1983).
  * Geva  _et al._ (1984) M. Geva, M. Krishnan, and J. L. Hirshfield, J. Appl. Phys. 56, 1398 (1984).
  * Bittencourt and Ludwig (1987) J. A. Bittencourt and G. O. Ludwig, Plasma Phys. Control. Fusion 29, 601 (1987).
  * Dolgolenko and Yu. A. Muromkin (2017) D. A. Dolgolenko and Yu. A. Muromkin, Phys.-Usp. 60, 994 (2017).
  * Ochs  _et al._ (2017) I. E. Ochs, R. Gueroult, N. J. Fisch, and S. J. Zweben, Phys. Plasmas 24, 043503 (2017).
  * Gueroult  _et al._ (2018a) R. Gueroult, J.-M. Rax, S. J. Zweben, and N. J. Fisch, Plasma Phys. and Control. Fusion 60, 014018 (2018a).
  * Gueroult and Fisch (2014) R. Gueroult and N. J. Fisch, Plasma Sources Sci. Technol. 23, 035002 (2014).
  * Gueroult  _et al._ (2015) R. Gueroult, D. T. Hobbs, and N. J. Fisch, J. Hazard. Mater. 297, 153 (2015).
  * Gueroult  _et al._ (2018b) R. Gueroult, J.-M. Rax, and N. J. Fisch, J. Clean. Prod. 182, 1060 (2018b).
  * Gueroult  _et al._ (2019) R. Gueroult, S. J. Zweben, N. J. Fisch, and J.-M. Rax, Phys. Plasmas 26, 043511 (2019).
  * Zweben  _et al._ (2018) S. J. Zweben, R. Gueroult, and N. J. Fisch, Physics of Plasmas 25, 090901 (2018).
  * Rax and Gueroult (2016) J.-M. Rax and R. Gueroult, J. Plasma Phys. 82, 595820504 (2016).
  * Gregorian  _et al._ (2005a) L. Gregorian, E. Kroupp, G. Davara, A. Starobinets, V. I. Fisher, V. A. Bernshtam, Y. V. Ralchenko, Y. Maron, A. Fisher, and D. H. H. Hoffmann, Phys. Rev. E 71, 056402 (2005a).
  * Gregorian  _et al._ (2005b) L. Gregorian, E. Kroupp, G. Davara, V. I. Fisher, A. Starobinets, V. A. Bernshtam, A. Fisher, and Y. Maron, Phys. Plasmas 12, 092704 (2005b).
  * Foord  _et al._ (1994) M. E. Foord, Y. Maron, G. Davara, L. Gregorian, and A. Fisher, Phys. Rev. Lett. 72, 3827 (1994).
  * Mikitchuk  _et al._ (2019) D. Mikitchuk, M. Cvejić, R. Doron, E. Kroupp, C. Stollberg, Y. Maron, A. L. Velikovich, N. D. Ouart, J. L. Giuliani, T. A. Mehlhorn, and E. Yu, Phys. Rev. Lett. 122, 045001 (2019).
  * Rosenzweig  _et al._ (2020) G. Rosenzweig, E. Kroupp, T. Queller, A. Starobinets, Y. Maron, V. Tangri, J. L. Giuliani, and A. Fruchtman, Phys. Plasmas 27, 022705 (2020).
  * Davara  _et al._ (1998) G. Davara, L. Gregorian, E. Kroupp, and Y. Maron, Phys. Plasmas 5, 1068 (1998).
  * Gregorian  _et al._ (2003) L. Gregorian, V. A. Bernshtam, E. Kroupp, G. Davara, and Y. Maron, Phys. Rev. E 67, 016404 (2003).
  * Alumot  _et al._ (2019) D. Alumot, E. Kroupp, E. Stambulchik, A. Starobinets, I. Uschmann, and Y. Maron, Phys. Rev. Lett. 122, 095001 (2019).
  * Doron  _et al._ (2016) R. Doron, B. Rubinstein, J. Citrin, R. Arad, Y. Maron, A. Fruchtman, H. R. Strauss, and T. A. Mehlhorn, Phys. Plasmas 23, 122126 (2016).
  * Mlodik  _et al._ (2020) M. E. Mlodik, E. J. Kolmes, I. E. Ochs, and N. J. Fisch, Phys. Rev. E 102, 013212 (2020).
  * Taylor (1961) J. B. Taylor, Phys. Fluids 4, 1142 (1961).
  * Spitzer (1952) L. Spitzer, Astrophys. J. 116, 299 (1952).
  * Braginskii (1965) S. I. Braginskii, “Transport processes in a plasma,” in  _Reviews of Plasma Physics_ , Vol. 1, edited by M. A. Leontovich (Consultants Bureau, New York, 1965) p. 205.
  * Mlodik  _et al._ (2021) M. E. Mlodik, E. J. Kolmes, I. E. Ochs, and N. J. Fisch, Physics of Plasmas 28, 052702 (2021).
  * Kolmes  _et al._ (2020) E. J. Kolmes, I. E. Ochs, M. E. Mlodik, and N. J. Fisch, Phys. Lett. A , 126262 (2020).
  * Kolmes  _et al._ (2021) E. J. Kolmes, I. E. Ochs, and N. J. Fisch, Comp. Phys. Communications 258, 107511 (2021).
  * Kolmes  _et al._ (2018) E. J. Kolmes, I. E. Ochs, and N. J. Fisch, Phys. Plasmas 25, 032508 (2018).
  * Richardson (2019) A. S. Richardson, “2019 NRL plasma formulary,” (2019).
  * Chung  _et al._ (2005) H.-K. Chung, M. Chen, W. Morgan, Y. Ralchenko, and R. Lee, High Energy Density Physics 1, 3 (2005).
  * Krstic and Schultz (1999) P. S. Krstic and D. R. Schultz, Journal of Physics B: Atomic, Molecular and Optical Physics 32, 3485 (1999).