# Statistics of Turbulence in the Solar Wind. I.
What is the Reynolds Number of the Solar Wind?

Daniel Wrench Victoria University of Wellington
Kelburn, Wellington 6012, New Zealand  Tulasi N. Parashar Victoria University
of Wellington
Kelburn, Wellington 6012, New Zealand  Sean Oughton University of Waikato,
Hamilton 3240, New Zealand  Kevin de Lange  Victoria University of Wellington
Kelburn, Wellington 6012, New Zealand  Marcus Frean Victoria University of
Wellington
Kelburn, Wellington 6012, New Zealand

###### Abstract

The Reynolds number, $Re$, is an important quantity for describing a turbulent
flow. It tells us about the bandwidth over which energy can cascade from large
scales to smaller ones, prior to the onset of dissipation. However,
calculating it for nearly collisionless plasmas like the solar wind is
challenging. Previous studies have used “effective” Reynolds number
formulations, expressing $Re$ as a function of the correlation scale and
either the Taylor scale or a proxy for the dissipation scale. We find that the
Taylor scale definition of the Reynolds number has a sizeable prefactor of
approximately 27, which has not been employed in previous works. Drawing from
18 years of data from the Wind spacecraft at 1 au, we calculate the magnetic
Taylor scale directly and use both the ion inertial length and the magnetic
spectrum break scale as approximations for the dissipation scale, yielding
three distinct $Re$ estimates for each 12-hour interval. Average values of
$Re$ range between 116,000 and 3,406,000, within the general distribution of
past work. We also find considerable disagreement between the methods, with
linear associations of between 0.38 and 0.72. Although the Taylor scale method
is arguably more physically motivated, due to its dependence on the energy
cascade rate, more theoretical work is needed in order to identify the most
appropriate way of calculating effective Reynolds numbers for kinetic plasmas.
As a summary of our observational analysis, we make available a data product
of 28 years of 1 au solar wind and magnetospheric plasma measurements from
Wind.

##  1 Introduction

Most naturally-occurring plasmas are either observed to be, or believed to be,
in a turbulent state. There is significant variation in the parameters of
these systems, including the length and time scales, the plasma $\beta$, the
turbulent Mach numbers, and the relative size of the system compared to
kinetic scales. Many of these systems are in what is called a “kinetic” state,
where the dynamical length and time scales of interest are comparable to or
smaller than the collisional time scales of interest. Astrophysical examples
include the solar wind (e.g., Bruno & Carbone, 2013), accretion disks (e.g.,
Balbus & Hawley, 1998), and the intracluster medium (e.g., Mohapatra et al.,
2020). For these systems, the collisional closures associated with fluid
models are no longer applicable (or at least not obviously so). This means one
has to resort to higher-order closures for the fluid models, or, in most
cases, to a kinetic description of the plasma (e.g., Marsch, 2006).

Turbulence theories utilize dimensionless parameters to categorize various
flow regimes. For homogeneous incompressible Navier–Stokes turbulence, the
most important of these is the Reynolds number $Re$, defined as the ratio of
the characteristic magnitudes of the non-linear inertial term and the viscous
term of the Navier–Stokes momentum equation (e.g., Pope, 2000). Herein we
define it by

$Re=\frac{U\lambda_{C}}{\nu},$ (1)

where $U=\sqrt{\langle\boldsymbol{v}\cdot\boldsymbol{v}\rangle}$ is the
characteristic root-mean-square (rms) speed of the fluctuations, $\lambda_{C}$
is the correlation scale (aka outer scale), and $\nu$ is the (kinematic)
viscosity; $\boldsymbol{v}(\boldsymbol{x},t)$ is the velocity field. Loosely,
$\lambda_{C}$ corresponds to the largest separation at which turbulent
fluctuations remain correlated, which in a hydrodynamic context can be thought
of as the size of the energy-containing eddies. (It is also often written as
$L$ and called the “characteristic length” scale.) Small $Re$ implies that the
viscous effects are significant and hence the nonlinear term is weak and will
not introduce significant nonlinearities into the system’s evolution.
Conversely, a large value of $Re$ implies that the nonlinear term plays a
significant role in the dynamics of the fluid.

This dynamic can be appreciated more clearly when $Re$ is expressed solely in
terms of length scales. One way this can be done is to introduce the
Kolmogorov dissipation scale (aka inner scale)
$\eta=(\nu^{3}/\epsilon)^{1/4}$, where
$\epsilon=\nu\langle\left(\nabla\times\boldsymbol{v}\right)^{2}\rangle$ is the
mean rate of kinetic energy dissipation (Kolmogorov, 1941; Tennekes & Lumley,
1972). A physical interpretation is that the Kolmogorov scale is where the
smallest eddies in the fluid become critically damped, due to their nonlinear
(aka turnover) timescale being equal to their dissipation timescale. Recall
also that the dissipation rate can be phenomenologically modelled as

$\displaystyle\epsilon_{\text{phenom}}=C_{\epsilon}\frac{U^{3}}{\lambda_{C}},$
(2)

where $C_{\epsilon}$ is treated as a fitting constant (e.g., Batchelor, 1970;
Tennekes & Lumley, 1972). Employing this in the definition of $\eta$ yields
the form

$Re\equiv
Re_{\eta}=\frac{1}{C_{\epsilon}^{1/3}}\left(\frac{\lambda_{C}}{\eta}\right)^{\frac{4}{3}},$
(3)

revealing that $Re$ is a measure of the bandwidth of the turbulent energy
cascade. A large $Re$ indicates there is a large separation between the outer
and inner scales. This larger bandwidth implies there are more scales where
the nonlinear term is strong enough to create turbulent structures and thereby
increase the intermittency of the flow (see, e.g., Matthaeus et al., 2015;
Parashar et al., 2019; Cuesta et al., 2022a). A small bandwidth, and hence a
small $Re$, implies that dissipation occurs very quickly and damps any
turbulent structures that the nonlinear term might try to create. Such
low-$Re$ situations are sometimes seen in planetary magnetosheaths (Czaykowska
et al., 2001; Hadid et al., 2015; Huang et al., 2017; Chhiber et al., 2018).

Estimating $Re$ for hydrodynamical systems, using Eq. (1), is straightforward
as all the required quantities are well defined and often readily determined
in experiments. For kinetic plasmas such as the solar wind, however, it is not
possible to write a Chapman-Enskog-like closure to define a viscosity (Chapman
& Cowling, 1990; Huang, 2008). (Some attempts have been made to estimate the
viscosity of kinetic systems; see e.g., Verma (1996); Zhuravleva et al.
(2019); Bandyopadhyay et al. (2023); Yang et al. (2023). This lack of a well-
defined viscosity also precludes using Eq. (3), as it means we cannot define
$\eta$. Typically, in kinetic plasmas, one must therefore resort to defining
an _effective_ Reynolds number. Some hydrodynamic studies have investigated
estimating the energy input into the system (Zhou et al., 2014), as well as
using more precise boundaries of the inertial range (Zhou, 2007; Zhou &
Thornber, 2016), in order to get around this lack of a clearly-defined inner
scale. Herein we describe two approaches to formulating an effective Reynolds
number.

The first approach is to apply Eq. (3) and use a different small scale — one
that is observationally calculable — as a signifier of the termination of the
inertial range. There are several reasonable options to choose from. For
example, in the solar wind, the _spectral break scale_ , $f_{b}$, the point at
which the power spectrum of the inertial range steepens, is thought to be a
good indicator of the onset of dissipation (Leamon et al., 1998; Yang et al.,
2022). Additionally, the _ion inertial length_ , $d_{\mathrm{i}}$, and also
the ion gyroradius, are frequently found to be in proximity to the break
scale, motivating their use as indicators of the onset of the kinetic range
(Chen et al., 2014; Franci et al., 2016; Wang et al., 2018; Woodham et al.,
2018; Parashar et al., 2019; Cuesta et al., 2022a; Lotz et al., 2023). We note
that $d_{\mathrm{i}}$ has the advantage that it only requires ion density to
calculate, rather than the high-resolution magnetic field data needed to
resolve spectral-steepening scales and calculate $f_{b}$. Its disadvantage is
that it does not capture the size of the turbulence amplitudes.

For example, consider the two different intervals in Fig. 1, each with very
similar $d_{\mathrm{i}}$ and outer scales $\lambda_{C}$ but with different
turbulence amplitudes. The use of $d_{\mathrm{i}}$ as an inner scale in Eq.
(3) consequently yields very similar $Re_{d_{\mathrm{i}}}$ for both cases
because it does not capture the different dynamics induced by the varying
turbulence strengths.

Figure 1: Power spectra (energy $E(k)$ vs. wavenumber $k$) for two spacecraft
data intervals with very similar outer and inner scales but different power
levels. Vertical lines indicate the respective correlation scales
$\lambda_{C}$ and ion inertial lengths $d_{\mathrm{i}}$. On the basis of the
power levels one would expect different turbulent behaviour from these
intervals. However, using $d_{\mathrm{i}}$ as the inner scale in Eq. (3)
implies they have almost the same $Re$. Using Eq. (5) we do capture this
difference, giving very different values of $Re_{\lambda_{T}}$.

Fortunately, there is a length scale that typically does depend on the energy
of the turbulent fluctuations due to their effect on the shape of the power
spectrum. This is the Taylor microscale $\lambda_{T}$ (Taylor, 1935;
Batchelor, 1970; Matthaeus et al., 2008), hereafter referred to simply as the
Taylor scale. See Fig. 1 for one such example. By employing it in a further
reformulation of $Re$ we can capture this strength-of-the-turbulence aspect.
The Taylor scale has multiple definitions and can be estimated in several
different ways which differ by order unity factors (denoted below by
$\gamma$). These can all be written, for the velocity field $\boldsymbol{v}$,
as

$\lambda_{T}^{2}=\gamma\frac{\langle\boldsymbol{v}^{2}\rangle}{\langle(\nabla\times\boldsymbol{v})^{2}\rangle},$
(4)

where the value of $\gamma$ depends on the specific definition of
$\lambda_{T}$ employed. For example, the traditional hydrodynamics usage is
that $\lambda_{T}$ is the curvature (at the origin) of the longitudinal
autocorrelation function so that $\gamma=5$ (e.g., Batchelor, 1970; Tennekes &
Lumley, 1972, p. 211). Herein we employ $\gamma=3$ because it corresponds to
the curvature of the traced correlation function, which is relatively simple
to calculate using spacecraft time series data; see Eqs. (7) and (8).

The inertial range comprises the scales $\ell$ which satisfy
$\eta\ll\ell\ll\lambda_{C}$.

Moreover, in hydrodynamics $\lambda_{T}$ lies between $\lambda_{C}$ and the
Kolmogorov scale (e.g., Pope, 2000). Eq. (4) makes it clear that $\lambda_{T}$
is related to the mean square spatial derivatives of the turbulent flow. It
can also be interpreted as the “single-wavenumber equivalent dissipation
scale” (Hinze, 1975). In plasma systems, the Taylor scale represents small-
scale turbulence physics that is not yet well understood, including its
relationship to other plasma parameters and the correlation length.

Re-expressed in terms of the Taylor scale, the exact hydrodynamic viscous
dissipation rate is
$\epsilon=\nu\gamma\langle\boldsymbol{v}^{2}\rangle/\lambda_{T}^{2}$. Equating
this to the $\epsilon_{\text{phenom}}$ relation, Eq. (2), yields another form
for the Reynolds number (Batchelor (1970, p. 118); Tennekes & Lumley (1972, p.
67)):

$Re\equiv
Re_{\lambda_{T}}=\frac{\gamma}{C_{\epsilon}}\left(\frac{\lambda_{C}}{\lambda_{T}}\right)^{2},$
(5)

The ratio of Taylor scale to the spectral break scale has been shown to have a
direct correlation with the decay rate (Matthaeus et al., 2008). Hence, one
would expect this definition of $Re$ to show variation with changing
turbulence amplitude and decay rates (as can be seen by the very different
values of $Re_{\lambda_{T}}$ in Fig. 1).

Note that $C_{\epsilon}$ is significantly less than unity. Hydrodynamic
simulations and experiments (Sreenivasan, 1998; Pearson et al., 2004) indicate
that

$\displaystyle C_{\epsilon}$ $\displaystyle\approx
0.5\frac{2}{9\sqrt{3}}\approx\frac{1}{15.6},$ (6)

where in the middle term the 0.5 value is empirical and the other values are
associated with “unit conversion” from a variant of Eq. (2) commonly used in
the hydrodynamic literature, namely
$\epsilon_{\text{phenom}}=Au_{1}^{3}/\ell_{f}$; here $U^{2}=3u_{1}^{2}$ and
$\ell_{f}=3\lambda_{C}/2$ is the correlation length for the longitudinal
velocity correlation function, all assuming isotropy (see, e.g., Batchelor
(1970); Tennekes & Lumley (1972); Pearson et al. (2004)). Thus, in
hydrodynamics, with $\gamma=3$, the prefactor in Eq. (5) is
$\gamma/C_{\epsilon}\approx 50$, and in Eq. (3) it is
$C_{\epsilon}^{-1/3}\approx 3$. The values in MHD, for solar wind-like
conditions, are $\gamma/C_{\epsilon}\approx 27$ and
$C_{\epsilon}^{-1/3}\approx 2$ (see Appendix A). These are the values we use
in the data analysis reported on below. However, one should keep in mind that
these values pertain to collisional MHD fluid models. The solar wind is an
almost collisionless plasma that can, in some circumstances, be well
approximated as an MHD fluid.

For a system like the solar wind, most velocity measurements have a time
cadence that is significantly longer than kinetic time scales (with the
exception of measurements from the MMS mission). Because of this, one cannot
reliably compute $\lambda_{T}$ for the velocity field. On the other hand,
magnetic field measurements have a significantly higher time cadence, allowing
one to explore kinetic scale physics. Hence most studies in the solar wind
compute the Taylor scale for the magnetic field. Given these constraints, we
also work (primarily) with magnetic field data in this study and compute
several types of _effective_ Reynolds numbers.

A history of estimating magnetic $Re$ in the solar wind is provided in the
introduction to Cartagena-Sanchez et al. (2022). Prior estimations have used
Eq. (5) and applied it to multi-spacecraft measurements, beginning with
Matthaeus et al. (2005) and continuing with Weygand et al. (2007, 2009, 2011)
and Zhou et al. (2020). Note that these studies use $\gamma/C_{\epsilon}=1$,
and thus essentially ignore this prefactor. The average values of
$\lambda_{C}$, $\lambda_{T}$, and $Re_{\lambda_{T}}$ from these studies are
summarised in Table 2, where we also indicate an appropriate value of
$\gamma/C_{\epsilon}$ to be used for comparison with the results we obtain
herein. All these studies used data from a combination of spacecraft at 1 au,
including ACE, Wind, and Cluster, and most investigated the relationship
between $Re$ and variables such as magnetic field orientation, wind speed, and
solar activity. Going beyond 1 au, this formulation has also been used to
estimate $Re$ at Mars (Cheng & Wang, 2022), and Voyager data has been used to
calculate it at very large distances from the Sun (Parashar et al., 2019).
Voyager data lacks sufficient resolution to calculate $\lambda_{T}$ and thus
$d_{\mathrm{i}}$ was used in the Eq. (3) formulation to estimate $Re$. Cuesta
et al. (2022a) supplemented this work with data from Parker Solar Probe and
Helios in a survey of variation in $Re$ throughout the heliosphere.

It is clear from the studies cited above that the Reynolds number plays a
pivotal role in understanding solar wind turbulence. Accurate estimation of
$Re$ can be used to validate theoretical predictions such as the enhanced
intermittency with increasing $Re$ (e.g., Van Atta & Antonia, 1980; Parashar
et al., 2015, 2019; Cuesta et al., 2022a), or its correlation with solar
activity (Zhou et al., 2020; Cheng & Wang, 2022). Different formulations need
to be compared to bolster these conclusions further. Additionally, a firmer
estimate of $Re$ will help refine the minimum scale separation required by an
experiment or simulation to faithfully capture the dynamics of such high $Re$
astrophysical systems; this is the so-called “minimum state” (Zhou, 2007,
2017). Therefore, to obtain reliable estimates of the solar wind’s (effective)
Reynolds number, a thorough comparison of computational techniques and their
implications is necessary.

This is the purpose of the present study. A large dataset of measurements from
the Wind spacecraft is compiled, allowing us to calculate $Re$ for nearly two
decades of data in three different ways: using either $f_{b}$ (obtained from
the magnetic energy spectrum) or $d_{\mathrm{i}}$ in Eq. (3), and using
$\lambda_{T}$, obtained from the autocorrelation function for
$\boldsymbol{b}$, in Eq. (5).

The structure of this paper is as follows. The dataset and its initial
cleaning are described in Sect. 2. Sect. 3 provides the methods for estimating
each of the scales; we calculate $Re_{\lambda_{T}}$ after first applying the
correction to $\lambda_{T}$ developed by Chuychai et al. (2014). In Sect. 4,
the three estimators are compared to each other and to the values obtained by
the aforementioned studies. Implications and limitations of these results are
discussed in Sect. 5.

##  2 Data

We use roughly 18 years (2004–2022) of data from NASA’s Wind spacecraft to
estimate $Re$ at 1 au. We process $\approx$ 12,000 12-hour intervals in the
solar wind. High-resolution (0.092 s) vector magnetic field data were obtained
from the Magnetic Field Investigation (MFI) (Lepping et al., 1995). Wind was
launched in 1994 and has operated at the Lagrangian point 1 (L1) since June
2004 in order to study plasma processes occurring in the near-Earth solar
wind. This mission has significantly contributed to understanding many aspects
of the solar wind, including electromagnetic turbulence (Wilson III et al.,
2021).

After downloading the data from NASA/GSFC’s Space Physics Data Facility
(SPDF), we split it into 12-hour intervals. This interval size is large enough
to contain a few correlation lengths but small enough to not average over
large-scale variations. Isaacs et al. (2015) demonstrated that (1 au)
intervals of 10-20 hours have “special significance” as they represent a range
where sufficient correlation times are sampled, making single-spacecraft
results coincide with those of multiple spacecraft.

Data gaps are linearly interpolated unless they comprise more than 10% of the
interval, in which case the interval is discarded. (This affected about 4% of
the intervals.) We initially processed 28 years of data, from 1995-01-01 to
2022-12-31, to compute various average quantities as well as turbulence
parameters such as the spectral slopes in the inertial and kinetic ranges, rms
amplitudes of the magnetic field and velocity, the Taylor scale, and the
correlation scale. The complete dataset comprises all available magnetic field
data, i.e., intervals containing shocks or from within the Earth’s
magnetosphere are _not_ removed. However, our analysis in the subsequent
sections of this paper focuses only on data from June 2004 and later, a period
when Wind was positioned at L1, away from the magnetosphere.

Given that it is also of interest to future analysis how quantities like the
Taylor scale relate to other properties of the turbulent plasma system — such
as electron density, cross-helicity, and solar activity — measurements of
electron and proton properties from Wind’s 3D Plasma (3DP) instrument were
also obtained, along with sunspot numbers from the World Data Center SILSO.

We note that the ion density from Wind has periods of anomalously small values
for a few months. To avoid issues associated with this we therefore always use
the electron density as a proxy for the proton density when calculating all
ion inertial lengths, ion plasma betas, and Alfvén speeds. Across the 28 years
of data, we obtained between 18,000 and 20,000 points for each variable,
depending on the amount of missing data. A full list of the variables in the
processed (and publicly available) data set can be found in Appendix B.

##  3 Method

Figure 2: Power spectral density (PSD) of a solar wind magnetic field
interval, raw (grey) and smoothed (black). Dashed power-law fits to estimates
of the inertial and kinetic ranges return spectral indices $q_{i}$ and
$q_{k}$. The left-hand vertical line indicates the intersection of these two
fits, denoted as the spectral break $f_{b}$; the ion inertial frequency
$f_{d_{\mathrm{i}}}$ is indicated by the right-hand vertical line.

We begin the analysis by determining several slopes for each of the magnetic
power frequency spectra obtained from the 12-hour intervals. Specifically, we
perform power-law fits in the inertial and kinetic ranges, denoting the power-
law exponents as $q_{i}$ and $q_{k}$, respectively. Nominal frequency
intervals for the inertial range (0.005–0.2 Hz) and kinetic range (0.5–1.4 Hz)
were chosen, consistent with those used by Wang et al. (2018). We then
identify the frequency at which (the extrapolation of) these powerlaws
intersect, calling this the spectral break frequency $f_{b}$. An example is
shown in Fig. 2. Any outliers, mostly in the form of anomalously large values
of $q_{k}$, are not included in the subsequent analysis, as described in Sect.
4. In the following, we will use the time scale associated with the break
frequency, i.e., $t_{b}=1/(2\pi f_{b})$, as a proxy for the inner (time)
scale.

Estimates for the Taylor scale $\lambda_{T}$ and the correlation scale
$\lambda_{C}$ are also needed and these are both computed using the
autocorrelation functions (see Fig. 26 in Bruno & Carbone (2013)). The
(normalized) temporal autocorrelation of the magnetic field fluctuations is
given by

$R(\tau)=\frac{\langle\boldsymbol{b}(t)\cdot\boldsymbol{b}(t+\tau)\rangle}{\langle\boldsymbol{b}^{2}\rangle},$
(7)

where $\boldsymbol{b}(t)=\boldsymbol{B}(t)-\langle\boldsymbol{B}(t)\rangle$ is
the magnetic field fluctuation at time $t$. The angle brackets denote a
suitable time ensemble average, implemented as a time average in this study.
Using Taylor’s frozen-in-flow hypothesis, we can convert time separations
$\tau$ into length separations $r$. (See Sect. 5 for a discussion of the
limitations of this hypothesis.)

Measurement of $\lambda_{C}$ requires a computation of the autocorrelation
function out to very large lags. On the other hand, measurement of
$\lambda_{T}$ requires iterative fitting at very small lags. It would quickly
become computationally expensive to use the high-time-cadence data to obtain
both quantities. Hence, for each 12-hour interval, the correlation length
$\lambda_{C}$ is computed from a down-sampled low-resolution (5 s) magnetic
field time-series out to roughly 10,000 s. We use the high-time-cadence (0.092
s) magnetic field data to compute autocorrelation functions only up to a lag
of 9.2 s; this is used to compute the magnetic Taylor scale $\lambda_{T}$.

Figure 3: A demonstration of the three methods used to calculate $\lambda_{C}$
using an interval comprising the second half of 2016-01-02. These include the
1/e (“e-folding”) method, giving $\lambda_{C}^{1/\mathrm{e}}$; the exponential
fit method, giving $\lambda_{C}^{\text{fit}}$; and the integral method, giving
$\lambda_{C}^{\text{int}}$.

The correlation scale $\lambda_{C}$ for $\boldsymbol{b}$ can be estimated from
$R(\tau)$ in three different ways, as shown in Fig. 3. We can perform an
exponential fit, we can find the separation at which the function falls to
1/e, or we can take the integral of the function
($\lambda_{C}=\int^{\infty}_{0}R(\tau)\,\mathrm{d}\tau$). The exponential fit
method is frequently used in the literature (Matthaeus et al., 2005; Zhou et
al., 2020; Bandyopadhyay et al., 2020; Phillips et al., 2022); multiple
exponential fits and a third-order polynomial have also been used (Weygand et
al., 2009, 2011; Cheng & Wang, 2022). In any case, this requires a decision
about how much of the autocorrelation to fit to. In this work, we fit a single
exponential to a range that extends to twice the value of the correlation
scale as obtained by the $1/\mathrm{e}$ method. We compute $\lambda_{C}$ from
the low-resolution autocorrelation using each of these three methods to
evaluate their consistency.

While it is straightforward to compute the Taylor scale in simulations, where
one has access to the full three-dimensional information, when working with
time series data from experiments we need to resort to an approximation. Since
$\lambda_{T}$ can be defined as the radius of curvature of the autocorrelation
function at the origin, we may use this definition to estimate it. (We do not
yet need to convert to spatial lags, so we work with the time-domain
equivalent, $\tau_{\text{TS}}$). This follows from the Taylor expansion of the
autocorrelation for $\tau\rightarrow 0$ (Batchelor, 1970; Tennekes & Lumley,
1972):

$R(\tau)=1-\frac{\tau^{2}}{2\tau_{\text{TS}}^{2}}+...$ (8)

Figure 4: An example of the process of refining the estimate of the Taylor
scale $\tau_{\text{TS}}$, using an interval comprising the second half of
2016-01-02. The three horizontal scales show the separations in units of lag,
time, and (Taylor frozen-flow equivalent) distance. Panel (a): Firstly, a
parabola is fit to the autocorrelation from the origin up to various values of
$\tau_{\text{fit}}$. In this example the fit is for lags less than
$\tau_{\text{fit}}=15$. The $x$-intercept of each parabola (which is off-scale
for this plot) produces an initial estimate $\tau_{\text{TS}}^{\text{est}}$.
Panel (b): Next, each of these estimates is plotted against
$\tau_{\text{fit}}$. A straight line is fit to these points and extrapolated
back to $\tau=0$, returning the Richardson extrapolation (R.E.) estimate,
$\tau_{\text{TS}}^{\text{ext}}$. Finally, the Chuychai correction (C.C.) is
applied using the kinetic range slope ($q_{k}=-3.03$ in this case) in Eq. (9).
This yields our final estimate, $\tau_{\text{TS}}$; we obtain a value of 9 s,
or approximately 4,000 km, for this particular interval.

In practice, this means fitting a parabola to $R(\tau)$ at the origin and
requires the high-resolution data provided by Wind so that we have enough
observations at small separations. It also requires an important decision: how
much of this high-resolution autocorrelation do we fit to? (Larger ranges
result in systematically larger estimates.) In order to reduce the
subjectivity of this decision, the Richardson extrapolation technique was
introduced in this context by Weygand et al. (2007): by fitting to a range of
values of maximum lag $\tau_{\text{fit}}$, then extrapolating back to 0 lags,
we obtain a refined estimate, $\tau_{\text{TS}}^{\text{ext}}$. In the
aforementioned work, the authors showed an apparent convergence of the final
estimate given by this technique as $\tau_{\text{fit}}$ increases. However,
Chuychai et al. (2014) showed with simulated data that, in fact, this
convergence depends upon the slope of the power spectrum at high frequencies.
In light of this, they produced a multiplicative correction factor, $r(|q|)$,
that is a function of this slope, given as

$r(|q|)=\begin{cases}-0.64(\frac{1}{|q_{k}|})+0.72,&\text{when }|q_{k}|<2\\\
-2.61(\frac{1}{|q_{k}|})+1.70,&\text{when }2\leq|q_{k}|<4.5\\\
-0.16(\frac{1}{|q_{k}|})+1.16,&\text{when }|q_{k}|\geq 4.5.\end{cases}$ (9)

We also apply this correction to our estimates, with the procedure we follow
depicted in Fig. 4. This gives us a final estimate $\tau_{\text{TS}}$. We fit
from a minimum lag of 1, equal to the time cadence (0.092 s), up to a maximum
lag $\tau_{\text{fit}}$ which was varied between 10 and 50 lags.

Finally, using the various (magnetic) scales, determined as outlined above, we
calculate estimates for effective Reynolds numbers in three distinct ways.
Specifically, we use $\lambda_{C}^{\text{fit}}$ (or its time scale analog) as
the outer scale and

  1. (i)

Eq. (3) with $d_{\mathrm{i}}$ as the inner scale and $C_{\epsilon}^{-1/3}=2$;

  2. (ii)

Eq. (3) with $t_{b}=1/(2\pi f_{b})$ as the inner (time)scale and
$C_{\epsilon}^{-1/3}=2$;

  3. (iii)

Eq. (5) with $\gamma/C_{\epsilon}=27$.

##  4 Results

Our analysis uses data from the period June 2004 to December 2022 when Wind
was always situated in the solar wind at L1. In about 6% of the intervals the
slope of the kinetic range, $q_{k}$, was unusually shallow (meaning
$|q_{k}|<1.7$) and therefore the final (corrected) estimate of the Taylor
scale came out to be negative. These outlier intervals were removed from the
following analysis but will be investigated in future work.

###  4.1 Correlation scale

Method | Mean (km) | Median (km) | SE (km)
---|---|---|---
Exponential fit | 899,000 | 769,000 | 5,000
$1/\mathrm{e}$ | 942,000 | 797,000 | 5,000
Integral | 880,000 | 808,000 | 4,000
Table 1: Statistical summary of estimates of the magnetic field correlation
length $\lambda_{C}$ by different methods. The standard error (SE) gives the
expected variation of the mean between samples of this size. Figure 5: Joint
(2D) histograms of the three $\lambda_{C}$ estimates with Pearson (linear)
correlation and Spearman (rank) correlation values and a dashed line of
equality. Marginal (1D) histograms of the $x$-variable are shown above each
plot. The $x$\- and $y$-axis limits have been set so as to include the bulk of
the data but exclude outliers. $\lambda_{C}^{\text{fit}}$: exponential fit
method, $\lambda_{C}^{1/\mathrm{e}}$: 1/e method, $\lambda_{C}^{\text{int}}$:
integral method.

Table 1 gives summary statistics of each of the three estimates of the
correlation length of the magnetic field, $\lambda_{C}$, and Fig. 5 shows
their marginal and joint distributions. Given the wide distribution of values,
all values are in line with those previously reported in the literature at 1
au, i.e., approximately $10^{6}$ km (see Table 2). Noting the logarithmic
scaling of the axes in this figure, we qualitatively find that the probability
distribution function of each estimator is log-normal. This is consistent with
the results of Ruiz et al. (2014) as well as the distribution of many other
solar wind quantities such as proton temperature, plasma beta, and Alfvén
speed (e.g., Hundhausen et al., 1970; Burlaga & Lazarus, 2000; Mullan & Smith,
2006; Veselovsky et al., 2010). In particular, the correlation scales are
positively skewed, with means larger than the corresponding medians. Looking
at the joint distributions, we see that the exponential fit and 1/e methods
agree very well with each other, with very high values of 0.99 for both the
Pearson and Spearman correlations, and most of the points lying close to the
equality line. (The Spearman correlation uses ranks to measure the
monotonicity of the relationship between two variables, rather than measure
their linear association.) This agreement is not surprising given the large-
scale statistical homogeneity of the solar wind. The autocorrelation functions
typically show approximately exponential fall-off (see, e.g., Fig. 3), with
deviations from $\lambda_{C}^{1/\mathrm{e}}\approx\lambda_{C}^{\text{fit}}$
only occurring for intervals that do not show steady turbulence (Ruiz et al.,
2014).

In contrast, the integral scale $\lambda_{C}^{\text{int}}$ shows a moderate
degree of scatter against either of the other two estimates, with correlations
of between 0.88 and 0.94. The greatest degree of scatter is present for values
of $\lambda_{C}^{\text{int}}$ less than about $10^{6}$ km. This disagreement
is likely due to occasional numerical issues with calculating the integral of
the autocorrelation. Ideally, the integral is computed out to infinity as $R$
asymptotically decays to 0. However, the finite size of the intervals and the
slight departures from “textbook-like” homogeneity and isotropy in some
intervals could introduce discrepancies between this and the exponential
estimates. Nonetheless, we conclude that, to a reasonable approximation, all
three methods give equivalent estimates for $\lambda_{C}$.

###  4.2 Taylor scale

Figure 6: Distributions of the uncorrected ($\lambda_{T}^{\text{ext}}$) and
corrected ($\lambda_{T}$) versions of the Taylor scale. Dashed vertical lines
indicate mean values for each distribution. The $x$-axis limits have been set
so as to include the bulk of the data but exclude outliers.  Figure 7: Joint
(2D) histograms of the three $Re$ estimates with Pearson correlation and
Spearman (rank) correlation values, log-space regression line fits, and a
dashed line of equality; and (above top axis) marginal (1D) histograms of the
$x$-variable in each plot. The $x$\- and $y$-axis limits have been set so as
to include the bulk of the data but exclude particularly extreme outliers
beyond these limits. Remaining outliers are shaded grey. Correlation
coefficients and marginal histograms are for all data values, whereas
regression lines are fitted to only the majority subsets of the data shown in
blue (see text).

Fig. 6 shows marginal distributions of both the uncorrected and corrected
versions of $\lambda_{T}$ for the magnetic field. Both have quasi-Gaussian
distributions, with a few large outliers. The distribution of $\lambda_{T}$
computed after applying the Chuychai correction factor is shifted to the left
because the (multiplicative) correction factor is almost always less than 1,
except for the 1% of intervals with particularly steep kinetic range slopes
($q_{k}<-3.7$). The mean $q_{k}$ is $-2.64$, resulting in an average
correction factor of $-2.61/2.64+1.7=0.71$, following Eq. (9). We therefore
end up with a mean of $\lambda_{T}$ that is about two thirds that of
$\lambda_{T}^{\text{ext}}$. We find that this final mean of 3,225 km is in
good agreement with the literature (see Table 2). Prior estimates of
$\lambda_{T}$ in the solar wind at 1 au vary between $\sim 1,000$ km and $\sim
7,000$ km, values which lie within the distributions of $\lambda_{T}$ (both
extrapolated or Chuychai-corrected) shown in Fig. 6.

Authors (Year) | Spacecraft |  $\lambda_{C}$ ($10^{6}$ km) |  $\lambda_{T}$ (km) | $Re_{\lambda_{T}}$ |
---|---|---|---|---|---
Matthaeus et al. (2005) | ACE-Wind-Cluster | 1.2 | 2,478 $\pm$ 702 | 230,000 | ($\times$27)
Weygand et al. (2007) | Cluster | 1.2 (from above) | 2,400 $\pm$ 100 | 260,000 | ($\times$27)
Weygand et al. (2009) | ACE-Wind-Cluster + 6 others | 2.92 | 1,000 $\pm$ 200 | 12,600,000 | ($\times$27)
Weygand et al. (2011) | ACE-Wind-Cluster + 8 others | 1–2.8 | 1,200–3,500 | 4,000,000 | ($\times$27)
Zhou et al. (2020) | ACE-Wind-Cluster | 1.14 | 2,459 | 300,000* | ($\times$27)
Bandyopadhyay et al. (2020) | MMS | 0.32 | 6,933 | 2,000** | ($\times$27)
This work | Wind | 0.899 | 3,220 | 3,406,000 |
Table 2: Average estimates of $\lambda_{C}$, $\lambda_{T}$ and an effective
Reynolds number in the solar wind at 1 au, with the latter calculated using
Eq. (5). Shown are the values determined in this work (given in bold) and in
some previous studies. Note that for direct comparison with this work, the
$Re_{\lambda_{T}}$ values from these earlier studies should be multiplied by
the previously neglected prefactor of $\gamma/C_{\epsilon}=27$, as indicated
by the ‘$(\times 27)$’ in the final column. This factor is already included in
our estimates. When calculating $\lambda_{T}$ all studies listed employed
$\gamma=3$, sometimes without explicitly stating so. All studies used at least
one exponential fit to compute $\lambda_{C}$. All except Matthaeus et al.
(2005) and Bandyopadhyay et al. (2020) used Richardson extrapolation to
compute $\lambda_{T}$; none, other than this work, used the Chuychai
correction. Values are expressed as ranges when the study grouped scales by
other variables such as magnetic field orientation.
* This mean value was reported in a follow-up article (Zhou & He, 2021).
** $Re$ was not calculated explicitly in this article.

###  4.3 Reynolds number

Having obtained estimates for the correlation length and Taylor scale of the
magnetic field fluctuations (and also for $d_{\mathrm{i}}$ and $t_{b}$) we may
use the procedures detailed at the end of Sect. 3 to calculate three distinct
effective Reynolds numbers. Fig. 7 shows the marginal and joint distributions
of these different estimates, as well as regression line fits. After applying
a logarithmic transformation each distribution appears approximately Gaussian,
suggesting a log-normal distribution. Comparing these marginal distributions
and the summary statistics given in Table 3, we can see that the three
estimates span multiple orders of magnitude, with, very roughly,
$Re_{\lambda_{T}}\sim 10Re_{d_{\mathrm{i}}}\sim 30Re_{t_{b}}$, for the mean
values.

The joint distributions show considerably more scatter than those of the
$\lambda_{C}$ estimates. The strongest linear association between any two
estimates is that between $Re_{d_{\mathrm{i}}}$ and $Re_{t_{b}}$ (Pearson
correlation = 0.72). This is shown by the majority of points lying in a
relatively thin linear band close to the equality line. We can also see that
the $Re_{t_{b}}$ estimates tend to be smaller than $Re_{d_{\mathrm{i}}}$. This
is an indication that the break scale is typically larger than
$d_{\mathrm{i}}$ by a factor of 2-3 in the solar wind (Leamon et al., 1998). A
dependence of the break scale on plasma $\beta$ is also well known (Chen et
al., 2014; Franci et al., 2016). The statistical details of any such potential
correlations will be explored in a follow-up study.

$Re_{\lambda_{T}}$ shows a much weaker linear association with the other two
methods of only 0.38 (with $Re_{d_{\mathrm{i}}}$) and 0.43 (with
$Re_{t_{b}}$). In addition to having the lowest Pearson correlation,
$Re_{\lambda_{T}}$ and $Re_{t_{b}}$ also have the lowest Spearman correlation,
showing that even after accounting for outliers, which have less influence on
this latter metric, it still remains a rather weakly positive association. On
the other hand, outliers do have a clear influence on the linear association
of $Re_{\lambda_{T}}$ vs. $Re_{d_{\mathrm{i}}}$, shown by the substantial
increase in the Spearman correlation (0.74) over the Pearson correlation
(0.38).

Despite these only moderately strong associations between the estimates, it is
important however to note the density of the points. All these distributions
show significant scatter of a small population in which the estimates differ
by up to an order of magnitude. Notably, the joint distributions of
$Re_{t_{b}}$ have a roughly triangular sub-population of points that shows
little to no relationship with the other estimates. This is seen in the upper
left of the plot of $Re_{\lambda_{T}}$ vs. $Re_{t_{b}}$, and the lower right
of $Re_{t_{b}}$ vs. $Re_{d_{\mathrm{i}}}$. This population (identified as
$Re_{\lambda_{T}}/Re_{t_{b}}>50$) represents about 27% of all observations,
and is shown as the grey points in Fig. 7. After removing this population, all
correlations increase to at least 0.68. The potential reasons for
significantly larger $t_{b}$ and hence a smaller $Re_{t_{b}}$ could include
errors in automated fitting and extreme intervals with atypical power spectra.
As with the other outliers, a detailed investigation of these is deferred to a
follow-up study. In cases where the power spectrum is well-behaved, with
typical inertial and kinetic range slopes ($q_{i}$ and $q_{k}$) and a well-
defined breakpoint, it might be safe to estimate $Re_{tb}$ and multiply it by
30 to estimate $Re_{\lambda_{T}}$.

$Re$ | Mean | Median | SE
---|---|---|---
$Re_{t_{b}}$ | 116,000 | 64,000 | 2,000
$Re_{d_{\mathrm{i}}}$ | 330,000 | 226,000 | 3,000
$Re_{\lambda_{T}}$ | 3,406,000 | 1,686,000 | 68,000
Table 3: Statistical summary of the estimates of effective $Re$ obtained by
the different methods. (SE = standard error of the mean.)

As well as the agreement between methods, it is also of interest how our
estimates of $Re$ match up with those previously reported. In particular,
given the prevalence of the Taylor scale method, we compare values of
$Re_{\lambda_{T}}$ in Table 2. The values for $Re$ in this table vary by a
factor of $\approx 6000$, from 54,000 to 340,000,000 after multiplying by the
prefactor. The mean value of $Re_{\lambda_{T}}=3,406,000$ from the present
study is on the same order of magnitude as the results from three of the
previous works. The much larger values given in Weygand et al. (2009) and
Weygand et al. (2011) were mainly attributed to the smaller values obtained
for $\lambda_{T}$. Conversely, Bandyopadhyay et al. (2020) noted that their
value of $\lambda_{T}$ calculated from a single 5-hour interval of MMS data
was about three times larger than previous estimates, while their estimation
of $\lambda_{C}$ was smaller than other estimates. Hence they computed a much
smaller value of $Re$. Three reasons were suggested for this: 1) interval
length, separation and mixing effects, 2) intrinsic variability in the solar
wind, and 3) differences between the geometric formation of the Cluster (to
which they were comparing their results) and MMS spacecraft. Our work herein
emphasizes that point 2) is indeed pertinent. In particular, our results show
the considerable intrinsic variability of the solar wind properties
(particularly $\lambda_{C}$ and $\lambda_{T}$), giving rise to variability in
the values of effective $Re$. On the plus side, this sampling variability
suggests that the results of all the cited studies may in fact be consistent
with each other, as they lie within the distribution of values found in our
study.

##  5 Conclusion

We present a thorough investigation and review of calculating estimates of
(effective) Reynolds numbers for the solar wind at 1 au, using 18 years of
data from NASA’s Wind spacecraft. As this dataset lacks high-time-cadence
velocity measurements, we employ magnetic field data to estimate $\lambda_{C}$
and $\lambda_{T}$ for the magnetic field. These are assumed to be comparable
to their velocity field equivalents, in line with previously published
results. More precisely, in using the magnetic length scales in (5) we are
assuming that
$\lambda_{C}^{b}/\lambda_{T}^{b}\approx\lambda_{C}^{v}/\lambda_{T}^{v}$.

We first compare three different ways of calculating the correlation scale and
find good agreement between all methods, albeit with a greater scatter for the
integral method. The mean values obtained for $\lambda_{C}$, between 880,000
km and 942,000 km, are consistent with previously reported values of about
$10^{6}$ km.

We then apply the correction factor developed by Chuychai et al. (2014) to our
estimates of the Taylor scale in order to reduce any remaining bias after
using the Richardson extrapolation technique. This correction factor typically
reduces the estimate of the Taylor scale, significantly shifting the
distribution to smaller values. In particular, the mean reduces to 3,225 km,
roughly 2/3 of the uncorrected mean value of 4,772 km. Both values are
consistent with previous estimates, given the wide spread of the distribution.

Finally, we compute effective Reynolds numbers using three distinct methods.
It should be noted that we include proportionality factors in our
calculations. In particular, we highlight that the factor in
$Re_{\lambda_{T}}$ of $C_{\epsilon}\approx 27$ was not included in many
previously published estimates (see Table 2 and Eqn. (5)).

While very strong correlations were observed for the three different methods
of estimating $\lambda_{C}$, the correlations between the associated estimates
of the effective Reynolds number were only moderate to strong, with a
considerable amount of scatter. The mean values determined by these methods
ranged from 116,000 for $Re_{t_{b}}$ to 3,406,000 for $Re_{\lambda_{T}}$.
Putting these into perspective, previously reported values of $Re$ at 1 au
exhibit substantial variability, ranging from approximately $10^{6}$ to
$10^{8}$. Most of our estimates of $Re$ comfortably fit within this
distribution, though an outlier population of small values of $Re_{t_{b}}$
warrants future investigation.

Ultimately, we conclude that more theoretical work is needed to better
understand which definition of an effective Reynolds number of the solar wind
is most appropriate. The key task is to identify scales that have a physical
meaning. For the $t_{b}$ or $d_{\mathrm{i}}$ approximations of the inner
scale, the implication is that ion-scale physics plays the most significant
role in energy dissipation and terminates the inertial range. This, however,
discounts the role of a sub-ion-scale cascade and its implications for
electron physics (Matthaeus et al., 2008; Alexandrova et al., 2009; Sahraoui
et al., 2009; Schekochihin et al., 2009; Boldyrev et al., 2013). Moreover,
these estimates are insensitive to the variability of the power input at large
scales and hence the cascade rate. On the other hand, the $\lambda_{T}$-based
estimate of $Re$ indirectly folds in the cascade rate through its dependence
(empirically in the solar wind but directly in hydrodynamics) (Pope, 2000;
Matthaeus et al., 2008). This makes $Re_{\lambda_{T}}$ a more physically
motivated estimate amongst the three considered.

Moreover, having obtained statistical relationships between different
estimates, these can be leveraged in situations where only one estimate is
calculable. The decision on which estimator to use rests on the assumptions
one elects to make and on the resolution of the available data. These
considerations are summarized below.

  * •

$Re_{t_{b}}$ requires calculation of the spectral break scale. This process is
subject to varying methods and numerical challenges, including spectra that do
not always show clear breaks. In our work, we calculated $t_{b}$ as the
intersection of (the extrapolation of) two power-law fits to magnetic field
spectra, which requires decisions on what intervals to choose for the inertial
and dissipation ranges.

  * •

Alternatively, one can simply use the ion inertial length $d_{\mathrm{i}}$ to
approximate the break scale and calculate $Re_{d_{\mathrm{i}}}$ (Parashar et
al., 2019; Cuesta et al., 2022a). This requires only the ion density (and
correlation length). However, it appears that changing solar wind conditions
affect which scale is best associated with the spectral break. Specifically,
$d_{\mathrm{i}}$ is the best approximation at low plasma $\beta$ values, the
ion gyroradius $\rho_{i}$ is best at high $\beta$ values (Chen et al., 2014),
and for typical solar wind values of $\beta\approx 1$, the ion cyclotron
resonance scale is the best (Woodham et al., 2018). Under conditions where one
might not have high-time-cadence measurements of the desired variables, it is
likely that one could still easily obtain reasonable estimates for both
$d_{\mathrm{i}}$ and the outer scale (e.g., $\lambda_{C}$) and employ these to
estimate an effective Reynolds number.

  * •

Using the Taylor scale-based Reynolds number, $Re_{\lambda_{T}}$, is a more
robust formulation for estimating $Re$ than the two listed above because of
its empirical dependence on the cascade rate (Matthaeus et al., 2008). This
benefit is shown by the prevalence of this formulation in the literature
(Matthaeus et al., 2005; Weygand et al., 2007, 2009, 2011; Zhou et al., 2020;
Phillips et al., 2022). We show here that the use of a correction factor
(Chuychai et al., 2014) makes a significant difference in the estimates of
$\lambda_{T}$ and hence $Re_{\lambda_{T}}$. However, as discussed above, the
$\lambda_{T}$ definition of $Re$ has a proportionality factor that is, in
part, determined by a phenomenological or empirical fitting for the (kinetic
energy) dissipation rate. Calculating $\lambda_{T}$ also requires high-
resolution data, which is not always available; for example, outer heliosphere
Voyager observations are so restricted (Parashar et al., 2019). Furthermore,
while this does not affect the validity of this formulation, we note that weak
cascade rates have been shown to result in $\lambda_{T}$ being smaller than
the break scale, inverting the hydrodynamic ordering (Matthaeus et al., 2008).
This is believed to be due to greater relevance of electron dissipation
(relative to proton) in these weak cascade circumstances (Matthaeus et al.,
2016).

A limitation of this work is that it relies on the Taylor hypothesis to
convert from single-spacecraft time separations to length separations. This
assumes that the bulk flow is sufficiently fast that local variations in time
can be effectively ignored (see Verma, 2022, for solar wind context). The
Taylor hypothesis relates to the well-studied “sweeping” hypothesis, whereby
large-scale fluctuations sweep (i.e., advect) smaller-scale fluctuations
(Kraichnan, 1965; Tennekes, 1975; Zhou, 2021). Although invoking Taylor’s
hypothesis at kinetic scales might introduce substantial inaccuracies, it has
nonetheless been shown, numerically and from observations, to be a reasonably
good approximation up to second-order statistics (Perri et al., 2017; Chhiber
et al., 2018; Roberts et al., 2022). This is also true under a model that
incorporates sweeping phenomenology (Bourouaine & Perez, 2019; Perez et al.,
2021). Furthermore, for the present analysis, we note that this assumption
does not affect the results for $Re_{\lambda_{T}}$, because both of the scales
involved are in fact left as time scales for this calculation. Another aspect
that we did not address in this study is the issue of anisotropy in
$\lambda_{C}$ and $\lambda_{T}$ (e.g., Weygand et al., 2009, 2011; Cuesta et
al., 2022b; Roy et al., 2022). We reiterate that no data filtering was
conducted, except to remove intervals with significant missing data, limit the
analysis time period to June 2004 onward, and remove outliers where
$|q_{k}|<1.7$.

Finally, we envisage that the full 28-year dataset and the accompanying code
that we have provided as a data product will be useful to the scientific
community for future large-scale statistical analysis and data mining, for
Wind and other missions. Future work will start investigating correlations,
dimensionality reduction, and machine learning models.

##  6 Acknowledgments

We are grateful to the Wind MFI and 3DP instrument teams for the data and NASA
GSFC’s Space Physics Data Facility for providing access to it
111https://spdf.gsfc.nasa.gov/. We acknowledge the World Data Center SILSO at
the Royal Observatory of Belgium in Brussels for providing the sunspot data
222http://www.sidc.be/silso/. We also thank Bill Matthaeus for providing
feedback on the manuscript.

##  7 Author contributions

  * •

TNP: Conceptualized and supervised the project and edited the draft
manuscript.

  * •

DW: Refined and extended the pipeline, created the final dataset, and wrote
the draft manuscript.

  * •

SO: Identified and calculated the prefactors in the Reynolds number equations.

  * •

KDL: Conducted preliminary analysis with guidance from TNP and MF.

  * •

All authors discussed results and implications and helped edit the final
manuscript.

##  Appendix A Determination of the $C_{\epsilon}$ prefactor for MHD

A standard phenomenological estimate for the kinetic energy dissipation rate
($\epsilon^{v}$) in Navier–Stokes turbulence is

$\epsilon_{\text{phenom}}=A\frac{u_{1}^{3}}{\ell_{f}^{v}}=C_{\epsilon}\frac{U^{3}}{\lambda_{C}}$
(A1)

(e.g., Batchelor, 1970; Tennekes & Lumley, 1972; Pope, 2000). Here,
$\langle\boldsymbol{v}\cdot\boldsymbol{v}\rangle=U^{2}$ and $u_{1}$ is the rms
velocity for one component of $\boldsymbol{v}$. Also $\ell_{f}^{v}$ is the
correlation length associated with the longitudinal correlation function
(Batchelor, 1970), whereas $\lambda_{C}$ is that for the traced correlation
function
$R^{v}(r)=\langle\boldsymbol{v}(\boldsymbol{x})\cdot\boldsymbol{v}(\boldsymbol{x}+\boldsymbol{r})\rangle$,
equivalent to our $\lambda_{C}$ definition in the main text. The dimensionless
coefficients $A$ and $C_{\epsilon}$ are treated as constants that may be
determined using experiments and/or simulations (Sreenivasan, 1998; Pearson et
al., 2004). For isotropic turbulence, the relations $U^{2}=3u_{1}^{2}$,
$\lambda_{C}=2\ell_{f}^{v}/3$, and $C_{\epsilon}=2A/(9\sqrt{3})$ hold. As the
middle ‘component-based’ version is founded on the assumption of isotropy, in
this work we instead employ the rightmost ‘trace-based’ variant which does not
assume isotropy; this is given as Eq. (2) above.

In the literature a variety of notations are in use for what we have called
$A$ and $C_{\epsilon}$, and indeed some works use $C_{\epsilon}$ for the $A$
in Eq. (A1) (e.g., Pearson et al., 2004); clearly this should not be confused
with the $C_{\epsilon}$ we employ herein. For clarity, and in line with the
notation of Batchelor (1970, eq. 6.1.1), we always use $A$ to denote a
component-based fitting value.

We wish to determine a value for $C_{\epsilon}$ that is applicable in MHD.
This requires taking into account the dissipation of magnetic as well as
kinetic energy. With superscripts $v$ and $b$ denoting velocity and magnetic
quantities, respectively, we may write the total energy decay rate as
$\epsilon^{\text{MHD}}=\epsilon^{v}+\epsilon^{b}$.

Using an Elsasser variable
($\boldsymbol{z}_{\pm}=\boldsymbol{v}\pm\boldsymbol{b}$) von Kármán–Howarth
equation analysis for incompressible MHD, Linkmann et al. (2015, 2017)
developed a theory for $A^{\text{MHD}}$ (denoted $C_{\epsilon,\infty}$
therein). For simplicity here we restrict attention to situations with low
cross helicity, i.e., $\langle\boldsymbol{v}\cdot\boldsymbol{b}\rangle\approx
0$. Consequently,
$\langle|\boldsymbol{z}_{+}|^{2}\rangle\approx\langle|\boldsymbol{z}_{-}|^{2}\rangle=Z^{2}=3W^{2}$,
and the two longitudinal Elsasser correlation lengths, $\ell_{f}^{\pm}$, are
approximately equal. Assuming further that the longitudinal correlation
lengths for $\boldsymbol{v}$ and $\boldsymbol{b}$ are also approximately
equal, the (low cross helicity) Linkmann et al. result is equivalent to

$\epsilon^{\text{MHD}}_{\text{phenom}}=A\frac{W^{3}}{\ell_{f}^{+}}=\frac{2A}{9\sqrt{3}}\frac{Z^{3}}{\lambda_{C}^{v}},$
(A2)

where $A\equiv A^{\text{MHD}}$ and we have made use of the isotropic relation
$\ell_{f}^{v}=\frac{3}{2}\lambda_{C}^{v}$. Eq. (A2) applies to the total
(viscous plus resistive) dissipation and is the MHD analogue of Eq. (2) above.
Although, formally, it only applies for low cross helicity cases, it is likely
to be approximately valid under somewhat more general circumstances (Linkmann
et al., 2017; Bandyopadhyay et al., 2018).

Next, we make use of the Alfvén ratio,
$r_{A}=U^{2}/\langle\boldsymbol{b}^{2}\rangle$, to write
$\epsilon^{\text{MHD}}$ in terms of $\epsilon^{v}$. In terms of $U$ and
$r_{A}$, the zero cross helicity form for $Z^{2}$ is

$Z^{2}=U^{2}+\langle\boldsymbol{b}^{2}\rangle=U^{2}\left(1+\frac{1}{r_{A}}\right)=U^{2}F(r_{A}),$
(A3)

which defines $F(r_{A})$.

Recall also that
$\epsilon^{\text{MHD}}=\nu\langle\boldsymbol{\omega}^{2}\rangle+\mu\langle\boldsymbol{j}^{2}\rangle$,
with $\nu$ the kinematic viscosity, $\mu$ the magnetic resistivity,
$\boldsymbol{\omega}=\nabla\times\boldsymbol{v}$ the vorticity, and
$\boldsymbol{j}=\nabla\times\boldsymbol{b}$ the electric current density.
Assuming a Prandtl number of order unity ($\nu\approx\mu$) and that
$\langle\boldsymbol{j}^{2}\rangle/\langle\boldsymbol{\omega}^{2}\rangle\approx
1/r_{A}$, as is commonly seen in MHD simulations, $\epsilon^{\text{MHD}}$ can
be re-expressed without explicit reference to the dissipation coefficients:

$\epsilon^{\text{MHD}}\approx\epsilon^{v}F(r_{A}).$ (A4)

Finally, using Eqs. (A2) and (A4), we can write a (small cross helicity)
approximation for the kinetic energy dissipation rate in MHD:

$\epsilon^{v}\approx\left(\frac{2A}{9\sqrt{3}}\sqrt{F}\right)\frac{U^{3}}{\lambda_{C}}.$
(A5)

The bracketed factor might be called $C_{\epsilon}^{\text{MHD},v}$ and can be
identified with $C_{\epsilon}$ in our Eq. (2). Observationally, for the solar
wind, $r_{A}\approx 1/2$ yielding $F\approx 3$ (e.g., Perri & Balogh, 2010).
Results from MHD simulations (Linkmann et al., 2017; Bandyopadhyay et al.,
2018)333In both these works $A$ is denoted as $C_{\epsilon,\infty}$ and here
we employ double their numerical value for $A$ because of a definitional
difference between their $L_{\pm}$ and our $\ell_{f}^{\pm}$. indicate that
$A\approx 0.5$ for situations with zero or moderate mean magnetic field and
low or moderate cross helicity, as is relevant to the solar wind. Using these
values we obtain $C_{\epsilon}^{\text{MHD},v}\approx 0.11$, which is about
twice the hydrodynamics estimate of $C_{\epsilon}^{\text{hydro}}\approx 0.064$
(via $A\approx 0.5)$; see Sreenivasan (1998) and Pearson et al. (2004). This
gives us values for the prefactor of Eq. (5) of
$\frac{\gamma}{C_{\epsilon}}=\frac{3}{0.11}\approx 27$, and of Eq. (3) of
$C_{\epsilon}^{1/3}\approx 2$.

The results obtained in this appendix are most relevant for systems governed
by the incompressible collisional MHD equations. Thus, application of these
results to the nearly collisionless solar wind needs to be undertaken with
caution.

##  Appendix B Data Product

Averages of each variable in our dataset are given in Table 4. The dataset (in
CSV form), along with metadata describing the variables and the code used to
extract and process the data, are available on
GitHub444reynolds_scales_project codebase: https://github.com/daniel-
wrench/reynolds_scales_project. under a 2-Clause BSD License and are archived
in Zenodo (Wrench, 2023). The code has been designed so as to make it
relatively simple to apply to data from other missions available in CDAWeb.
That is, it should be straightforward to adapt for projects interested in
calculating these variables for different heliophysics and space weather
environments.

Symbol | Name | Mean value | Unit
---|---|---|---
$SN$ | Sunspot number | 56.3 | -
$M_{A}$ | Alfvén Mach number | 7.36 | -
$M_{s}$ | Sonic Mach number | 15.31 | -
$\beta_{e}$ | Electron plasma beta | 0.82 | -
$\beta_{p}$ | Proton plasma beta | 0.53 | -
$\sigma_{c}$ | Cross helicity | 0.01 | -
$\sigma_{R}$ | Residual energy | -0.44 | -
$R_{A}$ | Alfvén ratio | 0.46 | -
$\cos(A)$ | Alignment cosine | 0.01 | -
$q_{i}$ | Inertial range slope | -1.68 | -
$q_{k}$ | Kinetic range slope | -2.64 | -
$Re_{\lambda_{T}}$ | Reynolds number ($\lambda_{T}$) | 3,406,000 | -
$Re_{d_{\mathrm{i}}}$ | Reynolds number ($d_{\mathrm{i}}$) | 330,000 | -
$Re_{t_{b}}$ | Reynolds number ($t_{b}$) | 116,000 | -
$f_{b}$ | Spectral break frequency | 0.25 | Hz
$t_{b}$ | Spectral break time scale | 14.3 | s
$B_{\text{0}}$ | Magnetic field magnitude (rms) | 5.49 | nT
$\delta b$ | Magnetic field fluctuations (rms) | 3.83 | nT
$\delta b/B_{0}$ | Normalized magnetic field fluctuations | 0.71 | nT
$n_{e}$ | Electron density | 4.18 | cm-3
$n_{\alpha}$ | Alpha density | 0.14 | cm-3
$T_{e}$ | Electron temperature | 12.9 | eV
$T_{p}$ | Proton temperature | 11.0 | eV
$T_{\alpha}$ | Alpha temperature | 63.8 | eV
$\rho_{e}$ | Electron gyroradius | 1.78 | km
$\rho_{p}$ | Proton gyroradius | 63.9 | km
$d_{e}$ | Electron inertial length | 3.12 | km
$d_{\mathrm{i}}$ | Proton inertial length | 134 | km
$l_{d}$ | Debye length | 0.02 | km
$\lambda_{C}^{\text{fit}}$ | Correlation length scale (exp. fit) | 899,000 | km
$\lambda_{C}^{\text{exp}}$ | Correlation length scale (1/e) | 942,000 | km
$\lambda_{C}^{\text{int}}$ | Correlation length scale (integral) | 880,000 | km
$\lambda_{T}^{\text{ext}}$ | Taylor length scale (raw) | 4,770 | km
$\lambda_{T}$ | Taylor length scale (corrected) | 3,220 | km
$V_{0}$ | Velocity magnitude (rms) | 439 | km/s
$V_{r}$ | Radial velocity | 438 | km/s
$\delta v$ | Velocity fluctuations (rms) | 26.2 | km/s
$v_{A}$ | Alfvén speed | 65.5 | km/s
$v_{T_{e}}$ | Electron thermal velocity | 1490 | km/s
$v_{T_{p}}$ | Proton thermal velocity | 30.5 | km/s
$\delta b_{A}$ | Magnetic field fluctuations (Alfven units, rms) | 42.4 | km/s
$z^{+}$ | Positive Elsasser variable (rms) | 48.9 | km/s
$z^{-}$ | Negative Elsasser variable (rms) | 48.4 | km/s
Table 4: List of the key variables in our Wind data product, comprised of
statistics for every 12-hours from 1995-2022. The mean values are for the
cleaned 18-year dataset at L1 used in this study. While not shown here, we
also provide a few additional variables such the time-scale versions of the
length scales, the uncertainty of the Taylor scale, and the amount of missing
data for each raw interval. The complete metadata, including the equations
used to derive secondary variables such as gyroradii and cross-helicity, can
be found in the GitHub README.

## References

  * Alexandrova et al. (2009) Alexandrova, O., Saur, J., Lacombe, C., et al. 2009, Phys. Rev. Lett., 103, 165003
  * Balbus & Hawley (1998) Balbus, S. A., & Hawley, J. F. 1998, Reviews of modern physics, 70, 1
  * Bandyopadhyay et al. (2018) Bandyopadhyay, R., Oughton, S., Wan, M., et al. 2018, Phys. Rev. X, 8, 041052, doi: 10.1103/PhysRevX.8.041052
  * Bandyopadhyay et al. (2020) Bandyopadhyay, R., Matthaeus, W. H., Chasapis, A., et al. 2020, The Astrophysical Journal, 899, 63, doi: 10.3847/1538-4357/ab9ebe
  * Bandyopadhyay et al. (2023) Bandyopadhyay, R., Yang, Y., Matthaeus, W. H., et al. 2023, Phys. Plasmas, 30, 080702, doi: 10.1063/5.0146986
  * Batchelor (1970) Batchelor, G. K. 1970, The Theory of Homogeneous Turbulence (Cambridge, England: Cambridge University Press)
  * Boldyrev et al. (2013) Boldyrev, S., Horaites, K., Xia, Q., & Perez, J. C. 2013, The Astrophysical Journal, 777, 41
  * Bourouaine & Perez (2019) Bourouaine, S., & Perez, J. C. 2019, The Astrophysical Journal Letters, 879, L16, doi: 10.3847/2041-8213/ab288a
  * Bruno & Carbone (2013) Bruno, R., & Carbone, V. 2013, Living Rev. Solar Phys, 10, doi: 10.12942/lrsp-2013-2
  * Burlaga & Lazarus (2000) Burlaga, L. F., & Lazarus, A. J. 2000, Journal of Geophysical Research: Space Physics, 105, 2357, doi: https://doi.org/10.1029/1999JA900442
  * Cartagena-Sanchez et al. (2022) Cartagena-Sanchez, C. A., Carlson, J. M., & Schaffner, D. A. 2022, Physics of Plasmas, 29, 032305, doi: 10.1063/5.0073207
  * Chapman & Cowling (1990) Chapman, S., & Cowling, T. G. 1990, The mathematical theory of non-uniform gases: an account of the kinetic theory of viscosity, thermal conduction and diffusion in gases (Cambridge University Press)
  * Chen et al. (2014) Chen, C. H. K., Leung, L., Boldyrev, S., Maruca, B. A., & Bale, S. D. 2014, Geophysical Research Letters, 41, 8081, doi: 10.1002/2014GL062009
  * Cheng & Wang (2022) Cheng, L., & Wang, Y. 2022, The Astrophysical Journal, 941, 37, doi: 10.3847/1538-4357/aca0f2
  * Chhiber et al. (2018) Chhiber, R., Chasapis, A., Bandyopadhyay, R., et al. 2018, Journal of Geophysical Research: Space Physics, 123, 2018JA025768, doi: 10.1029/2018JA025768
  * Chuychai et al. (2014) Chuychai, P., Weygand, J. M., Matthaeus, W. H., et al. 2014, Journal of Geophysical Research: Space Physics, 119, 4256, doi: 10.1002/2013JA019641
  * Cuesta et al. (2022a) Cuesta, M. E., Parashar, T. N., Chhiber, R., & Matthaeus, W. H. 2022a, The Astrophysical Journal Supplement Series, 259, 23, doi: 10.3847/1538-4365/ac45fa
  * Cuesta et al. (2022b) Cuesta, M. E., Chhiber, R., Roy, S., et al. 2022b, The Astrophysical Journal Letters, 932, L11, doi: 10.3847/2041-8213/ac73fd
  * Czaykowska et al. (2001) Czaykowska, A., Bauer, T. M., Treumann, R. A., & Baumjohann, W. 2001, Annales Geophysicae, 19, 275, doi: 10.5194/angeo-19-275-2001
  * Franci et al. (2016) Franci, L., Landi, S., Matteini, L., Verdini, A., & Hellinger, P. 2016, The Astrophysical Journal, 833, 91
  * Hadid et al. (2015) Hadid, L., Sahraoui, F., Kiyani, K. H., et al. 2015, The Astrophysical Journal Letters, 813, L29
  * Hinze (1975) Hinze, J. 1975, Turbulence, McGraw-Hill classic textbook reissue (McGraw-Hill).  https://books.google.co.nz/books?id=DfRQAAAAMAAJ
  * Huang (2008) Huang, K. 2008, Statistical mechanics (John Wiley & Sons)
  * Huang et al. (2017) Huang, S., Hadid, L., Sahraoui, F., Yuan, Z., & Deng, X. 2017, The Astrophysical Journal Letters, 836, L10
  * Hundhausen et al. (1970) Hundhausen, A., Bame, S., Asbridge, J., & Sydoriak, S. 1970, Journal of Geophysical Research, 75, 4643
  * Isaacs et al. (2015) Isaacs, J. J., Tessein, J. A., & Matthaeus, W. H. 2015, Journal of Geophysical Research: Space Physics, 120, 868, doi: 10.1002/2014JA020661
  * Kolmogorov (1941) Kolmogorov, A. N. 1941, Dokl. Akad. Nauk SSSR, 30, 301, doi: 10.1098/rspa.1991.0075
  * Kraichnan (1965) Kraichnan, R. H. 1965, The Physics of Fluids, 8, 1385
  * Leamon et al. (1998) Leamon, R. J., Smith, C. W., Ness, N. F., Matthaeus, W. H., & Wong, H. K. 1998, Journal of Geophysical Research: Space Physics, 103, 4775, doi: https://doi.org/10.1029/97JA03394
  * Lepping et al. (1995) Lepping, R. P., Slavin, J. A., Ness, N. E., et al. 1995, Space Science Reviews, 71, 207
  * Linkmann et al. (2017) Linkmann, M., Berera, A., & Goldstraw, E. E. 2017, Phys. Rev. E, 95, 013102, doi: 10.1103/PhysRevE.95.013102
  * Linkmann et al. (2015) Linkmann, M. F., Berera, A., McComb, W. D., & McKay, M. E. 2015, Phys. Rev. Lett., 114, 235001, doi: 10.1103/PhysRevLett.114.235001
  * Lotz et al. (2023) Lotz, S., Nel, A. E., Wicks, R. T., et al. 2023, Astrophys. J., 942, 93, doi: 10.3847/1538-4357/aca903
  * Marsch (2006) Marsch, E. 2006, Living Reviews in Solar Physics, 3, 1
  * Matthaeus et al. (2005) Matthaeus, W. H., Dasso, S., Weygand, J. M., et al. 2005, Physical Review Letters, 95, 231101, doi: 10.1103/PhysRevLett.95.231101
  * Matthaeus et al. (2016) Matthaeus, W. H., Parashar, T. N., & Wan, M. 2016, The Astrophysical Journal Letters, 827, L7
  * Matthaeus et al. (2015) Matthaeus, W. H., Wan, M., Servidio, S., et al. 2015, Phil. Trans. R. Soc. A, 373, doi: 10.1098/rsta.2014.0154
  * Matthaeus et al. (2008) Matthaeus, W. H., Weygand, J. M., Chuychai, P., et al. 2008, The Astrophysical Journal, 678, 141
  * Mohapatra et al. (2020) Mohapatra, R., Federrath, C., & Sharma, P. 2020, Monthly Notices of the Royal Astronomical Society, 493, 5838
  * Mullan & Smith (2006) Mullan, D., & Smith, C. 2006, Solar Physics, 234, 325
  * Parashar et al. (2019) Parashar, T. N., Cuesta, M., & Matthaeus, W. H. 2019, The Astrophysical Journal, 884, L57, doi: 10.3847/2041-8213/ab4a82
  * Parashar et al. (2015) Parashar, T. N., Matthaeus, W. H., Shay, M. A., & Wan, M. 2015, The Astrophysical Journal, 811, 112, doi: 10.1088/0004-637X/811/2/112
  * Pearson et al. (2004) Pearson, B. R., Yousef, T. A., Haugen, N. E. L., Brandenburg, A., & Krogstad, P.-A. 2004, Phys. Rev. E, 70, doi: 10.1103/PhysRevE.70.056301
  * Perez et al. (2021) Perez, J. C., Bourouaine, S., Chen, C. H. K., & Raouafi, N. E. 2021, A&A, 650, A22, doi: 10.1051/0004-6361/202039879
  * Perri & Balogh (2010) Perri, S., & Balogh, A. 2010, Geophys. Res. Lett., 37, doi: 10.1029/2010GL044570
  * Perri et al. (2017) Perri, S., Servidio, S., Vaivads, A., & Valentini, F. 2017, The Astrophysical Journal Supplement Series, 231, 4, doi: 10.3847/1538-4365/aa755a
  * Phillips et al. (2022) Phillips, C., Bandyopadhyay, R., & McComas, D. J. 2022, The Astrophysical Journal, 933, 33, doi: 10.3847/1538-4357/ac713f
  * Pope (2000) Pope, S. B. 2000, Turbulent Flows (Cambridge University Press), doi: 10.1017/CBO9780511840531
  * Roberts et al. (2022) Roberts, O. W., Alexandrova, O., Sorriso‐Valvo, L., et al. 2022, Journal of Geophysical Research: Space Physics, 127, doi: 10.1029/2021JA029483
  * Roy et al. (2022) Roy, S., Chhiber, R., Dasso, S., Ruiz, M. E., & Matthaeus, W. H. 2022, Phys. Rev. E, 105, 045204, doi: 10.1103/PhysRevE.105.045204
  * Ruiz et al. (2014) Ruiz, M. E., Dasso, S., Matthaeus, W. H., & Weygand, J. M. 2014, Solar Physics, 289, 3917, doi: 10.1007/s11207-014-0531-9
  * Sahraoui et al. (2009) Sahraoui, F., Goldstein, M., Robert, P., & Khotyaintsev, Y. V. 2009, Physical review letters, 102, 231102
  * Schekochihin et al. (2009) Schekochihin, A., Cowley, S., Dorland, W., et al. 2009, The Astrophysical Journal Supplement Series, 182, 310
  * Sreenivasan (1998) Sreenivasan, K. R. 1998, Phys. Fluids, 10, 528
  * Taylor (1935) Taylor, G. I. 1935, Proc. Roy. Soc. Lond. A, 151, 421, doi: 10.1098/rspa.1935.0158
  * Tennekes (1975) Tennekes, H. 1975, Journal of Fluid Mechanics, 67, 561–567, doi: 10.1017/S0022112075000468
  * Tennekes & Lumley (1972) Tennekes, H., & Lumley, J. L. 1972, A first course in turbulence (Cambridge, MA: MIT press)
  * Van Atta & Antonia (1980) Van Atta, C., & Antonia, R. 1980, Physics of Fluids (1958-1988), 23, 252
  * Verma (1996) Verma, M. K. 1996, Journal of Geophysical Research: Space Physics, 101, 27543
  * Verma (2022) —. 2022, Physics of Plasmas, 29, 082902, doi: 10.1063/5.0096743
  * Veselovsky et al. (2010) Veselovsky, I., Dmitriev, A., & Suvorova, A. 2010, Cosmic Research, 48, 113
  * Wang et al. (2018) Wang, X., Tu, C., He, J., & Wang, L. 2018, Journal of Geophysical Research: Space Physics, 123, 68, doi: 10.1002/2017JA024813
  * Weygand et al. (2011) Weygand, J. M., Matthaeus, W. H., Dasso, S., & Kivelson, M. G. 2011, Journal of Geophysical Research: Space Physics, 116, doi: 10.1029/2011JA016621
  * Weygand et al. (2009) Weygand, J. M., Matthaeus, W. H., Dasso, S., et al. 2009, Journal of Geophysical Research: Space Physics, 114, doi: 10.1029/2008JA013766
  * Weygand et al. (2007) Weygand, J. M., Matthaeus, W. H., Dasso, S., Kivelson, M. G., & Walker, R. J. 2007, Journal of Geophysical Research: Space Physics, 112, doi: 10.1029/2007JA012486
  * Wilson III et al. (2021) Wilson III, L. B., Brosius, A. L., Gopalswamy, N., et al. 2021, Reviews of Geophysics, 59, e2020RG000714, doi: https://doi.org/10.1029/2020RG000714
  * Woodham et al. (2018) Woodham, L. D., Wicks, R. T., Verscharen, D., & Owen, C. J. 2018, The Astrophysical Journal, 856, 49, doi: 10.3847/1538-4357/aab03d
  * Wrench (2023) Wrench, D. 2023, reynolds_scales_project, Zenodo, doi: 10.5281/ZENODO.8352767
  * Yang et al. (2022) Yang, Y., Matthaeus, W. H., Roy, S., et al. 2022, The Astrophysical Journal, 929, 142
  * Yang et al. (2023) Yang, Y., Matthaeus, W. H., Oughton, S., et al. 2023, submitted.  https://arxiv.org/abs/2309.02663
  * Zhou & He (2021) Zhou, G., & He, H.-Q. 2021, The Astrophysical Journal Letters, 911, L2, doi: 10.3847/2041-8213/abef00
  * Zhou et al. (2020) Zhou, G., He, H.-Q., & Wan, W. 2020, The Astrophysical Journal Letters, 899, L32, doi: 10.3847/2041-8213/abaaa9
  * Zhou (2007) Zhou, Y. 2007, Physics of Plasmas, 14, 082701, doi: 10.1063/1.2739439
  * Zhou (2017) —. 2017, Physics Reports, 720-722, 1, doi: https://doi.org/10.1016/j.physrep.2017.07.005
  * Zhou (2021) —. 2021, Physics Reports, 935, 1, doi: https://doi.org/10.1016/j.physrep.2021.07.001
  * Zhou et al. (2014) Zhou, Y., Grinstein, F. F., Wachtor, A. J., & Haines, B. M. 2014, Phys. Rev. E, 89, 013303, doi: 10.1103/PhysRevE.89.013303
  * Zhou & Thornber (2016) Zhou, Y., & Thornber, B. 2016, J. Fluids Eng., 138, 070905, doi: 10.1115/1.4032532
  * Zhuravleva et al. (2019) Zhuravleva, I., Churazov, E., Schekochihin, A., et al. 2019, Nature Astronomy, 3, 832