1 Introduction
An outstanding problem in the study of astrophysical plasmas is the origin and sustenance of dynamically important magnetic fields found throughout the Universe. It is believed that, for systems without any large-scale features such as rotational shear or net helicity, weak seed fields are amplified and sustained by the fluctuation (or small-scale) dynamo (Kulsrud & Zweibel Reference Kulsrud and Zweibel2008), in which a random turbulent flow amplifies a weak seed magnetic field to near-equipartition energies (Batchelor Reference Batchelor1950; Zel'dovich et al. Reference Zel'dovich, Ruzmaikin, Molchanov and Sokolov1984; Childress & Gilbert Reference Childress and Gilbert1995; Rincon Reference Rincon2019; Tobias Reference Tobias2019). This process is thought to be particularly important in the intracluster medium (ICM), a hot (temperatures 
$T \sim 1\text {--}10$ keV) and diffuse (densities 
$n \sim 10^{-4}\text {--}10^{-2}\ \textrm {cm}^{-3}$) ionized plasma that sits in the dark-matter-dominated gravitational potential of a galaxy cluster (Fabian Reference Fabian1994). Observations of Faraday rotation and synchrotron emission (Carilli & Taylor Reference Carilli and Taylor2002; Bonafede et al. Reference Bonafede, Feretti, Murgia, Govoni, Giovannini, Dallacasa, Dolag and Taylor2010; Govoni & Feretti Reference Govoni and Feretti2004) show that intracluster magnetic fields have strengths 
${\sim }\mathrm {\mu }\textrm {G}$, with an energy density that is comparable to that of the bulk turbulent motions. The latter have recently been directly measured in the Perseus cluster by the Hitomi Collaboration (2016). As the ICM is not known to contain a significant mean field or net helicity, it is generally assumed that the present intracluster magnetic field is a result of the fluctuation dynamo.
Previous theoretical work on the fluctuation dynamo has typically adopted a collisional magnetohydrodynamic (MHD) description in the 
${Pm}\gtrsim {1}$ (‘high-Pm’) limit, where Pm is the ratio of the plasma viscosity and magnetic diffusivity. In these situations, the shearing motions of the turbulent background organize the magnetic field into folded sheets, whose lengths are comparable to the viscous scale and whose reversal scale is controlled by the magnetic diffusivity (Maron, Cowley & McWilliams Reference Maron, Cowley and McWilliams2004; Schekochihin et al. Reference Schekochihin, Cowley, Taylor, Maron and McWilliams2004b). This results in a magnetic field whose energy is concentrated near the resistive scale 
$k_\eta ^{-1}$ (Kazantsev Reference Kazantsev1968; Kulsrud & Anderson Reference Kulsrud and Anderson1992), which is 
${\sim }10^4\ \textrm {km}$ in the ICM assuming Spitzer values (see table 1 of Schekochihin & Cowley Reference Schekochihin and Cowley2006b). By contrast, observations suggest that the energy spectra of intracluster magnetic fields peak on scales 
${\sim }1\text {--}10\ \textrm {kpc}$ (Murgia et al. Reference Murgia, Govoni, Feretti, Giovannini, Dallacasa, Fanti, Taylor and Dolag2004; Vogt & Enßlin Reference Vogt and Enßlin2005; Vacca et al. Reference Vacca, Murgia, Govoni, Feretti, Giovannini, Perley and Taylor2012; Govoni et al. Reference Govoni, Murgia, Vacca, Loi, Girardi, Gastaldello, Giovannini, Feretti, Paladino and Carretti2017). While it has been proposed that, during the nonlinear phase of the dynamo, the magnetic energy would migrate to larger scales (Kulsrud et al. Reference Kulsrud, Cen, Ostriker and Ryu1997b; Subramanian Reference Subramanian1998) – thus potentially becoming compatible with ICM observations – this scenario has yet to be borne out in simulations.
Table 1. Index of runs, sorted into those using isotropic MHD (prefix ‘MHD’), unlimited Braginskii MHD (prefix ‘U’) and limited Braginskii MHD (prefix ‘L’). Run names adorned with an ‘S’ employ viscosities approaching and entering the ‘Stokes-flow’ regime (§ 4.4). Note: ‘Res.’ denotes the number of collocation points in each simulation, with the effective resolution reduced by a factor of 
$(2/3)^3$ due to de-aliasing. Viscosity and diffusivity values with an ‘(H)’ indicate simulations with hyper-diffusion with value 
$\mu _{{\rm H}}^{-1}$, 
$\eta _{{\rm H}}^{-1} =1.8\times 10^7$; ‘mirror’ denotes a simulation with a hard-wall limiter at the mirror threshold and no limiter at the firehose threshold. Time averages for 
$\langle u_\mathrm {rms}^2 \rangle ^{1/2}_t$ are taken over the kinematic stage.
The solution to this discrepancy may be related the fact that the ICM is not rigorously a MHD fluid (Schekochihin et al. Reference Schekochihin, Cowley, Kulsrud, Hammett and Sharma2005; Kulsrud & Zweibel Reference Kulsrud and Zweibel2008). First, the ion–ion Coulomb collision frequency
\begin{equation} \nu_{\textrm{i}} \approx 0.2 \bigg(\frac{n}{10^{-3}\ \textrm{cm}^{-3}}\bigg) \left(\frac{T}{5\ \textrm{keV}}\right)^{-3/2}\ \textrm{Myr}^{-1} \end{equation}
is only a factor of 
${\sim }100$ larger than the inverse dynamical time of the turbulent fluid motions at the largest scales,Footnote 1
\begin{equation} t^{-1}_\textrm{dyn} \approx 0.002 \bigg(\frac{u_0}{200~\textrm{km}\,\textrm{s}^{-1}}\bigg) \left(\frac{\ell_0}{100\ \textrm{kpc}}\right)^{-1}\ \textrm{Myr}^{-1} . \end{equation}
Thus, 
${\sim }1\,\%$ deviations in the particle distribution function from local thermodynamic equilibrium are to be expected – i.e. the ICM plasma is weakly collisional. That the energy density of these deviations is comparable to that stored in the observed turbulent motions and magnetic-field fluctuations indicates that 
${\sim }1\,\%$, while small, is nevertheless enough to be dynamically important. Secondly, already magnetic-field strengths as small as
\begin{equation} B \gtrsim 10^{-18} \bigg(\frac{n}{10^{-3}\ \textrm{cm}^{-3}}\bigg) \left(\frac{T}{5\ \textrm{keV}}\right)^{-3/2}\ \textrm{G} \end{equation}
are sufficient to ensure that the ICM plasma is magnetized – i.e. that the ion gyrofrequency 
$\varOmega _\textrm{i} \doteq eB/m_{\textrm{i}}c$ is larger than 
$\nu _{\textrm{i}}$, e is the elementary charge, m i is the ion mass and c is the speed of light. As seed magnetic fields are thought to be produced by various processes in the era preceding galaxy formation with magnitudes 
${\sim }10^{-22}\text {--}10^{-19}\ \textrm {G}$ (e.g. Biermann Reference Biermann1950; Kulsrud & Zweibel Reference Kulsrud and Zweibel2008), the amplification of the intracluster magnetic field via the fluctuation dynamo occurs almost exclusively in the weakly collisional, magnetized regime and is thus not a priori appropriately described by MHD with isotropic transport.
At magnetic-field strengths larger than (1.3), departures of the plasma from local thermodynamic equilibrium are biased with respect to the magnetic-field direction (Chew, Goldberger & Low Reference Chew, Goldberger and Low1956), and the transport of momentum and energy across magnetic-field lines becomes stifled by the smallness of the particles’ Larmor radii. In weakly collisional plasmas like the ICM, these departures manifest themselves inter alia as a pressure anisotropy 
${\rm \Delta} p \doteq p_\perp - p_\parallel$, where 
$p_\perp$ (
$p_\parallel$) is the thermal pressure across (along) the magnetic-field direction. This anisotropy exerts an anisotropic viscous stress on the fluid motion (Braginskii Reference Braginskii1965), which in an incompressible plasma primarily damps the parallel component of the rate of strain, 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol {:}}\, {\boldsymbol \nabla } {\boldsymbol {u}}$, where 
${\boldsymbol {\hat {b}}} \doteq {\boldsymbol {B}}/B$ is the unit vector in the direction of the magnetic field 
${\boldsymbol {B}}$ and 
${\boldsymbol {u}}$ is the flow velocity. In the limit where 
$\nu _{\textrm{i}}\gg {\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol {:}}\, {\boldsymbol \nabla } {\boldsymbol {u}}$, this pressure anisotropy can be expressed as a balance between adiabatic production (via changes in the magnetic-field strength, since the Lagrangian time rate of change 
$\textrm {d}\ln B/\textrm {d}t = {\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol {:}}\, {\boldsymbol \nabla } {\boldsymbol {u}}$ in a perfectly conducting, incompressible plasma) and collisional isotropization of the pressure tensor, viz.
\begin{equation} {\rm \Delta} p = 3 \mu_{{\rm B}} {\boldsymbol{\hat{b}}}{\boldsymbol{\hat{b}}}\, {\boldsymbol{:}}\, {\boldsymbol \nabla} {\boldsymbol{u}}, \end{equation}
where 
$\mu _{{\rm B}} = p/\nu _{\textrm{i}}$ is the field-aligned Braginskii viscosity (Braginskii Reference Braginskii1965).Footnote 2
In practice, the Braginskii parallel viscous stress
\begin{equation} {\boldsymbol {\mathsf {\Pi}}}_\parallel \doteq -\left({\boldsymbol{\hat{b}}}{\boldsymbol{\hat{b}}} - \frac{1}{3} {\boldsymbol{\mathsf {I}}} \right) {\rm \Delta} p \end{equation}
with the pressure anisotropy 
${\rm \Delta} p$ given by (1.4) is only suitable for plasmas with small to order-unity values of 
$\beta \doteq 8 {\rm \pi} p/ B^2$, the ratio of the thermal and magnetic pressures. The reason is that high-
$\beta$, weakly collisional plasmas are susceptible to several kinetic instabilities when 
$|{\rm \Delta} p/p| \gtrsim 1/\beta$, such as the firehose (Rosenbluth Reference Rosenbluth1956; Chandrasekhar, Kaufman & Watson Reference Chandrasekhar, Kaufman and Watson1958; Hellinger & Matsumoto Reference Hellinger and Matsumoto2000) and mirror (Barnes Reference Barnes1966; Hasegawa Reference Hasegawa1969; Southwood & Kivelson Reference Southwood and Kivelson1993; Hellinger Reference Hellinger2007) instabilities. These instabilities regulate the pressure anisotropy to values near the instability thresholds, an effect that has been diagnosed in various kinetic particle-in-cell simulations (Kunz, Schekochihin & Stone Reference Kunz, Schekochihin and Stone2014; Hellinger & Trávníček Reference Hellinger and Trávníček2015; Riquelme, Quataert & Verscharen Reference Riquelme, Quataert and Verscharen2015; Melville, Schekochihin & Kunz Reference Melville, Schekochihin and Kunz2016) and directly observed using in situ measurements of particle distribution functions and magnetic fluctuations in the solar wind (Kasper, Lazarus & Gary Reference Kasper, Lazarus and Gary2002; Hellinger et al. Reference Hellinger, Trávníček, Kasper and Lazarus2006; Bale et al. Reference Bale, Kasper, Howes, Quataert, Salem and Sundkvist2009; Chen Reference Chen2016; Chen et al. Reference Chen, Matteini, Schekochihin, Stevens, Salem, Maruca, Kunz and Bale2016). Thus, because the magnitudes of pressure anisotropy (and thus parallel viscosity) specified by (1.4) are often unphysically large in weakly collisional, high-
$\beta$ plasmas, any fluid model of the fluctuation dynamo must adopt some form of microphysical closure to account for the enhanced (but not necessarily complete – see § 1.1.3) regulation of the pressure anisotropy due to microscale kinetic instabilities.
1.1 Theoretical considerations on pressure anisotropy and dynamo
1.1.1 Hard-wall limiters
One popular choice of closures in fluid simulations of weakly collisional, high-
$\beta$ plasmas is to limit, by hand, the pressure anisotropy to remain within the firehose- and mirror-instability thresholds
\begin{equation} -\frac{B^2}{4{\rm \pi}}\le {\rm \Delta} p \le \frac{B^2}{8{\rm \pi}}. \end{equation}The resulting ‘hard-wall’ limiters, which have their origin in pioneering work on the kinetic magnetorotational instability by Sharma et al. (Reference Sharma, Hammett, Quataert and Stone2006) and have since been used in Braginskii-MHD simulations of magnetothermal (Kunz et al. Reference Kunz, Bogdanović, Reynolds and Stone2012) and magnetorotational (Kempski et al. Reference Kempski, Quataert, Squire and Kunz2019) turbulence, take the form
\begin{equation} {\rm \Delta} p = {\min}\left(\frac{B^2}{8{\rm \pi}}, 3\mu_{{\rm B}}{\boldsymbol{\hat{b}}}{\boldsymbol{\hat{b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla} {\boldsymbol u}\right) \end{equation}
on the mirror (
${\rm \Delta} p>0$) side and
\begin{equation} {\rm \Delta} p = {\max}\left(-\frac{B^2}{4{\rm \pi}}, 3\mu_{{\rm B}}{\boldsymbol{\hat{b}}}{\boldsymbol{\hat{b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla} {\boldsymbol u}\right) \end{equation}
on the firehose (
${\rm \Delta} p<0)$ side (again, assuming incompressibility). Similarly effective limiters, in the form of a large anomalous collision frequency enacted in regions of firehose/mirror instability, were employed by Santos-Lima et al. (Reference Santos-Lima, de Gouveia Dal Pino, Kowal, Falceta-Gonçalves, Lazarian and Nakwacki2014) in simulations of turbulent dynamo using the double-adiabatic Chew et al. (Reference Chew, Goldberger and Low1956) equations.
1.1.2 Effective Reynolds number
Such pressure-anisotropy limiters effectively imply an enhanced collisionality in the unstable regions given by
\begin{equation} \nu_\mathrm{eff} \sim \beta{\boldsymbol{\hat{b}}}{\boldsymbol{\hat{b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla} {\boldsymbol u} . \end{equation}
Following the reasoning presented in § 4.2.2 of Melville et al. (Reference Melville, Schekochihin and Kunz2016), we estimate the effective parallel viscosity 
$\mu _{\parallel \textrm {eff}} \doteq v_{\textrm {th{i}}}^2/\nu _\mathrm {eff}$, and thus the effective parallel-viscous Reynolds number 
${Re}_\parallel$ associated with the enhanced collisionality (1.9), as follows. (Here, 
$v_\mathrm {thi} \doteq \sqrt {2 T_{{\rm i}} / m_{{\rm i}}}$ is the ion thermal velocity and 
$m_{{\rm i}}$ is the ion mass.) The Kolmogorov scaling 
$u_{\ell _\parallel } \propto \ell _\parallel ^{1/3}$ for the field-stretching turbulent velocity at parallel scale 
$\ell _\parallel$ implies that the magnitude of the field-parallel rate of strain 
$|{\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}| \sim u_{\ell _\parallel }/{\ell _\parallel } \propto {\ell _\parallel }^{-2/3}$ is largest at the effective parallel-viscous scale 
$\ell _{\mu _{\parallel \textrm {eff}}}$, where such motions are dissipated. The value of 
${Re}_{\parallel \mathrm {eff}}$ corresponding to (1.9) is then
\begin{equation} {Re}_{\parallel \mathrm{eff}} \doteq \frac{u_0\ell_0}{\mu_{\parallel\textrm{eff}}} \sim u_0 \ell_0 \frac{\nu_\mathrm{eff}}{v_\mathrm{thi}^2}\sim \beta M^2 {Re}_{\parallel \mathrm{eff}}^{1/2} \quad \Longrightarrow\quad {Re}_{\parallel \mathrm{eff}} \sim \beta^2 M^4 , \end{equation}
where 
$u_0$ is the characteristic speed of the outer-scale bulk fluid motions, 
$\ell _0$ is the energy injection (outer) scale and 
$M \doteq u_0/v_\mathrm {thi}$ is the Mach number. The resulting effective particle mean free path is
\begin{equation} \lambda_{\textrm{mfp},\textrm{eff}} \doteq \frac{v_{\textrm{th{i}}}}{\nu_\textrm{eff}} \sim \ell_0 \beta^{-2} M^{-3} . \end{equation}The result is a parallel-viscous cutoff
\begin{equation} \ell_{\mu_{\parallel\textrm{eff}}} \sim \ell_0 \,{Re}_{\parallel \mathrm{eff}}^{-3/4} \sim \ell_0 \beta^{-3/2} M^{-3} \end{equation}
that is smaller than the cutoff caused by Coulomb collisions, 
$\ell _\textrm {Coulomb} \sim \ell _0 \,{Re}^{-3/4}_\parallel$ with 
${Re}_{\parallel } \doteq u_0 \ell _0 / ({v_{\textrm {thi}}}^2/\nu _{{\rm i}})$, if the magnetic-field strength satisfies
\begin{equation} B \lesssim 2 \bigg(\frac{n}{10^{-3}\ \textrm{cm}^{-3}}\bigg)^{1/4} \left(\frac{T}{5\ \textrm{keV}}\right)\left(\frac{M}{0.2}\right)^{3/4} \left(\frac{\ell_0}{100\ \textrm{kpc}}\right)^{-1/4}\ \mathrm{\mu}\mathrm{G} . \end{equation}(In (1.13), we have chosen a normalization based on typical values found in the bulk ICM.) In this case, a larger maximal shear rate can be achieved:
\begin{equation} S_{\mu_{\parallel\textrm{eff}}} \sim \frac{u_0}{\ell_0} \, {Re}_{\parallel \mathrm{eff}}^{1/2} \sim \frac{u_0}{\ell_0} \beta M^2 . \end{equation}
Note that 
$\ell _{\mu _{\parallel \textrm {eff}}}\propto B^3$ increases and 
$S_{\mu _{\parallel \textrm {eff}}} \propto B^{-2}$ decreases with increasing magnetic-field strength. This suggests that the dynamo subject to ‘hard-wall’ limiters is likely faster than its unlimited counterpart, but it slows down on average as the magnetic energy grows.
1.1.3 Suitability of the hard-wall limiters
The use of hard-wall limiters implicitly assumes that the enhanced collisionality given by (1.9) is achievable in the plasma, viz., that the collision frequency needed to pin the pressure anisotropy to either marginal-stability threshold, 
${\sim }\beta {\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$ (see (1.9)), is smaller than the ion-Larmor frequency 
$\varOmega _{{\rm i}}$. Using (1.14) as an estimate for the maximal value of 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$, this criterion becomes
\begin{equation} \beta \gtrsim \bigg(\frac{\ell_0}{d_{{\rm i}}}\bigg)^{2/5} M^{-6/5} , \end{equation}
where 
$d_{{\rm i}}$ is the (field-strength-independent) ion inertial scale. Equation (1.15a) is referred to as the ‘ultra-high-
$\beta$’ limit in Melville et al. (Reference Melville, Schekochihin and Kunz2016). For characteristic parameters of the ICM, (1.15a) implies
\begin{equation} B \lesssim 6 \bigg(\frac{n}{10^{-3}\ \textrm{cm}^{-3}}\bigg)^{2/5} \left(\frac{T}{5\ \textrm{keV}}\right)^{1/2}\left(\frac{M}{0.2}\right)^{3/5} \left(\frac{\ell_0}{100\ \textrm{kpc}}\right)^{-1/5}\ \mathrm{nG}. \end{equation}For magnetic-field strengths satisfying (1.15b), the hard-wall limiters are unsuitable and the parallel Reynolds number is smaller than that given by (1.10).
1.1.4 Three dynamo regimes
The conditions (1.3) and (1.15) suggest that there exist three distinct regimes of magnetization and collisionality:
(i) the unmagnetized regime, when
$B \lesssim 10^{-18}$ G;(ii) the magnetized ‘kinetic’ regime (ultra-high-
$\beta$), when $10^{-18}\ \textrm {G} \lesssim B \lesssim 6\ \textrm {nG}$ and for which the hard-wall limiters (1.7) and (1.8) are unsuitable; and(iii) the magnetized ‘fluid’ regime, when
$B \gtrsim 6\ \textrm {nG}$ and for which the pressure anisotropy can be well regulated by the instabilities (i.e. $\nu _\textrm {eff} \lesssim \varOmega _{\rm i}$).
The nomenclature of ‘kinetic’ and ‘fluid’ used here is a reflection of our result (see § 1.2) that, when the pressure anisotropy 
${\rm \Delta} p/p$ is limited by the firehose and mirror instabilities to values 
${\sim }\beta ^{-1}\ll {1}$ (regime (iii)), the plasma dynamo largely resembles its fluid (MHD) counterpart, particularly when the magnetic field is dynamically weak. If instead the pressure anisotropy is poorly regulated by the instabilities (regime (ii)), the associated parallel viscous stress is dynamically important throughout the entire evolution of the dynamo, a feature that qualitatively changes the character of the fluctuation dynamo.
The saturated state of the dynamo, in which the magnetic and kinetic energies are comparable, would be obtained when
\begin{equation} B_\textrm{sat} \sim 3 \bigg(\frac{n}{10^{-3}\ \textrm{cm}^{-3}}\bigg)^{1/2} \left(\frac{T}{5\ \textrm{keV}}\right)^{1/2}\left(\frac{M}{0.2}\right)\mathrm{\mu}\textrm{G} \end{equation}
and thus occurs in the magnetized fluid regime. Perhaps coincidentally, this value is close to the field strength at which the effective scattering frequency due to efficient firehose/mirror is comparable to that due to Coulomb collisions, (1.13). Namely, for 
$B\gtrsim 2\ \mathrm {\mu }\textrm {G}$, the effective viscous length scale given by (1.12) satisfies 
$\ell _{\mu _{\parallel \textrm {eff}}} \lesssim \ell _\mathrm {Coulomb}$, and so Coulomb collisions are frequent enough to maintain a kinetically stable pressure anisotropy. At smaller magnetic-field strengths, the mean free path is effectively set by the kinetic instabilities (see (1.11)).
An intriguing feature of these regimes is that, while the Reynolds number in the unmagnetized regime is set by Coulomb collisions, resulting in 
${Re}_\parallel \sim 1 \text {--} 100$ (Schekochihin et al. Reference Schekochihin, Cowley, Kulsrud, Hammett and Sharma2005), at the transition from the second (magnetized kinetic) regime to the third (magnetized fluid) regime we find 
${Re}_{\parallel \mathrm {eff}} \sim \beta ^2 M^4 \gg 1$. This suggests that 
${Re}_{\parallel \mathrm {eff}}$ must experience a large increase at some time between these two epochs. Since the viscous-scale rate of strain increases as 
${Re}_{\parallel \mathrm {eff}}^{1/2}$, the dynamo in this intermediate second regime should, at some point, be self-accelerating, with the field-stretching eddies becoming smaller and faster as the magnetic field is amplified. Similar scenarios have previously formed the basis for theories of explosive dynamo in collisionless plasmas (Schekochihin & Cowley Reference Schekochihin and Cowley2006a,Reference Schekochihin and Cowleyb; Mogavero & Schekochihin Reference Mogavero and Schekochihin2014; Melville et al. Reference Melville, Schekochihin and Kunz2016), a topic that will be the subject of a separate publication.Footnote 3
1.1.5 Previous simulations and numerical results
The first two dynamo regimes have been previously studied through the use of hybrid-kinetic numerical simulations by Rincon et al. (Reference Rincon, Califano, Schekochihin and Valentini2016) and St-Onge & Kunz (Reference St-Onge and Kunz2018, hereafter Paper I). These studies observed the generation of firehose and mirror instabilities as the dynamo entered the magnetized regime, aiding in breaking adiabatic invariance and regulating the pressure anisotropy. Santos-Lima et al. (Reference Santos-Lima, de Gouveia Dal Pino, Kowal, Falceta-Gonçalves, Lazarian and Nakwacki2014) studied the effects of pressure-anisotropy regulation during the fluctuation dynamo using a collisionless double-adiabatic closure to evolve 
$p_\perp$ and 
$p_\parallel$ (Chew et al. Reference Chew, Goldberger and Low1956) supplemented by a non-zero collision frequency 
$\nu _{\textrm {eff}}$ that was activated in spatial regions of kinetic instability. It was found that simulations with instantaneous pressure-anisotropy relaxation exhibited magnetic-field growth rates similar to those in isotropic MHD, with some minor details differing in the saturated state. As the collision frequency was lowered, the dynamo growth rates and the final value of the saturated magnetic energy decreased. In the entirely collisionless case, the pressure anisotropy was allowed to grow arbitrarily large, and no growth of the magnetic energy was observed. This is consistent with theoretical considerations of magnetic-field amplification occurring under adiabatic invariance in collisionless plasmas by Helander, Strumik & Schekochihin (Reference Helander, Strumik and Schekochihin2016), who concluded that dynamo action always requires collisions or some small-scale kinetic mechanism for breaking the adiabatic invariance of the magnetic moment (cf. Kulsrud et al. Reference Kulsrud, Cowley, Gruzinov and Sudan1997a).
At the moment, the third (magnetized fluid) regime of the dynamo is prohibitively expensive to investigate using kinetic simulation, except near the saturated state (which was the focus of § 3.3 in Reference St-Onge and KunzPaper I). To appreciate this difficulty, let us imagine that one wishes to resolve two orders of magnitude of magnetic-energy growth (equivalently, one order of magnitude in growth of the magnetic-field strength) in this regime in a single simulation. The constraints on the initial plasma beta, 
$\beta _0$, required to simulate this regime are, in terms of the controllable simulation parameters,
\begin{equation} \frac{A}{M^2} \lesssim \beta_0 \lesssim \frac{1}{M^{3/2}} \left(\frac{\ell_0}{{\rho_{{\rm i}0}}}\right)^{1/2}, \end{equation}
where 
$A$ is the desired magnetic-energy amplification factor and 
$ {\rho _{{\rm i}0}}$ is the ion-Larmor radius at the transition from the second to this third regime. The first inequality in (1.17) follows from the dynamo not yet being saturated (i.e. 
$\beta M^{-2} \lesssim 1$), while the second inequality follows from the requirement 
$\nu _\textrm {eff} \lesssim \varOmega _{{\rm i}0}$. For such a range of 
$\beta _0$ to exist,
\begin{equation} A \lesssim M^{1/2} \left(\frac{\ell_0}{{\rho_{{\rm i}0}}}\right)^{1/2} . \end{equation}
To allow for an appreciable range of field amplification in this regime, we must then maximize 
$M$ and 
$\ell _0/ {\rho _{{\rm i}0}}$. If we demand 
$M \lesssim 0.2$ (so as to maintain relevance to the sub-sonic turbulence observed in the ICM), then resolving two orders of magnitude of energy growth (
$A \sim 10^2$) requires 
$\ell _0/ {\rho _{{\rm i}0}} \gtrsim 50\,000$. If we then wish to resolve the ion-Larmor radius at the end of this growth period by a minimum of two cells, then the number of cells needed in each spatial direction is 
${\sim }10^6$! For a problem that is intrinsically three-dimensional (Cowling Reference Cowling1933; Zel'dovich Reference Zel'dovich1957), this requirement is well beyond current computational capabilities.
1.1.6 This work
The goal of the work presented here is twofold. First, we study the plasma dynamo in the third (magnetized fluid) regime using a reduced set of fluid equations that takes the effects of kinetic instabilities into account through the use of microphysical closures that limit the pressure anisotropy. For this, we use the incompressible MHD equations including the parallel component of the Braginskii viscosity tensor, the latter being hard-wall limited to kinetically stable values (see (1.6) and § 2.1). Secondly, in order to understand better the effects of this regulation and the nature of the dynamo in the second (magnetized kinetic) regime, we also perform a number of Braginskii-MHD simulations without any hard-wall regulation of the pressure anisotropy, relying instead on incomplete regulation by a constant (e.g. Coulomb) collisionality. This approach is justified by measurements of imperfect pressure-anisotropy regulation during this regime in the hybrid-kinetic simulations presented in Paper I.
1.2 Main results and organization of this paper
Given the rather technical nature of the presentation that follows, we first provide an overview of our main results as a preamble, with pointers to the relevant sections and figures in the main text. We arrive at three main conclusions (which are re-summarized in § 5 alongside a recapitulation of ancillary results and a preview of future work):
(i) With hard-wall limiters on the pressure anisotropy that prevent
${\rm \Delta} p$ from growing beyond its kinetically stable values, the Braginskii-MHD dynamo is in most respects identical to the standard high-${Pm}$ MHD dynamo. This is to be expected, because the limited pressure anisotropy becomes dynamically important only once the Lorentz force does, viz., as the dynamo starts saturating. Certain minor differences with isotropic MHD dynamo do appear in the saturated state.(ii) When no pressure-anisotropy limiters are used (as relevant to regimes in which an effective collision frequency
$\nu _\textrm {eff} \gtrsim \varOmega _{\rm i}$ would be required to keep ${\rm \Delta} p$ near marginal stability), the Braginskii-MHD dynamo behaves quite differently to the MHD dynamo. Nearly all of these differences can be understood by noting that, in its growing phase, the structure and statistics of the magnetic field are remarkably similar to those found in the saturated state of the (high-${Pm}$) MHD dynamo. This occurs because the form of the Braginskii viscous stress, ${\boldsymbol \nabla }\,{\boldsymbol {\cdot }}\,({\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}{\rm \Delta} p)$ (the third term on the right-hand side of (2.1a)), is mathematically similar to that of the Lorentz force, ${\boldsymbol {B}}\,{\boldsymbol {\cdot }}\,{\boldsymbol \nabla }{\boldsymbol {B}}={\boldsymbol \nabla }\,{\boldsymbol {\cdot }}\,({\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}} B^2/4{\rm \pi} )$ (the first term on the right-hand side of (2.1a)), if one makes the substitution $B^2/4{\rm \pi} \rightarrow {\rm \Delta} p \sim \mu _{{\rm B}} \, \textrm {d}\ln B/\textrm {d} t = \mu _{{\rm B}} \, {\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$ (neglecting resistivity).(iii) Without pressure-anisotropy limiters, Braginskii MHD no longer supports a dynamo if the ratio of the isotropic viscosity (
$\mu$) and the Braginskii viscosity ($\mu _{{\rm B}}$) is too small. This behaviour may be understood heuristically by noting that the Braginskii viscosity, by targeting only those fluid motions that stretch the magnetic field (${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}\ne 0$), curbs the growth of that field while simultaneously promoting its resistive decay by allowing motions that mix the field lines. Finite isotropic viscosity moderates the mixing motions, thereby allowing the field to grow if this viscosity is sufficiently large. In the limit where the Braginskii viscosity is so strong that the outer-scale fluid motions become two-dimensionalized with respect to the magnetic-field direction, the dynamo shuts down unless the isotropic viscosity is large enough to diffuse velocity gradients into the field-perpendicular direction, thus once again rendering the dynamics three-dimensional.
The outline of this paper is as follows. We present the Braginskii-MHD equations and some definitions of relevant dimensionless numbers in § 2.1. We then describe the numerical method of solution and several useful diagnostics in §§ 2.2 and 2.3. The full list of runs is given in table 1. We open the results section (§ 3) with a brief overview of the fluctuation dynamo, broken down into its four evolutionary stages (§ 3.1). We then present evidence from our simulations in favour of our first conclusion, that the limited Braginskii-MHD dynamo is similar to the isotropic-MHD dynamo in a 
$Re\gg {1}$, 
${Pm}\gtrsim {1}$ fluid (§ 3.2). This discussion covers the basic structure of the fields and flows in the kinematic and saturated stages (see figures 1, 3–6) and the dynamo growth rates (figure 2). Section 3.3 provides evidence for our second conclusion, that the structure and statistics of the unlimited Braginskii-MHD dynamo imitate those in the saturated state of the more standard 
${Pm}\gtrsim {1}$ MHD dynamo (see figures 7–13). Much of the rest of the paper is devoted to exploring the consequences of this similarity. This leads us to a treatment of the anisotropization of the fluid flow by the Braginskii viscosity (§ 4), in which we study the structure of the rate-of-strain tensor (figure 14) and the alignment of the eigenvectors of the rate-of-strain tensor with the magnetic-field direction 
${\boldsymbol {\hat {b}}}$ (figures 15–17).Footnote 4 These studies also reveal minor differences in the properties of the saturated state of the limited Braginskii-MHD dynamo compared to isotropic MHD, alluded to in conclusion (i). Motivated by these results, in § 4.3 we formulate an analytic model for the kinematic stage of the unlimited Braginskii-MHD dynamo, which is based on an extension of the Kazantsev–Kraichnan model to anisotropic magnetic-field statistics first proposed to describe the saturated MHD dynamo by Schekochihin et al. (Reference Schekochihin, Cowley, Taylor, Hammett, Maron and McWilliams2004a) (a detailed derivation of the model, which was omitted in their paper, is given in appendix C). The model, although necessarily simplistic, appears to give a reasonable match to the spectra (figure 18) seen in simulations, and predicts that the dynamo shuts off for sufficiently small isotropic viscosities (figure 19). It also provides a reasonable explanation for the behaviour of a separate set of simulations in the Braginskii ‘Stokes-flow’ regime (
${Re}_{\parallel \mathrm {eff}}\lesssim 1$), described in § 4.4, in which we see the dynamo shut off for sufficiently small 
${Re}_{\parallel \mathrm {eff}}\lesssim 1$ (unless the isotropic viscosity is also very large, in which case the dynamo can operate). These results support our conclusion (iii). We finish with a more thorough discussion of our conclusions, as well as suggestions for future work, in § 5.
Figure 1. A two-dimensional cross-section of the magnetic field and flow in a hard-wall-limited Braginskii-MHD simulation (run L1) and in a comparable MHD simulation (run MHD4). The left, centre and right panels display the magnetic-field strength, velocity magnitude and parallel rate of strain, respectively. Panels (a) and (b) display results from run L1 (parameters 
$\mu _{{\rm B}}^{-1}=20$, 
$\mu ^{-1}=\eta ^{-1}=1500$): (a) the kinematic stage; (b) the saturated state. Panels (c) and (d) display results from run MHD4 (parameters 
$\mu ^{-1}=\eta ^{-1}=1500$): (c) the kinematic stage; (d) the saturated state. All plots are on a linear scale, with brighter regions denoting higher magnitudes.
Figure 2. Evolution of the kinetic (dash-dotted lines) and magnetic (solid lines) energies for simulations employing hard-wall pressure-anisotropy limiters with varying Braginskii viscosity 
$\mu_{{\rm B}}$ and fixed isotropic viscosity 
$\mu$ (left) and varying isotropic viscosity 
$\mu$ and fixed Braginskii viscosity 
$\mu _{{\rm B}}$ (right). A Braginskii-MHD simulation employing hard-wall limiters only on the mirror side, marked with an ‘(
${\rm M}$)’, is included in the left panel.
2 Method of solution
2.1 Equations and free parameters
The incompressible MHD equations including the field-aligned component of the Braginskii viscosity are
\begin{align} \frac{\textrm{d} {\boldsymbol{u}}}{\textrm{d} t} &\doteq \left(\frac{\partial }{\partial t} + {\boldsymbol{u}}\,{\boldsymbol{\cdot}}\,{\boldsymbol \nabla}\right){\boldsymbol{u}} = {\boldsymbol{B}}\,{\boldsymbol{\cdot}}\,{\boldsymbol \nabla}{\boldsymbol{B}} - {\boldsymbol \nabla} p + {\boldsymbol \nabla} \,{\boldsymbol{\cdot}}\,({\boldsymbol{\hat{b}}}{\boldsymbol{\hat{b}}}{\rm \Delta} p) - (-1)^h \mu_h \nabla^{2h}{\boldsymbol{u}}+ \tilde{{\boldsymbol f}}, \end{align}
\begin{align}\frac{\textrm{d} {\boldsymbol{B}}}{\textrm{d} t} &\doteq \left(\frac{\partial }{\partial t} + {\boldsymbol{u}}\,{\boldsymbol{\cdot}}\,{\boldsymbol \nabla}\right){\boldsymbol{B}} = {\boldsymbol{B}}\,{\boldsymbol{\cdot}}\, {\boldsymbol \nabla} {\boldsymbol{u}} - (-1)^h \eta_h\nabla^{2h} {\boldsymbol{B}} , \end{align}
where the magnetic field 
${\boldsymbol {B}}$ is expressed in Alfvénic units and the mass density has been scaled out. The last term on the right-hand side of (2.1a), 
$\tilde {{\boldsymbol f}}$, is a random driving body force, further described in § 2.2. The additional diffusive terms in (2.1), featuring 
$\mu _h$ and 
$\eta _h$, are Laplacian (
$h=1$) or hyper (
$h=2$) viscosity and diffusivity, respectively; these are introduced to truncate the cascades of kinetic and magnetic energy near the smallest wavelengths captured in our simulations. The pressure anisotropy is given by (1.4) except when limited by heuristic micro-instability limiters, as per (1.7) and (1.8).Footnote 5 Because the ratio of initial magnetic energy to kinetic energy is extremely small, much of the dynamo evolution occurs with the majority of the plasma residing in the firehose- or mirror-unstable regions. Depending upon the material properties of the plasma, this ratio can remain small even in the saturated state. As a result, the pressure-anisotropy limiters are likely to play an important dynamical role throughout the entire evolution of the plasma dynamo.
Equations (2.1) have four free parameters: the two isotropic diffusivities 
$\mu _h$ and 
$\eta _h$, the anisotropic Braginskii viscosity 
$\mu _{{\rm B}}$ and the specifications of the random forcing 
$\tilde {{\boldsymbol {f}}}$, described in § 2.2. For our analysis, it is useful to define the following Reynolds and Prandtl numbers:
\begin{equation} {Re} \doteq \frac{u_0 \ell_0}{\mu}, \quad {Re}_\parallel \doteq \frac{u_0 \ell_0}{\mu_{{\rm B}}}, \quad {Rm} \doteq \frac{u_0 \ell_0}{\eta}, \quad {Pm} \doteq \frac{{Rm}}{{Re}}, \end{equation}
where 
$u_0$ is the root-mean-square fluid velocity, 
$\ell _0 \doteq L/ 2 {\rm \pi} \sim k_{f}^{-1}$, 
$L$ is the size of the simulation domain and 
$k_{f}$ is the forcing wavenumber. In the unlimited case, the effective parallel-viscous Reynolds number introduced in (1.10) satisfies 
${Re}_{\parallel \mathrm {eff}} = {Re}_\parallel$; in the limited case, 
${Re}_\parallel \le {Re}_{\parallel \mathrm {eff}}\le {Re}$. The definitions (2.2a–d) are suitable for Laplacian dissipation, but generalized Reynolds numbers 
${Re}_h$ (
${Rm}_h$) can be formulated for higher-order dissipation. To do so, we make the substitutions 
$\mu \rightarrow \ell _\mu ^{-2(h-1)}\mu _h$ and 
$\eta \rightarrow \ell _\eta ^{-2(h-1)}\eta _{h}$ in the standard definitions (2.2a–d). This replacement is done so that, for a given value of 
${Re}_h$ (
${Rm}_h$), 
$\ell _\mu /\ell _0$ (
$\ell _\eta /\ell _0$) is held fixed for all values of 
$h$. This allows one to compare two systems by directly comparing their generalized Reynolds numbers. Assuming Kolmogorov scalings (viz. 
$|{\boldsymbol \nabla } {\boldsymbol {u}}| \propto {Re}_h^{2/3}$) to compute the dissipation scales 
$k_\mu$ and 
$k_\eta$, we have
\begin{gather} {Re}_h \doteq \frac{u_0 \ell_0}{\ell_\mu^{-2(h-1)} \mu_{h}} =\left(\frac{u_0 \ell_0^{2h-1}}{\mu_h}\right)^{2/(3h-1)}, \end{gather}
\begin{gather}{Rm}_h \doteq \frac{u_0\ell_0}{\ell_\eta^{-2(h-1)} \eta_{h}} = \left(\frac{ u_0 \ell_0^{2h-1}}{\eta_h {Re}_{\parallel \mathrm{eff}}^{(h-1)/2}}\right)^{1/h}. \end{gather}While more general definitions can be devised without assuming Kolmogorov scalings, they will generically depend on some characteristic of the underlying fields that must be determined a posteriori.
2.2 Numerical simulations
To solve (2.1), we use a version of the pseudospectral incompressible-MHD code Snoopy (Lesur & Longaretti Reference Lesur and Longaretti2007) equipped with the parallel Braginskii viscosity. All simulations are run on a triply periodic, 
$224^3$ grid with 
$2/3$ de-aliasing (with the exception of the ‘Stokes-flow’ runs discussed in § 4.4, which use 
$112^3$ collocation points).Footnote 6 Turbulence is driven using an incompressible, zero-mean-helicity,Footnote 7 random forcing at wavenumbers 
$k_{f} \in [2 {\rm \pi} /L, 4{\rm \pi} /L]$; its power is distributed evenly among these Fourier modes. The forcing is time-correlated by an Ornstein–Uhlenbeck process. Code units are based on the box size 
$L=1$ and energy injection rate 
$\varepsilon =1$. The latter appears in the Ornstein–Uhlenbeck process via
\begin{equation} \tilde{{\boldsymbol f}}(t+{\rm \Delta} t) = \tilde{{\boldsymbol f}}(t) \exp({-{\rm \Delta} t/{t_{\textrm{corr},f}}}) + \left[\frac{\varepsilon}{{t_{\textrm{corr},f}}}\left(1-\exp({-2 {\rm \Delta} t /{t_{\textrm{corr},f}}})\right)\right]^{1/2}{\boldsymbol{\tilde{g}}}, \end{equation}
where 
${\rm \Delta} t$ is the simulation timestep, 
${t_{\textrm {corr},f}}$ is the correlation time of the forcing and 
${\boldsymbol {\tilde {g}}}$ is a Gaussian noise at wavenumber 
$k_{f}$ (Gillespie Reference Gillespie1996). These normalizations of 
$L$ and 
$\varepsilon$ lead to saturated turbulence amplitudes of order unity (
$u_\mathrm {rms} \sim 1$). The correlation time is chosen to be 
${t_{\textrm {corr},f}}=(2{\rm \pi} )^{-1} \approx \ell _0/u_\textrm {rms}$, which is comparable to the inverse decorrelation rate at the outer scale for 
${Re} \ge 1$ turbulence. Simulations are initiated with 
${\boldsymbol {u}} =0$ and a random zero-net-flux magnetic field with power at 
$k \in [2 {\rm \pi} /L, 4{\rm \pi} /L]$. All runs have an initial root-mean-square (r.m.s.) field strength 
$B_{0,\textrm {rms}} = 10^{-3}$. Because we have assumed incompressibility, the thermal velocity 
$v_{{\rm thi}}$, and thus 
$\beta$, have been eliminated from the equations. Accordingly, we formulate the stability thresholds (1.6) in terms of 
${\rm \Delta} p$ and 
$B^2$ directly and subsume 
$p$ and 
$\nu _{{\rm i}}$ into the definition of 
$\mu _{{\rm B}}$.
For comparisons, we have also performed simulations of the dynamo in standard incompressible MHD with isotropic diffusion (i.e. 
$\mu _{{\rm B}}=0$) for Reynolds numbers ranging from order unity to large values.
A list of all the simulations used in this work, along with some parameters of note, is given in table 1. We limit our study to isotropic 
${Pm} \gtrsim 1$ and modestly large 
${Rm}$; future work should examine whether our conclusions hold in other limits, e.g. 
${Re}\gg {Rm} \gg {Re}_\parallel$ or 
${Rm} \gg {Re} \gg {Re}_\parallel$.
2.3 Diagnostics
Before presenting our results, we define various diagnostics that will be used to study the structure and statistics of the turbulent velocity and magnetic fields. In what follows, 
$\langle \cdot \rangle$ (
$\langle \cdot \rangle _t$) denotes a volume (time) average. The r.m.s. value of a quantity 
$A$ is given by 
$A_\textrm {rms} \doteq \langle A^2 \rangle ^{1/2}$.
2.3.1 Characteristic wavenumbers
The following assortment of characteristic wavenumbers are useful for diagnosing the structure of the magnetic field (see Schekochihin et al. Reference Schekochihin, Cowley, Taylor, Maron and McWilliams2004b):
\begin{gather} k_\parallel \doteq \left( \frac{\left\langle|{\boldsymbol{B}}\,{\boldsymbol{\cdot}}\,{\boldsymbol \nabla}{\boldsymbol{B}}|^2\right\rangle}{\langle B^4\rangle}\right)^{1/2} , \quad k_{{\boldsymbol{B}}{\boldsymbol{\times}}{\boldsymbol{J}}} \doteq \left( \frac{\left\langle|{\boldsymbol{B}}\,{\boldsymbol \times}\,{\boldsymbol{J}}|^2\right\rangle}{\langle B^4\rangle}\right)^{1/2}, \end{gather}
\begin{gather}k_{{\boldsymbol{B}}\,{\boldsymbol{\cdot}}\,{\boldsymbol{J}}} \doteq \left( \frac{\left\langle|{\boldsymbol{B}} \,{\boldsymbol{\cdot}}\,{\boldsymbol{J}}|^2\right\rangle}{\langle B^4\rangle} \right)^{1/2}, \quad k_\textrm{rms} \doteq \left( \frac{\left\langle|{\boldsymbol \nabla}{\boldsymbol{B}}|^2\right\rangle}{\langle B^2\rangle} \right)^{1/2} . \end{gather}
These quantities have simple interpretations: 
$k_\parallel$ measures the variation of the magnetic field along itself, and is typically set by the smallest-scale field-aligned stretching motions; 
$k_{{\boldsymbol {B}}{\boldsymbol {\times }}{\boldsymbol {J}}}$ measures the variation of the magnetic field across itself, and captures the field reversals in folds, which are ultimately limited by resistive dissipation; 
$k_{{\boldsymbol {B}}\,{\boldsymbol {\cdot }}\,{\boldsymbol {J}}}$ measures the variation of the field in the direction orthogonal both to 
${\boldsymbol {B}}$ and to 
${\boldsymbol {B}}\,{\boldsymbol \times }\, {\boldsymbol {J}}$, which tends to be the direction of the greatest compression (Zel'dovich et al. Reference Zel'dovich, Ruzmaikin, Molchanov and Sokolov1984; Schekochihin et al. Reference Schekochihin, Cowley, Taylor, Maron and McWilliams2004b); and 
$k_\textrm {rms}$ provides a general measure of the overall variation of the magnetic field, but tends to follow 
$k_{{\boldsymbol {B}}{\boldsymbol {\times }}{\boldsymbol {J}}}$. For magnetic fields that are arranged in folded sheets – a typical realization during the kinematic stage of the 
${Pm}\gg {1}$ MHD dynamo – the relative ordering of these wavenumbers is 
$k_\parallel \lesssim k_{{\boldsymbol {B}}\,{\boldsymbol {\cdot }}\, {\boldsymbol {J}}} \ll k_{{\boldsymbol {B}}{\boldsymbol {\times }}{\boldsymbol {J}}} \sim k_\textrm {rms} \sim k_\eta$, where 
$k_\eta$ is the spectral cutoff due to resistivity. For a magnetic field arranged in folded ribbons – a typical realization during the saturated state of the 
${Pm}\gg {1}$ MHD dynamo – 
$k_\parallel \ll k_{{\boldsymbol {B}}\,{\boldsymbol {\cdot }}\, {\boldsymbol {J}}} \lesssim k_{{\boldsymbol {B}}{\boldsymbol {\times }}{\boldsymbol {J}}}\approx k_\textrm {rms} \sim k_\eta$.
2.3.2 Transfer functions
We define the shell-filtered kinetic-energy transfer function 
$\mathcal {T}_k[{\boldsymbol {A}}]$ of an arbitrary vector field 
${\boldsymbol {A}}$ via
\begin{equation} \mathcal{T}_k[{\boldsymbol{A}}] \doteq \frac{2}{u_\mathrm{rms}^{2}} \sum_{q \in (2^{-1/4}k, 2^{1/4}k]} {\boldsymbol{u}}_{\boldsymbol{q}}^* \,{\boldsymbol{\cdot}}\, {\boldsymbol{A}}_{\boldsymbol{q}}, \end{equation}
where the star denotes the complex conjugate and 
${\boldsymbol {A}}_{\boldsymbol {q}}$ denotes the Fourier amplitude of 
${\boldsymbol {A}}$ with wavenumber 
${\boldsymbol {q}}$. If 
${\boldsymbol {A}}$ has the units of length squared over time (such as any term from the momentum equation (2.1a)), then 
$\mathcal {T}_k[{\boldsymbol {A}}]$ has the units of inverse time, and 
$\mathcal {T}_k[{\boldsymbol {A}}]$ represents the rate of the kinetic energy flowing due to 
${\boldsymbol {A}}$ into the Fourier shell at 
$k$ of width 
$2^{1/2}k$. For example, 
$\mathcal {T}_k[{\boldsymbol {B}}\,{\boldsymbol {\cdot }}\, {\boldsymbol \nabla } {\boldsymbol {B}}]$ denotes the energy flowing into the velocity field at shell 
$k$ due to the Lorentz force. This diagnostic can be used to probe the scale-by-scale energy balance in (2.1). We also define the r.m.s. shell-filtered kinetic-energy transfer function,
\begin{equation} \mathcal{T}^{\mathrm{rms}}_k[{\boldsymbol{A}}] \doteq \frac{2}{\langle u^4\rangle^{1/2}} \left[ \sum_{q \in (2^{-1/4}k, 2^{1/4}k]}\left(\frac{1}{2}{\rm Re}( {\boldsymbol{u}}_{\boldsymbol{q}}^* \,{\boldsymbol{\cdot}}\, {\boldsymbol{A}}_{\boldsymbol{q}})\right)^2\right]^{1/2}. \end{equation}This diagnostic serves as an alternative to the fourth-order spectra previously used by Schekochihin et al. (Reference Schekochihin, Cowley, Taylor, Maron and McWilliams2004b) and features the added benefit of enabling quantitative comparison between the nonlinear terms in the momentum equation.
Finally, we define a shell-filtering procedure on a vector field 
${\boldsymbol {A}}$ as
\begin{equation} {\boldsymbol{A}}^{[\mathrm{range}]} \doteq \int_{k \in \mathrm{range}} \frac{\textrm{d}^3 {\boldsymbol{k}}}{(2{\rm \pi})^3} \,{\boldsymbol{A}}_{\boldsymbol{k}} \,\textrm{e}^{\textrm{i} {\boldsymbol{k\,{\boldsymbol{\cdot}}\, x}}}, \end{equation}
where the integration is taken over the specified range in 
$k$ space. We utilize three ranges: 
$[k]$ denotes modes in the shell of width 
$\sqrt {2}$ with range 
$(2^{-1/4}k,2^{1/4}k]$; 
$[{<}k]$ denotes all modes with wavenumber magnitude less than 
$2^{-1/4}k$; and 
$[{>}k]$ denotes all modes with wavenumber magnitude greater than 
$2^{1/4}k$. The shell-filtered quantity 
${\boldsymbol {A}}^{[\mathrm {range}]}$ can be used to determine the amount of energy transfer from one region of 
$k$ space to another. For instance, the quantity 
$\mathcal {T}_k[{\boldsymbol {u}}\,{\boldsymbol {\cdot }}\, {\boldsymbol \nabla } {\boldsymbol {u}}^{[< K]}]$ denotes the net transfer of kinetic energy from all modes with wavenumber magnitudes 
${<}2^{1/4}K$ to modes in the shell 
$k \in (2^{-1/4}K,2^{1/4}K]$. Such shell-filtered quantities have been used in analyses of spectral energy transfer in MHD guide-field turbulence (e.g. Alexakis, Mininni & Pouquet Reference Alexakis, Mininni and Pouquet2005; Grete et al. Reference Grete, O'Shea, Beckwith, Schmidt and Christlieb2017) and have also been used alongside Hölder's inequality to establish constraints on non-local transport in the fluctuation dynamo (Beresnyak Reference Beresnyak2012). In this paper, we use the shell-filtered quantity 
${\boldsymbol {u}}^{[\mathrm {range}]}$ to compare the relative strengths of the hydrodynamic nonlinearity and the viscous stresses (figures 10 and 22), as well as to determine how the motions at scale 
$k$ affect the growth of the magnetic energy (figure 24).
2.3.3 Structure functions
To probe the structure of the turbulent velocity field, it is useful to introduce structure functions, which relay information about the scale-by-scale structure and spatial anisotropy with respect to the local magnetic-field direction. In particular, we employ three-point, second-order structure functions for increment 
${\boldsymbol {\ell }}$, defined by
\begin{equation} \mathrm{SF}_2[{\boldsymbol{u}}]({\boldsymbol{\ell}}) \doteq \langle | {\boldsymbol{u}}({\boldsymbol{x}}+{\boldsymbol{\ell}}) - 2{\boldsymbol{u}}({\boldsymbol{x}}) + {\boldsymbol{u}}({\boldsymbol{x}}-{\boldsymbol{\ell}}) |^2\rangle. \end{equation}
These can be used to extract information about variations of a given field along and across the local magnetic-field direction by conditioning the box average on the alignment of the point separation vector 
${\boldsymbol {\ell }}$ with the local magnetic field, defined by 
${\boldsymbol {B}}_{\boldsymbol {\ell }} \doteq [{\boldsymbol {B}}({\boldsymbol {x}}+{\boldsymbol {\ell }})+{\boldsymbol {B}}({\boldsymbol {x}})+{\boldsymbol {B}}({\boldsymbol {x}}-{\boldsymbol {\ell }})]/3$. This conditioning yields the parallel and perpendicular structure functions (cf. Chen et al. Reference Chen, Mallet, Schekochihin, Horbury, Wicks and Bale2012)
\begin{align} \mathrm{SF}_2[{\boldsymbol{u}}](\ell_\parallel) &\doteq \langle | {\boldsymbol{u}}({\boldsymbol{x}}+{\boldsymbol{\ell}}) - 2{\boldsymbol{u}}({\boldsymbol{x}}) + {\boldsymbol{u}}({\boldsymbol{x}}-{\boldsymbol{\ell}}) |^2 ; \:0 \le \theta < {\rm \pi}/18 \rangle, \end{align}
\begin{align}\mathrm{SF}_2[{\boldsymbol{u}}](\ell_\perp) &\doteq \langle | {\boldsymbol{u}}({\boldsymbol{x}}+{\boldsymbol{\ell}}) - 2{\boldsymbol{u}}({\boldsymbol{x}}) + {\boldsymbol{u}}({\boldsymbol{x}}-{\boldsymbol{\ell}}) |^2 ;\: 8{\rm \pi}/18 < \theta \le {\rm \pi}/2\rangle, \end{align}
respectively, where 
$\theta \doteq \arccos | {\boldsymbol {B}}_{\boldsymbol {\ell }}\,{\boldsymbol {\cdot }}\,{\boldsymbol {\ell }}/B_\ell \ell |$ is the angle between the point separation vector and the local magnetic field.
The parallel and perpendicular structure functions may be combined to calculate the scale-dependent anisotropy of the fluctuations. For example, equating the two,
\begin{equation} \mathrm{SF}_2[{\boldsymbol{u}}](\ell_\parallel) = \mathrm{SF}_2[{\boldsymbol{u}}](\ell_\perp), \end{equation}
provides 
$\ell _{\parallel ,{\boldsymbol {u}}}$ as a function 
$\ell _{\perp ,{\boldsymbol {u}}}$, or vice versa.
Power laws with exponent 
$\alpha$ that appear in structure functions translate to Fourier spectra with spectral index 
$-\alpha -1$. The use of a three-point stencil allows one to resolve spectral indices that are less steep than 
$-5$, or structure functions with exponents less than 
$4$ (Lazarian & Pogosyan Reference Lazarian and Pogosyan2008).Footnote 8
3 Results
3.1 Stages of fluctuation dynamo
There are four stages in the typical evolution of the fluctuation dynamo:
(i) The diffusion-free stage, during which the magnetic diffusion due to resistivity has yet to become large enough to influence the growth of the magnetic field. This stage occurs only if the Prandtl number
${Pm}$ is sufficiently large and the scale of the initial magnetic field is much larger than the resistive scale. Once the magnetic field has become folded enough by the flow that the bulk of the magnetic energy reaches the resistive scale, magnetic diffusion becomes important and the dynamo enters…;(ii) …the kinematic stage, in which the magnetic energy continues to grow despite the resistivity (if the flow is indeed a dynamo and
${Rm}$ is above threshold) but the Lorentz force (${\boldsymbol {B}}\,{\boldsymbol {\cdot }}\, {\boldsymbol \nabla } {\boldsymbol {B}}$) is still too feeble to exert any dynamical feedback on the field-amplifying turbulence. In MHD, the kinematic dynamo is linear in ${\boldsymbol {B}}$ (though nonlinear in the random fields), resulting in exponential growth of $B_\mathrm {rms}$ (Kazantsev Reference Kazantsev1968; Moffatt Reference Moffatt1978; Kulsrud & Anderson Reference Kulsrud and Anderson1992). The Braginskii viscosity introduces a dependence in the velocity equation on the unit vector ${\boldsymbol {\hat {b}}}$. If the parallel viscous stress (which is just the Maxwell stress with $B^2/4{\rm \pi}$ replaced by ${\rm \Delta} p$) is sufficiently large, the ‘kinematic’ phase is then fundamentally nonlinear: even though the magnetic field is dynamically weak, its structure influences the properties of the flow, which in turn affects induction in a nonlinear way. On the other hand, if the Braginskii viscosity is subject to hard-wall limiters, then its efficacy is reduced to be comparable with that of the Maxwell stress, rendering this stage truly ‘kinematic’. Eventually, the magnetic field becomes strong enough for the Lorentz force to exert a back reaction on the smallest-scale eddies, suppressing their ability to amplify the magnetic field. The dynamo then enters…;(iii) …the nonlinear stage, in which the magnetic-field amplification is driven by progressively larger (and slower) eddies and the dynamo begins to slow down (e.g. Schekochihin et al. Reference Schekochihin, Cowley, Hammett, Maron and McWilliams2002a, Reference Schekochihin, Cowley, Taylor, Maron and McWilliams2004b; Maron et al. Reference Maron, Cowley and McWilliams2004; Cho et al. Reference Cho, Vishniac, Beresnyak, Lazarian and Ryu2009; Beresnyak Reference Beresnyak2012). An estimate of the growth rate in this stage may be obtained by positing that the Lorentz force disables all eddies with energy less than the total magnetic energy. If one assumes Kolmogorov scalings with eddies at scale
$\ell$ having energy $E(\ell ) \propto \ell ^{2/3}$, then the scale of the smallest eddies that appreciably stretch and grow the field is given by $E(\ell ) \sim \langle B^2\rangle$. This leads to secular evolution:
(3.1)Therefore, as the field-stretching scales shift, exponential growth gives way to a linear-in-time growth of magnetic energy; and\begin{equation} \frac{1}{2}\frac{\textrm{d} \langle B^2\rangle}{\textrm{d} t} = \langle B^2{\boldsymbol{\hat{b}}}{\boldsymbol{\hat{b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla} {\boldsymbol u} \rangle \sim \frac{\langle B^2 \rangle [E(\ell)]^{1/2}}{\ell} \sim \frac{[E(\ell)]^{3/2}}{\ell} \sim \mathrm{const}. \end{equation}(iv) the fourth and final stage of the dynamo is saturation, which is achieved when the magnetic and kinetic energies become comparable (though not necessarily scale by scale – see, e.g. Schekochihin et al. Reference Schekochihin, Cowley, Hammett, Maron and McWilliams2002a, Reference Schekochihin, Cowley, Taylor, Maron and McWilliams2004b).
In what follows, the evolution and characteristics of the dynamo in each of these stages are examined using results from the hard-wall-limited Braginskii-MHD, unlimited Braginskii-MHD, and isotropic-MHD simulations. We begin with a comparison of the limited Braginskii-MHD and 
${Pm} = 1$, isotropic-MHD dynamo, which we find to be similar to one another in almost every respect.
3.2 Limited Braginskii-MHD dynamo is similar to ${Re}\gg {1}$, ${Pm}\gtrsim {1}$ MHD dynamo
The first two stages of the dynamo take place while the magnetic field is dynamically weak. As a result, unless the Braginskii viscosity is negligibly small, a majority (by volume) of the plasma will have pressure anisotropies that exceed the firehose and mirror-instability thresholds, viz. 
$\mu _{{{\rm B}}}|{\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}| \gtrsim B^2$. Applying the hard-wall limiters then effectively disables the Braginskii viscosity in most of the plasma volume, effectively rendering viscous transport mostly isotropic, at least until the saturated state is reached and the magnetic field becomes dynamically strong. In what follows, we demonstrate this point both qualitatively and quantitatively through a series of diagnostics.
3.2.1 Visual appearance of the velocity and magnetic fields
We first demonstrate that most of the qualitative features of the Braginskii-MHD dynamo with hard-wall pressure-anisotropy limiters are similar to those found in isotropic MHD. Figure 1 displays two-dimensional cross-sections in the 
$x$–
$y$ plane of the magnetic-field strength, the velocity magnitude, and the parallel component of the rate-of-strain tensor from a limited Braginskii-MHD simulation (run L1) and from a comparable 
${Re}\gg {1}$, 
${Pm}=1$ simulation using isotropic MHD (run MHD4), both in the kinematic stage and in saturation. The only difference between these simulations is that 
$\mu ^{-1}_{{\rm B}} = 20$ in the former. Despite this difference, the cross-sections of all displayed quantities are difficult to distinguish between the two systems. In both runs, the magnetic field is dominated by small-scale fluctuations that grow to somewhat larger scales in the saturated state. The usual folded structure of the magnetic field, with direction reversals at the resistive scale and field lines curved at the scale of the flow (e.g. Schekochihin et al. Reference Schekochihin, Cowley, Taylor, Maron and McWilliams2004b), is manifest in both the kinematic and nonlinear stages. Both cases also feature a 
${Re}\gg {1}$ flow characterized by chaotic structures across multiple scales, with a tendency for the flow to shift to somewhat larger scales in the saturated state when the now dynamically important Lorentz force is able to exert an influence on the dynamics. In principle, the Braginskii viscous stress could also influence the flow structure and dynamics, but its regulation to values comparable to the Maxwell stress by the hard-wall limiters merely changes the Maxwell stress by a factor of order unity.
3.2.2 Evolution of magnetic energy
Figure 2 displays the evolution of the boxed-averaged magnetic energy 
$\langle B^2 \rangle /2$ for a series of Braginskii-MHD simulations with hard-wall limiters and isotropic-MHD simulations. In the left panel, various Braginskii viscosities 
$\mu _{{\rm B}}$ and fixed Laplacian diffusivities (
$\mu =\eta = 5\times 10^{-4}$) are used. When the magnetic field is very weak, the hard-wall limiters affect the majority of the plasma, thus effectively disabling the parallel viscous stress. As a result, the growth rate of the magnetic energy in the kinematic stage is largely independent of 
$\mu _{{\rm B}}$ and exhibits magnetic-field growth closely resembling its isotropic MHD counterpart (the orange line). In the right panel, the Braginskii viscosity is fixed at 
$\mu ^{-1}_{{\rm B}}=20$ while the isotropic viscosity 
$\mu$ is varied and the magnetic diffusivity is held fixed at 
$\eta ^{-1}=1500$. As in the left panel, the parallel viscosity is effectively disabled at 
$\langle B^2\rangle \ll 1$ by the hard-wall limiters, and so the scale of the fastest eddies is determined instead by the isotropic viscosity. Accordingly, as 
${Pm} \sim 1$ is approached the growth rate of the magnetic energy converges to a maximal value (the growth rate would begin to decrease once 
${Pm}$ decreased to be 
${<}1$; see Vincenzi Reference Vincenzi2002 and Iskakov et al. Reference Iskakov, Schekochihin, Cowley, McWilliams and Proctor2007).
3.2.3 Power spectra of the velocity and magnetic fields
The similarity between MHD and limited Braginskii-MHD is further illustrated in figure 3, which displays kinetic- and magnetic-energy spectra for the hard-wall-limited Braginskii simulation L1 and for the 
${Pm} = 1$ isotropic-MHD simulation MHD4, both in the kinematic stage (purple lines) and in saturation (green lines). For reference, both the Kasantsev 
$k^{3/2}$ scaling and a 
$k$ scaling are shown in the magnetic-energy spectra.Footnote 9 The kinetic-energy spectra from both runs have all the characteristics of a high-
$Re$ turbulent dynamo: a Kolmogorov 
$-5/3$ spectrum in the kinematic stage and a steeper 
$-7/3$ spectrum in the saturated state (as seen at these moderate Reynolds numbers but not explained by Schekochihin et al. Reference Schekochihin, Cowley, Taylor, Maron and McWilliams2004b). The magnetic spectrum in the kinematic stage is consistent with the Kazantsev (Reference Kazantsev1968) 
$k^{3/2}$ scaling at small wavenumbers, while in saturation it is shallower. The peaks of the magnetic spectra for both simulations move to larger scales as saturation is reached. In the isotropic-MHD dynamo, this migration to larger scales has been explained as a consequence of ‘selective decay’ (Schekochihin et al. Reference Schekochihin, Cowley, Hammett, Maron and McWilliams2002a, Reference Schekochihin, Cowley, Taylor, Maron and McWilliams2004b) – the increased importance of resistive dissipation on smaller-scale magnetic-field fluctuations as the Lorentz force begins to suppress field-stretching motions. Because limited Braginskii viscosity in regions of positive parallel rate of strain satisfying 
$3\mu _{{\rm B}}{\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u} \gtrsim 1/\beta$ effectively enhances the magnetic tension by only a factor of 3/2 (since 
$B^2/4{\rm \pi} \rightarrow B^2/4{\rm \pi} + {\rm \Delta} p = (3/2) B^2/4{\rm \pi}$ at the mirror threshold), the similarities between the spectra in the saturated state of the limited Braginskii and isotropic-MHD runs are not particularly surprising. The pressure anisotropy in regions of negative parallel rate of strain (
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u} \lesssim 0$) is instead adiabatically pushed towards the firehose threshold, at which the limited Braginskii stress exactly nullifies the magnetic tension. The effect of this cancellation on the magnetic spectrum in saturation appears to be minimal, however, likely because the volume-filling factor of regions whose pressure anisotropy lies beyond the firehose threshold is small – see figure 25(a). (The paucity of regions with 
$3\mu _{{\rm B}}{\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u} \lesssim -2/\beta$ does not necessarily imply that the magnetic energy is growing in the majority of the volume. Rather, resistive dissipation competes with suppressed field-stretching motions to give saturation, despite the majority of the volume having 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u} > 0$.)
Figure 3. Time-averaged wavenumber spectra of kinetic energy (a) and magnetic energy (b) from run L1 with hard-wall-limited pressure anisotropy (solid lines) and run MHD4 with 
${Pm}=1$ (dotted lines) in the exponential phase (purple lines) and in saturation (green lines).
3.2.4 Scale-by-scale energy transfer
The relative importance of the various terms of the momentum equation (2.1a) in shaping the kinetic-energy spectrum may be calculated using the energy-transfer functions introduced in § 2.3. Figure 4 shows the time-averaged kinetic-energy transfer function 
$\mathcal {T}_k$ (see (2.6)) for isotropic MHD with 
${Pm} \gg 1$ (run MHD1; panel a) and 
${Pm} = 1$ (run MHD4; panel b), as well as for Braginskii MHD employing hard-wall limiters (run L1; panel c). These plots quantify the flow of kinetic energy into and out of the Fourier shell centred geometrically around 
$k$, with solid (dotted) lines signifying energy flowing into (out of) each shell.
Figure 4. Shell-filtered kinetic-energy transfer functions 
$\mathcal {T}_k$ (see (2.6)) for MHD (a,b) and hard-wall-limited Braginskii MHD (c) in the kinematic (left) and saturated (right) stages. Solid (dotted) lines denote energy flowing into (out of) the shell centred at wavenumber 
$k$.
During the kinematic regime of the 
${Pm} = 1$ isotropic-MHD dynamo, the magnetic tension (purple lines) remains subdominant to the other terms in the momentum equation, and thus the energy balance is predominantly between the hydrodynamic nonlinearity (green lines) and energy injection at large scales (orange lines), and the hydrodynamic nonlinearity and viscous dissipation (blue lines) at small scales. In the case of 
${Pm} \gg 1$ (run MHD1), the kinetic-energy cascade is terminated by viscosity at large scales, and so the magnetic tension is balanced by viscous dissipation. This results in a steep sub-viscous velocity spectrum (previously explained by Schekochihin et al. (Reference Schekochihin, Cowley, Taylor, Maron and McWilliams2004b), see figure 18 in that reference). These motions play no dynamical role in the evolution of the turbulence or magnetic field.
As the isotropic-MHD dynamo approaches saturation, the magnetic field becomes dynamically important and its tension begins to influence the turbulent motions. Interestingly, the magnetic tension serves to redistribute turbulent energy from large to small scales, a feature not previously noted in the literature. It is shown in § 3.3 that the Braginskii-MHD dynamo in the unlimited regime also exhibits this behaviour throughout its entire evolution. We demonstrate in appendix B.2 that the small-scale turbulent motions driven by the magnetic tension (or Braginskii stress) serve to nullify the parallel rate of strain 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$, thereby counteracting the overall growth of the magnetic field and allowing the dynamo to saturate.
The hard-wall-limited Braginskii-MHD dynamo exhibits features that are largely similar to the 
${Pm} = 1$ isotropic-MHD dynamo. This is not surprising, as the limiters reduce the strength of the parallel viscous stress to be comparable to the magnetic tension. Importantly, the parallel viscous stress (red lines) mimics the magnetic tension in both the ‘kinematic’ and saturated regimes, whereby it serves to redistribute kinetic energy from large to small scales (see appendix B.2).
3.2.5 Structure functions of the velocity field
The energy spectra shown in figure 3 do not provide information about the relative contributions of field-parallel/perpendicular gradients of field-parallel/perpendicular quantities to the energetics. Because it is only the field-aligned component of the rate-of-strain tensor that is responsible for increasing the magnetic energy (
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$), it is worthwhile to separate out the field-parallel and -perpendicular components of the velocity field and study their individual contributions to the kinetic-energy spectra. To do so, we calculate the structure functions (2.9) and (2.10) of the velocity field in the hard-wall-limited, 
$\mu ^{-1}_{{\rm B}}=20$ Braginskii-MHD simulation (L1H), as well as for the MHD simulation employing hyper-diffusion (MHDH). These simulations use hyper-diffusion in order to maximize the inertial range; this facilitates cleaner measurements of scale-dependent anisotropy than those obtained using only Laplacian dissipation.Footnote 10Figure 5(a) displays the resulting curves in the kinematic stage (left) and in saturation (right). The structure functions for both systems during their kinematic stage are largely isotropic, exhibiting a 
$\ell ^{2/3}$ power law (corresponding to a spectral index of 
$-5/3$) in the inertial range. Beyond the viscous cutoff, all structure functions steepen to a slope close to 
$\textrm {SF}_2 \propto \ell ^4$, which is the maximal slope measurable by structure functions using a three-point stencil.Footnote 11 In the saturated state, both runs exhibit an anisotropization of the turbulence, with the parallel and perpendicular structure functions exhibiting different scalings. The perpendicular structure functions are steeper than Kolmogorov, being roughly proportional to 
$\ell ^{4/3}_\perp$ (corresponding to a 
$-7/3$ spectral index). The 
$\ell _\perp ^{4/3}$ scaling was previously observed by Yousef, Rincon & Schekochihin (Reference Yousef, Rincon and Schekochihin2007), who studied the effects of disparate-scale interactions between turbulence and dynamo-generated magnetic fields on the exact scaling laws of structure functions typically found in isotropic MHD turbulence. We find the slopes of the parallel structure functions to be even steeper, with a scaling of approximately 
$\ell ^{5/3}_\parallel$ (corresponding to a 
$-8/3$ spectral index). The structure functions of the limited-Braginskii and isotropic-MHD systems differ only in the saturated state, with those from run L1H being slightly steeper than those from run MHDH. This is likely due to addition of the parallel viscous stress to the magnetic tension force, resulting in a stronger influence of the magnetic field on the flow.
Figure 5. Various three-point, second-order structure functions for the Braginskii-MHD system with hard-wall-limited pressure anisotropy (L1H) as well as the MHD simulation employing hyper dissipation (MHDH) in the kinematic (left) and saturation (right) stages. Solid (dotted) line denotes the Braginskii-MHD (MHD) simulation. ‘
$\parallel$’ (‘
$\perp$’) refers to the direction parallel (perpendicular) to the local scale-dependent magnetic field (see (2.10)). Panel (b) shows the scale-dependent anisotropy scaling of the parallel coherence length 
$\ell _\parallel$, as defined by (2.11).
The difference in perpendicular and parallel scalings implies a scale-dependent anisotropy in the saturated state of the dynamo, which we quantify using (2.11) and display in the right side figure 5(b). Both runs exhibit a scaling close to 
$\ell _\parallel \sim \ell _\perp ^{3/4}$.Footnote 12 By contrast, in the kinematic stage (left side of figure 5b) the anisotropy scaling is linear in both systems, indicating isotropic turbulence.
3.2.6 Characteristic scales of the magnetic field
As a final point of contact with results from isotropic MHD, we present in figure 6 the characteristic wavenumbers (2.5) quantifying the geometry of the magnetic field in the runs L1 and MHD4. The similarity between the two systems is quite good, even in the saturated state, in which the limited Braginskii viscosity is dynamically important. As the magnetic field is stretched and folded by the flow, it is organized into long thin structures (folds). As a result, the wavenumbers in the kinematic phase satisfy the ordering 
$k_\parallel < k_{{\boldsymbol {B}}\boldsymbol {\cdot } \boldsymbol {J}} < k_{{\boldsymbol {B}}\boldsymbol {\times } {\boldsymbol {J}}}$, with each of these wavenumbers decreasing and the latter two becoming more comparable in the saturated state, as the magnetic folds become more ribbon-like.Footnote 13 The probability distribution function (p.d.f.) of 
$K$ (not shown) is nearly identical to that found in run MHD4 (the blue line in figure 13 below), having a peak concentrated near the viscous scale and a power-law tail 
${\propto }K^{-13/7}$ (as predicted for the high-Pm MHD dynamo by Schekochihin et al. Reference Schekochihin, Cowley, Maron and Malyshkin2002b). The shallowness of these tails implies that all non-zero moments of the curvature p.d.f. diverge in the absence of resistive regularization. As a result, the r.m.s. value of the magnetic-field-line curvature 
$K\doteq |{\boldsymbol {\hat {b}}}\,{\boldsymbol {\cdot }}\,{\boldsymbol \nabla }{\boldsymbol {\hat {b}}}|$ is comparable to 
$k_\textrm {rms} \sim k_\eta$ (Malyshkin Reference Malyshkin2001). The curvature p.d.f. is then representative of a three-dimensional field whose regions of large curvature occupy only a small fraction of the volume.
Figure 6. Evolution of the characteristic wavenumbers (2.5) for isotropic MHD (dashed) and Braginskii MHD with hard-wall pressure-anisotropy limiters (solid lines).
3.3 Unlimited Braginskii-MHD dynamo is similar to saturated MHD dynamo
Having provided several lines of evidence that pressure-anisotropy limiters revert the Braginskii-MHD dynamo to its more traditional 
${Re}\gg {1}$, 
${Pm}\gtrsim {1}$ counterpart, and motivated by the observation of imperfect regulation of pressure anisotropy in the hybrid-kinetic simulations presented in Paper I, we now turn off those limiters and let the full Braginskii viscous stress act unabated on the flow. In this case, the dynamo takes on a very different character. Left unchecked by limiters, and without a rapidly growing mirror instability in Braginskii MHD to reign it in (see appendix A), the pressure anisotropy is proportional to the parallel rate of strain and to 
$\mu _\textrm {B}$. For large values of 
$\mu _{{\rm B}}$, it spills over both the firehose and mirror thresholds. No longer bound to the relatively meagre Lorentz force in the kinematic stage, the large parallel viscous stress exerts a strong dynamical feedback on those fluid motions responsible for amplifying the magnetic field. The result is a strong suppression of 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$, leading to a viscous anisotropization of the fluid flow. Because (2.1b) implies 
$\textrm {d}\ln B/\textrm {d}t = {\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$ in the absence of resistivity, this anisotropization ultimately results in a less efficient dynamo. Since the Braginskii viscous stress is similar in form to the magnetic tension force 
${\boldsymbol {B}}\,{\boldsymbol {\cdot }}\, {\boldsymbol \nabla } {\boldsymbol {B}} = {\boldsymbol \nabla } \,{\boldsymbol {\cdot }}\, (B^2 {\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}})$, with the pressure anisotropy playing the role of the magnetic-field strength (viz. 
$B^2 \rightarrow {\rm \Delta} p \propto \textrm {d}_t B^2$), one may expect similarities between the unlimited Braginskii-MHD dynamo and the isotropic-MHD dynamo in its saturated state. As in § 3.2, we confirm these expectations by using a variety of diagnostics to compare the unlimited Braginskii-MHD runs with the isotropic-MHD runs. Some notable differences are also highlighted.
3.3.1 Visual appearance of the velocity and magnetic fields
Figure 7 displays the same quantities as in figure 1 but now for a Braginskii-MHD simulation without limiters (U1) and an isotropic-MHD simulation with 
$\mu =\mu _{{\rm B}}/5$ (MHD2). This choice of isotropic viscosity for comparison with Braginskii MHD was advocated by Malyshkin & Kulsrud (Reference Malyshkin and Kulsrud2002) as an effective closure for systems with a magnetic field that is isotropically tangled on sub-viscous scales. Indeed, the velocity fields in these two runs appear to be similar on scales larger than those on which the magnetic field is tangled. However, a major difference is that the velocity field in the unlimited Braginskii run features thin striations with sharp gradients across the local magnetic-field direction. The isotropic-MHD system instead exhibits only large-scale motions typical of order-unity Reynolds-number flow. As a result, the Braginskii turbulent state throughout the evolution of the dynamo resembles more closely the saturated state of the high-
$Re$ MHD system (cf. centre panel of figure 1d), in which large-scale features of the velocity field are accompanied by small-scale turbulence. Another notable difference between the unlimited-Braginskii and isotropic-MHD runs is that, in the unlimited-Braginskii run, the magnetic field, velocity, and parallel rate of strain all exhibit remarkably little change in going from the kinematic stage to the saturated state – a feature to which we return repeatedly throughout this section.
Figure 7. Same as figure 1 but for an unlimited Braginskii-MHD simulation (U1) and a comparable MHD run (MHD2). The left, central and right panels display the magnetic-field strength, velocity magnitude and parallel rate of strain, respectively. Panels (a) and (b) display results from the unlimited simulation U1 (parameters 
$\mu ^{-1}=1500$, 
$\mu _{{\rm B}}^{-1}=20$, 
$\eta ^{-1}=1500$): (a) the kinematic stage; (b) the second row in the saturated state. Panels (c) and (d) display results from the MHD simulation MHD2 (parameters 
$\mu ^{-1}=100$, 
$\eta ^{-1}=1500$): (a) the kinematic stage; (b) the second row in the saturated state. All plots are on a linear scale, with brighter regions denoting higher magnitudes.
A more subtle difference between runs U1 and MHD2 concerns the parallel rate of strain, displayed in the rightmost panels of figure 7. The MHD simulation features larger-scale patches of 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$, with more extreme values, than found in the Braginskii-MHD simulation. This is because, in the unlimited Braginskii system, there is a dynamical feedback whereby the full pressure anisotropy driven by the field-stretching motions (viz. 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}\ne 0$) dynamically suppresses those very same motions. There are also fewer regions that exhibit strong negative values of 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$ in the unlimited run, most likely because the act of decreasing the magnetic-field strength with 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}<0$ is unstable to the production of firehose fluctuations that grow small-scale magnetic fields (and thus contribute a positive 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$; see Schekochihin et al. (Reference Schekochihin, Cowley, Kulsrud, Rosin and Heinemann2008), Rosin et al. (Reference Rosin, Schekochihin, Rincon and Cowley2011), and Melville et al. (Reference Melville, Schekochihin and Kunz2016) for further discussion of this point in the context of the parallel firehose instability). Further discussion of the contributions of small and large scales to the total parallel rate of strain and the suppression of 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$ is deferred to appendix B.
3.3.2 Evolution of magnetic energy
The dynamical feedback of the full Braginskii viscosity on the flow affects the time evolution of the magnetic energy, shown in figure 8. As the Braginskii viscosity is increased, the viscous scale of the field-stretching motions becomes larger and the dynamo growth rate decreases accordingly (left side), indicating that in the unlimited regime the dynamo growth rate is ultimately controlled by the Braginskii viscosity. Indeed, the growth of the average magnetic energy in the 
$\mu ^{-1}_{\rm B} = 20$ and 
$100$ unlimited-Braginskii runs (purple and green solid lines, respectively) is almost identical to that in the 
$\mu ^{-1} = 20$ and 
$100$ isotropic-MHD runs (red and blue dashed lines, respectively, in the same plot). Note that the value 
$\mu =\mu _{{\rm B}}/5$ advocated by Malyshkin & Kulsrud (Reference Malyshkin and Kulsrud2002) leads to a much faster dynamo than the unlimited-Braginskii-MHD dynamo it is intended to model (e.g. compare the blue-dashed and purple-solid lines). The reason for this difference is that, while the two simulations exhibit similar stretching motions, the resistive dissipation is effectively much higher in the unlimited Braginskii system due to the increased fraction of field-mixing motions in the rate-of-strain tensor (relative to the field-stretching motions, which are viscously suppressed). This increased magnetic dissipation results in a slower dynamo (see § 4.3 for more on this subject). This indicates that the details of the sub-parallel-viscous range are important for the overall operation of the dynamo, and that the closure advocated by Malyshkin & Kulsrud (Reference Malyshkin and Kulsrud2002) is inappropriate for the case of large isotropic Reynolds number.
Figure 8. Evolution of the kinetic (dash-dotted lines) and magnetic (solid lines) energies for simulations that do not employ pressure-anisotropy limiters with varying Braginskii viscosity and fixed isotropic viscosity (left) and varying isotropic viscosity and fixed Braginskii viscosity (right). The evolution of magnetic energy in the MHD1 and MHD2 simulations (
$\mu ^{-1} = 20$ and 
$100$, respectively) is included for reference in the left panel (dashed lines). Solid lines in the left (right) panel have 
$\mu ^{-1} = 1500$ (
$\mu _{{\rm B}}^{-1} = 20$).
In the right side of figure 8, the Braginskii viscosity is fixed at 
$\mu ^{-1}_{{\rm B}} = 20$ while the isotropic viscosity is varied. Somewhat counter-intuitively, the growth rate appears to decrease as the isotropic viscosity is decreased. A plausible explanation of this result is that a small isotropic viscosity can allow a cascade of perpendicular (or ‘interchange’-like) motion to small scales. These motions, while not capable of growing the magnetic field, can bring field lines closer together and thereby accelerate the resistive reconnection of the field. The deleterious effect of these mixing motions is a central issue in the 
${Pm} < 1$ dynamo (Vincenzi Reference Vincenzi2002; Boldyrev & Cattaneo Reference Boldyrev and Cattaneo2004; Iskakov et al. Reference Iskakov, Schekochihin, Cowley, McWilliams and Proctor2007; Schekochihin et al. Reference Schekochihin, Iskakov, Cowley, McWilliams, Proctor and Yousef2007), where mixing from all kinetic-energy scales promotes resistive annihilation, while only motions larger than the resistive scale can aid in amplification of the magnetic energy via stretching. This idea is further developed in § 4.3. Another anti-dynamo factor is that the rate of strain of the smaller-scale motions allowed by the decreased isotropic viscosity could act to cancel partially the parallel rate of strain driven by the large scales, modifying the dynamo growth rate without significantly increasing resistive mixing. This effect has been recently seen in Braginskii-MHD simulations of the magnetorotational instability (Kempski et al. Reference Kempski, Quataert, Squire and Kunz2019). This partial cancellation, related to the concept of ‘magneto-immutability’, is investigated further in appendix B. It will be shown that, while both effects are present in our simulations, the latter is of only minor consequence, while the former – that of the efficiency of sub-parallel-viscous mixing motions – has a significant impact on whether or not unlimited Braginskii MHD can exhibit a dynamo.
3.3.3 Power spectra of the velocity and magnetic fields
A theme of the last two sub-sections is that, in Braginskii MHD, the perpendicular component of the fluid motions is allowed to cascade to sub-parallel-viscous scales. Figure 9(a) provides quantitative evidence that this is the case, with less isotropic viscosity allowing a shallower spectrum and thus stronger small-scale velocity fluctuations. At the smallest non-zero value of isotropic viscosity (purple line), the unlimited runs exhibit a kinetic-energy spectrum 
$E(k)\sim k^{-2.75}$, which appears to be asymptotic in 
$\mu \rightarrow 0$ at fixed 
$\mu _{{\rm B}}^{-1} = 20$ (based on preliminary studies at even higher Re). While steep, this power law still implies a local rate of strain that increases with 
$k$. Eventually the spectrum experiences a break at 
$kL/2{\rm \pi} \approx 30$, whereupon it exhibits a 
$k^{-4}$ power law down to the grid. This break can be taken as the effective perpendicular viscous scale, the slope beyond which is sufficiently steep (
$\text {spectral index} <-3$) for the rate of strain to decrease with 
$k$, i.e. for the fastest eddies to occur at the largest scales. The run with the next smallest value of isotropic viscosity (green line) also exhibits a similar spectrum, but the spectral break occurs around 
$kL/2{\rm \pi} \approx 8$. For 
$\mu ^{-1}=96$ (orange line), almost the entire kinetic-energy spectrum is proportional to 
$k^{-4}$.
Figure 9. Time-averaged wavenumber spectra of kinetic energy (a) and magnetic energy (b) from simulations with unlimited pressure anisotropy for various values of isotropic viscosity 
$\mu$ in the kinematic stage (top) and the saturated state (bottom).
Note that a small kink exists at the de-aliased grid-scale wavenumber 
$kL/2{\rm \pi} = 224/3 \approx 75$ in all of the kinetic-energy spectra, regardless of the value of the isotropic viscosity. This is a result of small-scale energy injection by the unregulated mirror and firehose instabilities that are present in unlimited Braginskii MHD, the latter of which exhibits an ultraviolet catastrophe when 
$\mu = \eta = 0$ (see appendix A). In our simulations, this small-scale energy injection is balanced by dissipation via the isotropic viscosity (this is shown later in figure 10).
Figure 10. Shell-filtered kinetic-energy transfer functions 
$\mathcal {T}_k$ (see (2.6)) for the unlimited Braginskii-MHD system in the (left) kinematic and (right) saturated stages. Top (bottom) row corresponds to 
$\mu ^{-1} = 1500$ (
$\mu ^{-1} = 96$). Solid (dotted) lines denote energy flowing into (out of) the shell centred at mode 
$k$.
The kinetic-energy spectra shown in figure 9 appear to be independent of whether the dynamo is in the kinematic stage (top) or in saturation (bottom) – a notable difference from the high-
$Re$ MHD dynamo, in which the kinetic-energy spectrum steepens from the Kolmogorovian 
$k^{-5/3}$ to 
$k^{-7/3}$ in the saturated state (at least at these limited resolutions; cf. Schekochihin et al. Reference Schekochihin, Cowley, Taylor, Maron and McWilliams2004b). The accompanying magnetic-energy spectra, shown in figure 9(b), also show little evolution from the kinematic stage to saturation (consistent with the visualization of the magnetic-field evolution in figure 7). In both the kinematic and saturated stages, these magnetic spectra are shallower than the Kazantsev 
$3/2$ scaling, and are instead closer to the scaling seen in the saturated state of the high-
$Re$ MHD simulation (cf. green dotted line in figure 3b).
3.3.4 Scale-by-scale energy transfer
The shell-filtered kinetic-energy transfer functions are again used here to probe the origin of the kinetic-energy spectra found in the previous subsection and are shown in figure 10 for the unlimited Braginskii-MHD simulations with 
$\mu ^{-1} = 1500$ (top) and 
$\mu ^{-1}= 96$ (bottom). For both systems, the 
$k^{-4}$ spectrum of sub-perpendicular-viscous motions (i.e. motions occurring at scales smaller than the viscous cutoff caused by the isotropic viscosity) is seen to be the result of a balance between the pressure-anisotropy stress (red lines) and the isotropic viscosity (blue lines). Similar behaviour was measured and explained in Schekochihin et al. (Reference Schekochihin, Cowley, Taylor, Maron and McWilliams2004b) by balancing the viscous dissipation 
$\mu \nabla ^2 {\boldsymbol {u}}$ and the magnetic tension 
$({\boldsymbol {B}}\,{\boldsymbol {\cdot }}\, {\boldsymbol \nabla } {\boldsymbol {B}})$, resulting in sub-viscous motions that, while initially small, grew along with the magnetic energy. In our simulations of unlimited Braginskii-MHD, the same scenario applies, but now the pressure-anisotropy stress 
${\boldsymbol \nabla } \,{\boldsymbol {\cdot }}\, ({\boldsymbol {\hat {b}}} {\boldsymbol {\hat {b}}} {\rm \Delta} p)$ takes on the role of the magnetic tension 
${\boldsymbol \nabla } \,{\boldsymbol {\cdot }}\, ({\boldsymbol {\hat {b}}} {\boldsymbol {\hat {b}}} B^2)$. In order for the balance 
${\boldsymbol \nabla } \,{\boldsymbol {\cdot }}\, ({\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}{\rm \Delta} p) \sim \mu \nabla ^2 {\boldsymbol {u}}$ to be valid, it must be the case that 
$|{\boldsymbol {u}}\,{\boldsymbol {\cdot }}\, {\boldsymbol \nabla } {\boldsymbol {u}}| \ll |\mu \nabla ^2 {\boldsymbol {u}}|\sim |{\boldsymbol \nabla } \,{\boldsymbol {\cdot }}\, ({\boldsymbol {\hat {b}}} {\boldsymbol {\hat {b}}} {\rm \Delta} p) |$, inequalities that are satisfied in the unlimited Braginskii-MHD simulation beyond the spectral break. Similar to the role of magnetic tension on small scales in the isotropic MHD dynamo, the Braginskii viscosity in this range gives energy to the velocity field instead of dissipating it, a feature that is discussed further in appendix B.2 in the context of magneto-immutability. Both the sub-viscous motions in the 
${Pm} \gg 1$ MHD dynamo (§ 3.2) and sub-perpendicular-viscous motions in unlimited Braginskii-MHD dynamo are passive and do not influence the evolution of either the turbulent motions above the viscous scale or the magnetic field.
The observed 
$k^{-2.75}$ slope in the kinetic-energy spectrum for the simulation with 
$\mu ^{-1}= 1500$ can be seen as a balance between the Braginskii viscosity (red line) and the hydrodynamic nonlinearity (solid and dotted green lines). In this range of 
$k$, the Braginskii viscosity removes kinetic energy at every scale, resulting in a spectrum that is steeper than 
$-5/3$. Much like the spectral index of 
$-7/3$ reported for the kinetic-energy spectra of saturated MHD dynamo (Schekochihin et al. Reference Schekochihin, Cowley, Taylor, Maron and McWilliams2004b), the exact nature of this index is difficult to pin down, although the difference in the two indices (2.75 vs. 
$7/3\simeq 2.33$) does represent one of the minor distinctions between the unlimited Braginskii-MHD and the saturated MHD dynamo. Because the unlimited Braginskii viscosity renders the turbulence anisotropic with respect to the magnetic field regardless of the strength of the latter, the unlimited Braginskii-MHD dynamo increasingly resembles the MHD dynamo in the saturated state as the parallel viscosity is increased. For a system with a sufficiently large parallel viscosity and small perpendicular viscosity (resulting in strongly anisotropic turbulence), one only needs a small amount of magnetic tension to anisotropize the velocity field further and push the dynamo into its saturated state. This is clearly seen in figure 10, in which the magnetic tension is shown to be sub-dominant to the Braginskii stress for the 
${Re} \gg 1$ unlimited Braginskii system, even in the saturated state. It is thus no surprise then that the spectrum of the underlying fluid motions remains relatively constant throughout the entire dynamo process when the pressure anisotropy is unlimited. The reduced magnetic tension required to saturate the dynamo also implies that the unlimited Braginskii system saturates with a smaller magnetic to kinetic-energy ratio than the MHD dynamo, a feature that can be seen by comparing the purple, green, blue and orange solid lines in figure 8 that show how the evolution of the magnetic energy changes with decreasing 
$\mu _{{\rm B}}$. However, the reduction in this ratio is minor (order unity) unless the parallel-viscous scale implied by 
${Re}_\parallel$, viz. 
$\ell _{\mu _\parallel } \sim \ell _0 \, {Re}_\parallel ^{-3/4}$, resides near the outer scale. In the next two sections, we show that the field-biased structure of the flow, the characteristic scales of the magnetic field, and the statistics of the magnetic-field curvature are all relatively constant as well.
3.3.5 Structure functions of the velocity field
We have shown that an unlimited anisotropic viscous stress strongly biases the properties of the flow with respect to the magnetic-field direction. To quantify this bias further, we calculate the structure functions (2.9) and (2.10) of the velocity field in the unlimited 
$\mu ^{-1}_{{\rm B}}=20$ Braginskii-MHD simulation (run U1) during the kinematic phase and in saturation. The result is shown in figure 11(a). We find 
$\mathrm {SF}_2 \propto \ell _\perp ^{1.75}$ for the perpendicular structure function (corresponding to the 
$-2.75$ slope seen in figure 9) and a much steeper 
$\mathrm {SF}_2 \propto \ell _\parallel ^{2.5}$ for the parallel structure function. The steepness of the parallel structure function confirms that small-scale stretching motions play no dynamical role in this Braginskii-MHD dynamo, and so the largest scales are those primarily responsible for growing the magnetic field. The corresponding spectral anisotropy in both the kinematic and saturated stages, shown in figure 11(b), scales roughly as 
$\ell _\parallel \sim \ell _\perp ^{3/4}$, the same as in the saturated states of the limited Braginskii-MHD and isotropic-MHD dynamos (cf. figure 5). Most notably, none of these properties change from the kinematic stage to saturation.
Figure 11. Perpendicular and parallel three-point, second-order structure functions (a) and corresponding spectral anisotropy (b) from the unlimited Braginskii-MHD run U1 in the kinematic (left) and saturated (right) stages.
3.3.6 Characteristic scales of the magnetic field
The geometry of the magnetic field in the unlimited Braginskii-MHD runs does not change much either as the dynamo saturates. This is diagnosed quantitatively by measuring the characteristic wavenumbers (2.5) of the magnetic field and the distributions of the field-line curvature. Figure 12 shows the evolution of the former from runs U1, U4 and MHD1, 2 and 4. The wavenumbers in both unlimited Braginskii-MHD runs (U1 and U4, denoted by the purple and green lines) are roughly constant in time, with those from run U1 holding values similar to those found in the saturated state of the 
$\mu ^{-1}=1500$ MHD run (blue lines) and those from run U4 matching closely its MHD counterpart (
$\mu ^{-1} = 100$, orange line). This is to be contrasted with the high-
${Re}$ MHD simulation (light blue line) whose characteristic scales of the magnetic field increase as the system saturates, which is due to selective decay of small-scale modes of the magnetic field as a result of the anisotropization of the turbulence by the magnetic tension.
Figure 12. Evolution of the characteristic wavenumbers (2.5) in the unlimited Braginskii-MHD runs U1 and U4 (purple and green lines, respectively) and in runs MHD1, 2, 4 (yellow, orange and blue lines, respectively).
A somewhat surprising result for the Braginskii-MHD system is that 
$k_\parallel$ is larger in run U1 than in the isotropic-MHD simulations with either 
$\mu ^{-1}= \mu _{{\rm B}}^{-1} = 20$ or 
$\mu ^{-1} = \mu _{{\rm B}}^{-1}/5 = 100$. This is likely because anisotropic viscosity allows a cascade of perpendicular energy to the perpendicular viscous scale (which is set by 
$\mu$), thereby allowing small-scale mixing motions that can bring field lines closer together and promote resistive annihilation. This causes characteristically shorter fold lengths. As the isotropic viscosity for the unlimited Braginskii simulations is increased, the saturated value of 
$k_\parallel$ approaches that found in the small-Re MHD runs (as demonstrated by U4).
Figure 12 also suggests that Braginskii viscosity does not promote magnetic fields with larger-scale structure than those produced in isotropic MHD. It has been conjectured by Malyshkin & Kulsrud (Reference Malyshkin and Kulsrud2002) that anisotropic viscosity might cause the turbulent dynamo to inverse cascade the saturated magnetic fields from resistive scales to the larger viscous scales (i.e. those independent of 
${Rm}$) by allowing small-scale perpendicular motions that might unfold the field. If this conjecture turned out to be true, then it might explain the relatively large scale of the observed magnetic fields in the weakly collisional ICM. At least at our modest resolutions, no additional unfolding seems to take place, as the magnetic fields generated in the Braginskii simulations exhibit similar 
$k_{{\boldsymbol {B}}{\boldsymbol {\times }} {\boldsymbol {J}}}$ and 
$k_\textrm {rms}$ to those found in the saturated state of the isotropic-MHD runs. Efforts to extend our work to larger Rm, and thus larger scale separation, could clarify what sets the peak of the magnetic-energy spectrum in the Braginskii-MHD dynamo.
The p.d.f. of the field-line curvature 
$K \doteq |{\boldsymbol {\hat {b}}}\,{\boldsymbol {\cdot }}\, {\boldsymbol \nabla } {\boldsymbol {\hat {b}}}|$ is displayed in figure 13 for the same runs as in figure 12. The asymptotic form of the curvature p.d.f. predicted by Schekochihin et al. (Reference Schekochihin, Cowley, Maron and Malyshkin2002b) (viz., a peak concentrated near the viscous scale and a power-law tail of 
$K^{-13/7}$ – see § 3.2.6) is manifest in all of our simulations. As with the other quantities that we have analysed, the curvature p.d.f.s from the unlimited Braginskii run in both its kinematic and saturated stages match closely those measured in the saturated state of isotropic MHD. The form of the p.d.f. changes somewhat for the 
${Pm} = 1$ case (blue line), in which the peak moves to smaller scales in accordance with the higher Reynolds number. By contrast, the unlimited Braginskii-MHD simulation exhibits slightly more curvature than comparable MHD simulations (runs MHD1 and MHD2; compare the purple line with the yellow and orange lines). This small increase likely occurs for the same reasons that the characteristic wavenumber 
$k_\parallel$ increases in the unlimited Braginskii case (see preceding paragraph): small-scale motions that mix field lines and promote resistive reconnection of the magnetic field are allowed by the Braginskii viscosity (but not by the isotropic viscosity).
Figure 13. Probability distribution function of the magnetic-field curvature 
$K \doteq |{\boldsymbol {\hat {b}}} \,{\boldsymbol {\cdot }}\, {\boldsymbol \nabla } {\boldsymbol {\hat {b}}}|$ in run U1 (purple line) and in runs MHD1, 2 and 4 (yellow, orange and blue lines, respectively) in the kinematic stage and in saturation.
4 Anisotropization of the rate of strain and its consequences
As was motivated at the start of § 3.3 and supported by the accompanying figures, the ‘kinematic’ stage of the unlimited Braginskii-MHD dynamo has many characteristics in common with the saturated state of the 
${Pm}\gtrsim {1}$ isotropic-MHD dynamo. This is because the Braginskii viscous stress has a mathematical form very similar to that of the magnetic tension, which, in saturation, biases the fluid velocity with respect to the magnetic-field direction to reduce the parallel rate of strain 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$. In this section, we explore this bias further. Namely, we focus on how the rate-of-strain tensor remains field-anisotropic throughout the evolution of the Braginskii-MHD dynamo and on the consequences of this anisotropy. We also provide further evidence for, and formulate a theoretical dynamo model that is based on, the observed similarity between the unlimited-Braginskii and saturated-MHD dynamos.
4.1 Constraints on stretching versus mixing
The starting point for this discussion is an assessment of the roles that field-line stretching and mixing play in the efficacy of the dynamo – roles that have already been glimpsed in the preceding sections. These two effects are quantified, respectively, by the parallel rate of strain, 
${\nabla }_\parallel u_\parallel \doteq {\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$, and by the perpendicular variations of the perpendicular velocity, 
${\boldsymbol \nabla }_\perp {\boldsymbol {u}}_\perp \doteq ({\boldsymbol {\mathsf {I}}}-{\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}})\,{\boldsymbol {\cdot }}\, {\boldsymbol \nabla } {\boldsymbol {u}}\,{\boldsymbol {\cdot }}\,({\boldsymbol {\mathsf {I}}}-{\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}})$. The key property of the fluctuation dynamo, exemplified by the Zel'dovich (Reference Zel'dovich1957) anti-dynamo theorem, is that amplification of the magnetic energy by field-line stretching can overcome resistive diffusion only in three dimensions. This is because a folded magnetic field undergoing stretching in a two-dimensional incompressible flow necessarily involves bringing together oppositely directed magnetic fields in the compression direction of the flow (see figure 10 of Schekochihin et al. Reference Schekochihin, Cowley, Taylor, Maron and McWilliams2004b). This causes very rapid annihilation of the field by resistivity. In three dimensions, this dynamo-adverse configuration may be avoided by orienting the directions of any magnetic-field reversals (which are parallel to 
${\boldsymbol {B}}\,{\boldsymbol \times }\, {\boldsymbol {J}}$) to be perpendicular to both the stretching and compression directions of the flow – the ‘null’ direction – thus ensuring that current sheets in the folded field are not subject to further compression and thus to enhanced resistive diffusion. Of course, this special configuration cannot be maintained everywhere in a chaotic flow, and so the triumph of field-line stretching over field-line mixing and enhanced diffusion is a statistical one: in the now-classic Zel'dovich et al. (Reference Zel'dovich, Ruzmaikin, Molchanov and Sokolov1984) model of the kinematic fluctuation dynamo in a random, linear velocity field, only an exponentially decreasing fraction of initial wavenumbers 
${\boldsymbol {k}}(t=0)$ characterizing the magnetic-field geometry – those corresponding to field reversals occurring only in the null direction – contribute to the growing field at any given time.
In the Braginskii-MHD dynamo, the relative sizes of 
$|{\nabla }_\parallel u_\parallel |$ and 
$|{\boldsymbol \nabla }_\perp {\boldsymbol {u}}_\perp |$ are controlled both by the dynamical influence of the magnetic field on the flow through the Lorentz force (in the saturated state) and by the disparate values of the parallel (
$\mu _{{\rm B}}$) and perpendicular (
$\mu$) viscosities (throughout the entire evolution, including in the ‘kinematic’ regime). The latter truncate the various components of the rate-of-strain tensor at disparate scales, thereby forcing an imbalance in the maximal rates of strain along and across the magnetic field and placing constraints on the alignment between the magnetic folds and the stretching direction of the flow. The principal theme of the following sub-sections is that it is this imbalance, controlled by the anisotropic viscosity, that determines whether or not the Braginskii-MHD dynamo is viable. In effect, the magnetic Prandtl number 
$\mu /\eta$, which governs the relative rates of viscous-scale stretching and resistive-scale diffusion in MHD, is replaced as the key parameter in Braginskii MHD by the pair of parameters 
$\mu /\eta$ and 
$\mu /\mu _{{\rm B}}$. The survival of the Braginskii-MHD dynamo depends on both of these ratios, the latter (when small) constraining the geometry of the flow and thus its ability to amplify magnetic fields.
4.2 Alignment of the rate-of-strain tensor with the magnetic-field direction
It is therefore of interest to use our numerical simulations to measure directly the alignment of the rate-of-strain tensor with the magnetic-field direction during the various stages of the dynamo and across a range of parallel and isotropic viscosities.
Figure 14 displays the time-averaged values of 
${\nabla }_\parallel u_\parallel$ (left) and 
${\boldsymbol \nabla }_\perp {\boldsymbol {u}}_\perp$ (right) for various simulations using isotropic MHD (a) and unlimited Braginskii MHD (b and c). For the high-
${Re}$ MHD simulations, the saturated state is characterized by a reduction in the magnitude of the parallel rate of strain 
$|{\nabla }_\parallel u_\parallel |$ (left side of panel a for 
$\mu ^{-1} = 500$, 1500) as the magnetic tension force becomes dynamically important. The mixing motions (figure 14a, right side) are partially suppressed as well. By contrast, the parallel rate of strain in the unlimited Braginskii run with 
$\mu ^{-1}_{{\rm B}}=20$ is nearly constant in time and independent of isotropic viscosity (figure 14b, left side). The perpendicular rate of strain is nearly unchanged in the saturated state as well, but varies with 
$\mu ^{-1}$ in a predictable way: smaller viscosity allows smaller-scale mixing motions, as in isotropic MHD. When the parallel viscosity is decreased at fixed 
$\mu$ (figure 14c), the evolution of both the parallel and perpendicular rates of strain in the Braginskii-MHD system (rightmost points in figure 14c) converges to those seen in the analogous isotropic MHD system (rightmost points in figure 14a), as expected. Interestingly, the overall levels of the perpendicular mixing motions in the unlimited-Braginskii MHD simulations (figure 14b, right side) are comparable to those found in the isotropic-MHD systems in the saturated state (figure 14a, right side), indicating that these motions are partially suppressed relative to those in the kinematic stage of analogous isotropic MHD systems, even though the parallel viscosity does not affect them directly. It will be shown in § 4.3 that, in order for the dynamo to be viable, the relative size of the perpendicular (mixing) motions cannot greatly exceed that of the parallel (stretching) ones.
Figure 14. Time- and box-averaged values of the parallel (
${\nabla }_\parallel u_\parallel \doteq \boldsymbol {\hat {b}}\boldsymbol {\hat {b}}\boldsymbol{:}\boldsymbol {\nabla } \boldsymbol {u}$, left) and perpendicular [
$\boldsymbol {\nabla }_\perp {\boldsymbol {u}}_\perp \doteq ( \boldsymbol {\mathsf {I}} - {\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}) \,{\boldsymbol {\cdot }}\, \boldsymbol {\nabla } {\boldsymbol {u}} \,{\boldsymbol {\cdot }}\, ( \boldsymbol {\mathsf {I}} - {\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}})$, right] rates of strain for MHD (a) and unlimited Braginskii MHD with varying isotropic viscosity 
$\mu$ (b) and varying parallel viscosity 
$\mu _{{\rm B}}$ (c). Purple pluses (green crosses) are time-averages taken in the exponential growth phase (saturated state). All simulations have 
$\eta ^{-1}= 1500$.
An alternative way to quantify the anisotropization of the velocity caused by the Braginskii viscosity is to examine the relationship between the magnetic-field direction 
${\boldsymbol {\hat {b}}}$ and the eigenvectors of the symmetrized rate-of-strain tensor 
${\boldsymbol {\mathsf {S}}} \doteq ({\boldsymbol \nabla } {\boldsymbol {u}} + {\boldsymbol \nabla } {\boldsymbol {u}}^{T})/2$. As we are considering only incompressible turbulence, the trace of the tensor 
${\boldsymbol {\mathsf {S}}}$, and thus the sum of its eigenvalues 
$\lambda _i$, are zero. We order these eigenvalues so that 
$\lambda _1 > \lambda _2 > \lambda _3$. Because 
${\boldsymbol {\mathsf {S}}}$ is real and symmetric, its eigenvalues are real, and its eigenvectors 
$\hat {\boldsymbol {e}}_i$ are orthogonal, hence 
$|\boldsymbol {\hat {e}}_1\,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}|^2 + |\boldsymbol {\hat {e}}_2\,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}|^2 + |\boldsymbol {\hat {e}}_3\,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}|^2 = 1$. The eigenvectors 
$\boldsymbol {\hat {e}}_1$ and 
$\boldsymbol {\hat {e}}_3$ correspond to the directions of field-line stretching and compression, while 
$\boldsymbol {\hat {e}}_2$ points in the so-called ‘null’ direction (which can be either stretching or compressing). The alignment angles 
$\theta _i$ are defined by 
$|\boldsymbol {\hat {e}}_i \,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}| = \cos \theta _i$. For isotropic hydrodynamic turbulence, 
$\lambda _1 \sim -\lambda _3 \sim (3\text {--}5) \lambda _2 > 0$ (Ashurst et al. Reference Ashurst, Kerstein, Kerr and Gibson1987; Tsinober, Kit & Dracos Reference Tsinober, Kit and Dracos1992).
Figure 15 displays the p.d.f.s of the alignment cosines 
$|\boldsymbol {\hat {e}}_i \,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}|$ for 
${Re}\gg {1}$, 
${Pm} =1$ MHD (run MHD4; top) and unlimited Braginskii MHD (run U1; bottom) in the kinematic stage (solid lines) and saturated state (dashed lines).Footnote 14 For the MHD run in the kinematic stage, the most probable values of both 
$\boldsymbol {\hat {e}}_1 \,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}$ and 
$\boldsymbol {\hat {e}}_2 \,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}$ are 1, while for 
$\boldsymbol {\hat {e}}_3\,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}$, it is 0. Thus, the magnetic field aligns itself either with 
$\boldsymbol {\hat {e}}_1$ or 
$\boldsymbol {\hat {e}}_2$. In saturation, the p.d.f. of 
$\boldsymbol {\hat {e}}_2\,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}$ remains largely unchanged, while the part of the field that is not aligned with 
$\boldsymbol {\hat {e}}_2$ aligns itself predominantly between 
$\boldsymbol {\hat {e}}_1$ and 
$\boldsymbol {\hat {e}}_3$. In contrast, the unlimited Braginskii-MHD run mimics the statistics of the saturated MHD dynamo throughout its evolution, with the alignment angles hardly changing at all between the exponential-growth regime and saturation.
Figure 15. Probability distribution function of the rate-of-strain eigenvector alignment cosines 
$\cos \theta _i = |{\boldsymbol {\hat {e}}}_i \,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}|$ for simulations of isotropic MHD (top) and unlimited Braginskii MHD (bottom) in the kinematic stage (solid line) and saturated state (dotted line). Eigenvalues are ordered so that 
$\lambda _1 > \lambda _2 > \lambda _3$. The vertical dotted line marks 
$\theta _i = 45^\circ$.
These alignments can also be seen in figure 16, which shows time-averaged joint p.d.f.s of 
$|\boldsymbol {\hat {e}}_1\,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}|^2$ and 
$|\boldsymbol {\hat {e}}_2\,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}|^2$ both for MHD simulations (a and b) and for unlimited Braginskii MHD ones (c and d) in the kinematic stage (a and c) and in the saturated state (b and d). To highlight the alignment of dynamically important fields, these p.d.f.s are conditioned on regions of high field strength, 
$B > 2 B_\mathrm {rms}$. As the alignment of 
${\boldsymbol {\hat {b}}}$ with 
$\boldsymbol {\hat {e}}_3$ can be related to the other two alignments via the identity 
$|\boldsymbol {\hat {e}}_1\,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}|^2 + |\boldsymbol {\hat {e}}_2\,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}|^2 + |\boldsymbol {\hat {e}}_3\,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}|^2 = 1$, these p.d.f.s contain all the information about the field-alignment statistics (e.g. any probability density appearing in the bottom left corner of these plots denotes magnetic fields primarily aligned with 
$\boldsymbol {\hat {e}}_3$). It is clear that, in the MHD simulations, the magnetic field tends to align itself principally with either 
$\boldsymbol {\hat {e}}_1$ or 
$\boldsymbol {\hat {e}}_2$ in the kinematic stage, rather than with some combination of all three eigenvectors. The statistics in saturation vary significantly for different amounts of viscosity; for 
${Re}\sim 1$, the dynamo saturates by having the magnetic field align itself uniformly between 
$\boldsymbol {\hat {e}}_1$ and 
$\boldsymbol {\hat {e}}_2$ along the 
$|\boldsymbol {\hat {e}}_3 \,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}| = 0$ contour, while for 
${Re} \gg 1$, it aligns itself between 
$\boldsymbol {\hat {e}}_1$ and 
$\boldsymbol {\hat {e}}_3$ (namely, 
$\theta _1 \approx 40^\circ$ and 
$\theta _3 \approx 50^\circ$). The simulation with intermediate Reynolds number (centre panel) contains regions where both of these situations occur.
Figure 16. Time-averaged joint p.d.f.s of the cosines of the angles between the magnetic field and the first two ordered eigenvectors of the rate-of-strain tensor for the MHD simulations MHD1/2/4, and for the unlimited Braginskii-MHD simulations U1, U4 and U1b. The p.d.f.s are conditioned on high magnetic energy (
$B > 2 B_\mathrm {rms}$) to probe the strongest dynamo-generated fields. Dashed lines correspond to 
$|\boldsymbol {\hat {e}}_1\,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}| = |\boldsymbol {\hat {e}}_2\,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}|$; dotted lines correspond to solutions of 
$|\boldsymbol {\hat {e}}_1\,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}|^2 + |\boldsymbol {\hat {e}}_2\,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}|^2+ |\boldsymbol {\hat {e}}_3\,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}|^2 =1$ for fixed 
$|\boldsymbol {\hat {e}}_3\,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}|$ (i.e. 
$|\boldsymbol {\hat {e}}_3\,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}|$ is fixed along these lines). The simulations shown in figure 15 are marked by the thicker boxes.
For unlimited Braginskii MHD with large parallel viscosity 
$\mu _{{\rm B}}$ (left side of panels (c) and (d) in figure 16), the statistics of the field alignment hardly change between the kinematic stage and saturation, and both are similar to those found in the high-
$Re$ MHD dynamo (
$\mu ^{-1}=1500$) in its saturated state (figure16b, bottom right side). (The only notable difference is that the enhancement in probability density localized near 
$|\boldsymbol {\hat {e}}_1 \,{\boldsymbol {\cdot }}\, {\boldsymbol {\hat {b}}}|^2 \approx 0.6$ in the latter is slightly less pronounced in the Braginskii-MHD case.) This again emphasizes that the role of the parallel Braginskii viscosity is to render the statistics of the magnetic field nearly identical to those found in the saturated state of the MHD dynamo.Footnote 15
Finally, figure 17 displays the p.d.f.s of the eigenvalues 
$\lambda _1 > \lambda _2 > \lambda _3$ of the rate-of-strain tensor 
$\boldsymbol {\mathsf {S}}$ for 
${Re}\gg {1}$ MHD (run MHD4; figure 17a) and for unlimited Braginskii MHD (run U1; figure 17b) in the kinematic (solid line) and saturated states (dashed line). The statistics of the rate-of-strain tensor change noticeably between the two stages in the isotropic-MHD simulation, while the change is minimal (if any) in the unlimited Braginskii-MHD system. Once again, the unlimited Braginskii-MHD dynamo behaves throughout its evolution much like the isotropic-MHD dynamo in its saturated state.
Figure 17. Probability distribution function of the rate-of-strain eigenvalues for simulations of isotropic MHD (a) and unlimited Braginskii MHD (b) in the kinematic stage (solid lines) and saturated state (dotted lines). Eigenvalues are ordered so that 
$\lambda _1 > \lambda _2 > \lambda _3$.
4.3 Modified Kazantsev–Kraichnan model for the Braginskii-MHD dynamo
In light of the main conclusion of §§ 3.3–4.2 – that the turbulence statistics during the ‘kinematic’ phase of the unlimited Braginskii-MHD dynamo are strikingly similar to those found in the saturated state of the isotropic-MHD dynamo – we propose that the ‘kinematic’ stage of the unlimited Braginskii-MHD dynamo can be described by an analytic model originally developed to describe the saturated state of the MHD dynamo by Schekochihin et al. (Reference Schekochihin, Cowley, Taylor, Hammett, Maron and McWilliams2004a,Reference Schekochihin, Cowley, Taylor, Maron and McWilliamsb). This model uses a modified form of the Kazantsev–Kraichnan (Kazantsev Reference Kazantsev1968; Kraichnan Reference Kraichnan1974) velocity correlation tensor
\begin{equation} \overline{u^i(t,{\boldsymbol{x}})u^j(t',{\boldsymbol{x}}')} = \delta(t-t')\kappa^{ij}({\boldsymbol{x}}-{\boldsymbol{x}}'), \end{equation}
where an overline denotes the ensemble average, that allows anisotropic statistics with respect to the magnetic-field direction 
${\boldsymbol {\hat {b}}}$: the Fourier transform of 
$\kappa ^{ij}({\boldsymbol {x}}-{\boldsymbol {x}}')$ is
\begin{equation} \kappa^{ij}({\boldsymbol{k}}) = \kappa^{({\rm i})}(k,|\xi|) \left(\delta^{ij} - \hat{k}_i \hat{k}_j \right) + \kappa^{({\rm a})}(k,|\xi|) \left({b}^i {b}^j + \xi^2 \hat{k}_i \hat{k}_j - \xi {b}^i \hat{k}_j - \xi \hat{k}_i {b}^j \right), \end{equation}
where 
$\hat {{\boldsymbol {k}}} \doteq {\boldsymbol {k}}/k$, 
$\xi \doteq \hat {{\boldsymbol {k}}}\,{\boldsymbol {\cdot }}\,{\boldsymbol {\hat {b}}}$, and 
${b}^i$ denotes the 
$i$th component of the magnetic-field unit vector 
${\boldsymbol {\hat {b}}}$. The 
$k$- and 
$|\xi |$-dependent amplitudes 
$\kappa ^{({\rm i})}$ and 
$\kappa ^{({\rm a})}$ quantify the sizes of the isotropic and anisotropic components of the correlation tensor. Note that 
$k_i\kappa ^{ij} = k_j \kappa ^{ij} = 0$ as a result of incompressibility. When this model was put forward, the idea was that the dynamically important magnetic field in saturation exerts a feedback on the velocity by means of the Lorentz force, biasing the velocity's statistics with respect to the local field direction. Here, we propose that this model might accurately describe the impact of the Braginskii viscosity on the turbulent flow during the kinematic stage.
An outline and derivation of the model is given in § 4.3.1 and appendix C. The model is then solved in the general case of finite dynamo growth rate and non-negligible resistivity in § 4.3.2. In §§ 4.3.3–4.4, we use this solution and the results of numerical simulation to show that the dynamo fails in the limit 
${Re}_\parallel /{Re} \rightarrow 0$.
4.3.1 Formulation of the model
We first obtain equations for the evolution of magnetic-field fluctuations, the wavevector of these fluctuations, and the scale-dependent magnetic-field direction as functions of the fluctuating velocity field. These are used to determine the statistics of the magnetic-field strength, structure and direction, from which we then derive an evolution equation for the magnetic-energy spectrum.
To simplify the calculation, we focus only on the straight regions of the resulting folded magnetic field, so that spatial variations in 
${\boldsymbol {\hat {b}}}$ are limited to changes of sign. Because 
${\boldsymbol {\hat {b}}}$ arises in the momentum equation as part of the dyadic 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}$, these changes of sign cancel each other out, and so 
${\boldsymbol {\hat {b}}}$ can be taken to depend only on time. This approximation is not as drastic as it may appear, as the bends of the magnetic folds occupy a small fraction of the total volume (Schekochihin et al. Reference Schekochihin, Cowley, Maron and Malyshkin2002b, Reference Schekochihin, Cowley, Taylor, Maron and McWilliams2004b). The evolution equations then follow straightforwardly from the non-ideal induction equation after adopting the ansatz
\begin{equation} {\boldsymbol{B}}(t,{\boldsymbol{x}}) = {\boldsymbol{\hat{b}}}(t) \int \textrm{d}^3{\boldsymbol{k}}' \, {B}(t,{\boldsymbol{k}}') \exp({\textrm{i} {\boldsymbol{x}}\,{\boldsymbol{\cdot}}\, {{\boldsymbol{k}}}(t,{\boldsymbol{k}}')}) , \end{equation}
where 
${\boldsymbol {k}}(t,{\boldsymbol {k}}')$ is the wavevector that evolves in time from the initial value 
${\boldsymbol {k}}'$. Assuming statistical homogeneity and setting 
${\boldsymbol {x}} = 0$ results in the closed set of equations
\begin{gather} \partial_t {B} = {b}^i {b}^m (\partial_m u^i) {B} - \eta {k}^2 {B}, \end{gather}
\begin{gather}\partial_t {k}_m = - (\partial_m u^i) {k}_{i}, \end{gather}
\begin{gather}\partial_t {b}^i = {b}^m (\partial_m u^i) - {b}^l {b}^m (\partial_m u^l) {b}^i, \end{gather}
where we have used the Einstein summation convention for repeated indices. The upper and lower indices are used for convenience of keeping track of the summations. The joint probability density function of 
$B$, 
${\boldsymbol {k}}$ and 
${\boldsymbol {\hat {b}}}$ must then have the form
\begin{equation} \mathcal{P}({B},{{\boldsymbol{k}}},{\boldsymbol{\hat{b}}}) = \delta(|{\boldsymbol{\hat{b}}}|^2 -1)\delta({\boldsymbol{\hat{b}}}\,{\boldsymbol{\cdot}}\,{{\boldsymbol{k}}})(4{\rm \pi}^2 {k})^{-1}P({B},{k}), \end{equation}
and a closed equation for 
$\mathcal {P}(B,k)$ can be derived using (4.1). Several quantities of interest may then be calculated by taking the appropriate moments of 
$\mathcal {P}(B,k)$, e.g. the magnetic-energy spectrum
\begin{equation} M(k) \doteq \frac{1}{2} \int_0^\infty \textrm{d}{B}\, {B}^2 P({B},{k}) . \end{equation}
A detailed derivation of the evolution equation for 
$\mathcal {P}({B},k)$ is given in appendix C.
For 
$k\gg k_\mu$ (i.e. at sub-viscous scales), the evolution equation for the magnetic-energy spectrum obtained via (4.6) is
\begin{align} \frac{\partial M}{\partial t} = \frac{\gamma_\perp}{8} \frac{\partial }{\partial k} \left[ \left( 1 + 2\sigma_\parallel\right) k^2 \frac{\partial M}{\partial k} - \left( 1 + 4\sigma_\perp + 10\sigma_\parallel \right) kM \right] + 2 \left( \sigma_\perp + \sigma_\parallel\right)\gamma_\perp M - 2\eta k^2 M, \end{align}where
\begin{align} \gamma_\perp = \int\frac{\textrm{d}^3{\boldsymbol{k}}}{(2{\rm \pi})^3} \, k^2_\perp \kappa_\perp({\boldsymbol{k}}) , \quad \sigma_\perp = \frac{1}{\gamma_\perp} \int\frac{\textrm{d}^3{\boldsymbol{k}}}{(2{\rm \pi})^3} \, k^2_\parallel \kappa_\perp({\boldsymbol{k}}) , \quad \sigma_\parallel = \frac{1}{\gamma_\perp} \int\frac{\textrm{d}^3{\boldsymbol{k}}}{(2{\rm \pi})^3} \, k^2_\parallel \kappa_\parallel({\boldsymbol{k}}) . \end{align}
Here, 
$k_\parallel = k\xi$, 
$k_\perp = k(1-\xi ^2)^{1/2}$ and
\begin{align} \kappa_\perp({\boldsymbol{k}}) &= \frac{1}{2} \left( \delta^{ij} - {b}^i{b}^j\right)\kappa^{ij}({\boldsymbol{k}}) \nonumber\\ &= \frac{1}{2} \left[ \left(1+\xi^2\right)\kappa^{({\rm i})}(k,|\xi|) + \xi^2\left(1-\xi^2\right) \kappa^{({\rm a})}(k,|\xi|) \right] , \end{align}
\begin{align} \kappa_\parallel({\boldsymbol{k}}) &= \frac{1}{2} {b}^i {b}^j \kappa^{ij}({\boldsymbol{k}})\notag\\ &= \frac{1}{2} \left(1-\xi^2\right) \left[\kappa^{({\rm i})}(k,|\xi|) + \left(1-\xi^2\right) \kappa^{({\rm a})}(k,|\xi|) \right] \end{align}
are the correlation spectra of velocities respectively perpendicular and parallel to the local magnetic field, and 
$\kappa ^{ij}$ is defined by (4.1) and (4.2). The quantity 
$\gamma _\perp$ measures the perpendicular variation of the perpendicular velocities, and thus gives the mixing rate of the magnetic field. It is at this rate that interchange-like motions shuffle and bring direction-reversing magnetic fields close enough together for them to diffuse resistively. The quantities 
$\sigma _\perp$ and 
$\sigma _\parallel$ measure the relative strengths compared to this mixing rate of the field-aligned stretching rates of the perpendicular and parallel velocities. For reference, the isotropic case with 
$\kappa ^{({\rm i})} = \kappa ^{({\rm i})}(k)$ and 
$\kappa ^{({\rm a})} = 0$ has 
$\sigma _\perp = 2/3$, 
$\sigma _\parallel = 1/6$ and 
$\gamma _\perp = (6/5)\overline{\gamma }$, where
\begin{equation} \overline{\gamma} = \frac{1}{3}\left[\int_{k_0}^\infty \textrm{d} k \, k^2 E(k)\right]^{1/2}, \end{equation}
$k_0 \doteq \ell _0^{-1} = 2 {\rm \pi} / L$ is the wavenumber of the outer scale, 
$E(k) = \int \textrm {d} \varOmega _{{\boldsymbol {k}}}\, k^2\delta _{ij}I^{ij}({\boldsymbol {k}})$ is the kinetic-energy spectrum and 
$\textrm {d} \varOmega _{{\boldsymbol {k}}}$ is an element of solid angle in 
$k$ space.
4.3.2 Solution of the model and comparison with simulation results
Equation (4.7) is solved by assuming 
$M(k,t) \propto \textrm {e}^{\gamma t}$, where 
$\gamma$ is the magnetic-energy growth rate. Requiring 
$M(k) \rightarrow 0$ as 
$k \rightarrow \infty$ leads to the solution
\begin{equation} M(k) = k^s \,\textrm{e}^{\gamma t}K_{\textrm{i} \varpi}(k/k_\eta), \end{equation}
where 
$K_{\textrm {i} \varpi }$ is the Macdonald function, 
$k_\eta = [(1+2\sigma _\parallel )\gamma _\perp / 16 \eta ]^{1/2}$ is the (inverse) resistive scale,
\begin{equation} s = \frac{2(\sigma_\perp + 2 \sigma_\parallel)}{1+2 \sigma_\parallel} \end{equation}is the spectral index and
\begin{align} \varpi = \frac{1}{1 + 2\sigma_\parallel}\left[16(\sigma_\perp + \sigma_\parallel)(1+2\sigma_\parallel) + 8(\gamma/\gamma_\perp)(1+2\sigma_\parallel) - (1+2\sigma_\perp + 6 \sigma_\parallel)^2 \right]^{-1/2}. \end{align}
In order for this solution to be uniquely determined, the four coefficients 
$\gamma$, 
$\gamma _\perp$, 
$\sigma _\perp$ and 
$\sigma _\parallel$ must be determined. One constraint is obtained by imposing a zero-flux boundary condition at small 
$k$, i.e. the term in square brackets in (4.7) is set to zero:
\begin{equation} \frac{k_\ast}{k_\eta}K'_{\textrm{i} \varpi}\left(\frac{k_\ast}{k_\eta}\right) - \frac{1+2 \sigma_\perp + 6\sigma_\parallel}{1 + 2\sigma_\parallel}K_{\textrm{i} \varpi}\left(\frac{k_\ast}{k_\eta}\right)=0, \end{equation}
where 
$k_\ast$ is the infrared cutoff of the magnetic-energy spectrum. In practice, we let 
$k_\ast = k_0$. Three more constraints must then be provided.Footnote 16
We can gauge the validity of this model (and, therefore, of the idea of the local anisotropization) by comparing the magnetic-energy spectra obtained from our simulations to those predicted by (4.12) with (4.15). To do so, we specify the three remaining constraints as follows. First, we model the mixing rate by
\begin{equation} \gamma_\perp = c (\epsilon/\mu)^{1/2}, \end{equation}
where 
$\epsilon$ is the energy injection rate and 
$c$ is an adjustable parameter (we use 
$c=0.3$ for all calculations, as did Schekochihin et al. Reference Schekochihin, Cowley, Taylor, Maron and McWilliams2004b). Secondly, we assume 
$\sigma _\perp = 4\sigma _\parallel$, as in the isotropic case. We justify this simple closure by noting that the structure functions in figure 11 indicate that the parallel scales of velocity variations are predominantly on the outer scale (i.e. 
$k_\parallel \sim k_0$) and that, at these scales, we indeed measure 
$\sigma _\perp \approx 4 \sigma _\parallel$ (using (4.17) below; see figure 19). Finally, the magnetic-energy growth rate 
$\gamma$ is obtained by time averaging the measured growth rate during the exponential-growth phase of the simulated dynamo. With this information, (4.15) can be solved to obtain 
$\sigma _\parallel$ (and thus 
$\sigma _\perp$), which then completely determines the magnetic-energy spectra.
The results of this calculation and a comparison of the model and simulated spectra are shown in figure 18. Surprisingly good quantitative agreement is found between the time-averaged spectra obtained from the simulations and those predicted by (4.12).
Figure 18. Time-averaged magnetic-energy spectra from simulations of unlimited Braginskii MHD (solid line), along with predicted spectra calculated from solutions of (4.13) and (4.15) (dotted lines) using the assumption 
$\sigma _\perp = 4 \sigma _\parallel$. Curves with varying 
$\mu$ have been shifted vertically for clarity.
4.3.3 Is unlimited Braginskii-MHD dynamo possible when ${Re}\rightarrow \infty$?
Numerical investigation of solutions to (4.15) has shown that, in the limit 
${Rm} \rightarrow \infty$, the viability of the dynamo – that is, whether its growth rate is positive or negative – depends only on the relative size of the parallel (stretching) motions compared to the perpendicular (mixing) ones, rather than on their absolute magnitude. In particular, as 
$\gamma _\perp$ increases and the mixing motions become more important (thereby decreasing 
$\sigma _\parallel$ and 
$\sigma _\perp$), the dynamo growth rate will decrease if the parallel motions do not increase commensurately. Physically, this is due to an enhancement of small-scale mixing, which tends to bring field lines closer together and promotes their resistive annihilation. This suggests that there is no unlimited Braginskii-MHD dynamo in the limit 
$\mu \rightarrow 0$, a somewhat paradoxical result whose physical origin is further discussed at the end of § 4.4. This scenario can be inferred from figure 8, where the unlimited Braginskii-MHD system with the largest isotropic viscosity (and thus the weakest mixing motions) has the fastest growing dynamo. It is also consistent with the observation that the structure functions in figure 11 reveal that parallel variations of the velocity field reside mainly at the outer scale, while its perpendicular variations tend to reside at the viscous cutoff. This indicates that 
$\sigma _\parallel$, 
$\sigma _\perp \rightarrow 0$ as 
${Re} \rightarrow \infty$ if 
${Re}_\parallel$ is held fixed.
This asymptotic behaviour is inferred by determining the dependence of 
$\sigma _\perp$, 
$\sigma _\parallel$ on 
$\mu ^{-1}$ from our simulations and then extrapolating the trends. We do so using two separate methods, which fit the model outlined above to our simulation data in fundamentally different ways and are compared in figure 19. The first method is to use the model solution (4.12) and the constraints given by (4.15), (4.16) and 
$\sigma _\perp = 4 \sigma _\parallel$. The final undetermined variable, the magnetic-energy growth rate 
$\gamma$, is obtained directly by measuring the growth rate from simulation data. The second method is to extract 
$\gamma _\perp$, 
$\sigma _\parallel$, and 
$\sigma _\perp$ directly from the simulations by modifying (4.8) to account for the finite correlation time 
$\tau _{c}$ of the velocity field. In this case, the Fourier-space velocity correlation tensor 
$I^{ij}({\boldsymbol {k}}) \doteq \int \textrm {d}^3 {\boldsymbol {k}}\exp ({\textrm {i} {\boldsymbol {k}}\,{\boldsymbol {\cdot }}\,({\boldsymbol {x}}-{\boldsymbol {x}}')})\overline {u^i({\boldsymbol {x}})u^{j}({\boldsymbol {x}}')}$ is related to 
$\kappa ^{ij}({\boldsymbol {k}})$ via 
$I^{ij}({\boldsymbol {k}}) = \tau ^{-1}_{c} \kappa ^{ij}({\boldsymbol {k}})$ (cf. (4.1)). We then take the correlation time of a given type of motions to be proportional to the ‘turnover’ time of these motions; e.g. the correlation time of the perpendicular variations of velocity perpendicular to the magnetic field is 
$\tau _{c,\perp \perp } = C_{\perp \perp } \gamma ^{-1}_\perp$, where 
$C_{\perp \perp }$ is an undetermined coefficient. In principle, the correlation time may depend on both the alignment of the motion with respect to the magnetic field and the alignment of its variation, so we must also introduce 
$\tau _{c,\parallel \perp } = C_{\parallel \perp } (\sigma _\perp \gamma _\perp )^{-1}$ and 
$\tau _{c,\parallel \parallel } = C_{\parallel \parallel } (\sigma _\parallel \gamma _\perp )^{-1}$, where 
$C_{\parallel \perp }$ and 
$C_{\parallel \parallel }$ are further undetermined constraints. Then (4.8) become
\begin{align} \gamma_\perp &= \left[C_{\perp\perp}\int\frac{\textrm{d}^3{\boldsymbol{k}}}{(2{\rm \pi})^3} \, k^2_\perp I_\perp\right]^{1/2} = \left[C_{\perp\perp} \langle | (\boldsymbol{\mathsf {I}} - {\boldsymbol{\hat{b}}}{\boldsymbol{\hat{b}}})\,{\boldsymbol{\cdot}}\, {\boldsymbol \nabla}{\boldsymbol{u}}\,{\boldsymbol{\cdot}}\, (\boldsymbol{\mathsf {I}} - {\boldsymbol{\hat{b}}}{\boldsymbol{\hat{b}}})|^2 \rangle \right]^{1/2} , \end{align}
\begin{align}\sigma_\perp &= \frac{1}{\gamma_\perp}\left[C_{\parallel\perp} \int\frac{\textrm{d}^3{\boldsymbol{k}}}{(2{\rm \pi})^3} \,k^2_\parallel I_\perp\right]^{1/2} = \frac{1}{\gamma_\perp}\left[C_{\parallel \perp}\langle | {\boldsymbol{\hat{b}}}\,{\boldsymbol{\cdot}}\, {\boldsymbol \nabla}{\boldsymbol{u}}\,{\boldsymbol{\cdot}}\, (\boldsymbol{\mathsf {I}} - {\boldsymbol{\hat{b}}}{\boldsymbol{\hat{b}}})|^2 \rangle \right]^{1/2}, \end{align}
\begin{align}\sigma_\parallel &= \frac{1}{\gamma_\perp} \left[C_{\parallel\parallel}\int\frac{\textrm{d}^3{\boldsymbol{k}}}{(2{\rm \pi})^3} \, k^2_\parallel I_\parallel\right]^{1/2} = \frac{1}{\gamma_\perp}\left[C_{\parallel \parallel}\langle | {\boldsymbol{\hat{b}}}{\boldsymbol{\hat{b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla} {\boldsymbol u} |^2 \rangle \right]^{1/2}, \end{align}
where 
$I_\perp$ and 
$I_\parallel$ are defined in a similar way to (4.9) and (4.10), respectively. In order to eliminate some of the undetermined coefficients, we make the assumptions 
$C_{\parallel \parallel } = (1/6)C_{\perp \perp }$, 
$C_{\parallel \perp } = (2/3) C_{\perp \perp }$, chosen to ensure that 
$\sigma _\perp = 2/3$ and 
$\sigma _\parallel = 1/6$ in the isotropic case 
$I^{{({\rm a})}} =0$.
Figure 19. Modified Kazantsev–Kraichnan model parameters 
$\gamma _\perp$ (top of panel a), 
$\sigma _\perp$ (top of panel b) and 
$\sigma _\parallel$ (bottom of panel a), calculated using the velocity gradient from simulations of unlimited Braginskii MHD according to (4.17), assuming 
$C_{\perp \perp }= 1$ (purple crosses), as well as from (4.15) assuming 
$\sigma _\perp = 4 \sigma _\parallel$ (green ‘x’s). The former values are time averaged over the exponential-growth phase. Bottom right: magnetic-energy growth rate 
$\gamma$ normalized by 
$\gamma _\perp$ for simulations of unlimited Braginskii MHD, time averaged over the exponential-growth phase.
Figure 19 shows the results of these measurements. The assumption 
$c = 0.3$ in (4.16) appears to agree well with 
$\gamma _\perp$ estimated directly from the simulation data. The predicted value of 
$\sigma _\parallel$ also closely matches the estimated value obtained from the simulated 
${\boldsymbol \nabla } {\boldsymbol {u}}$. Finally, and most importantly for this discussion, the scaled magnetic-energy growth rate decreases along with the isotropic viscosity 
$\mu$, indicating that there exists a critical value of 
$\mu$ at which the dynamo ceases to operate.Footnote 17
4.4 Stokes flow: is an unlimited Braginskii-MHD dynamo possible when ${Re}_\parallel < 1?$
It was predicted in § 4.3.2 that the unlimited Braginskii-MHD dynamo ceases to operate in the limit 
${Re} \rightarrow \infty$, 
${Re}_\parallel = \mathrm {const}$. In this section, we ask whether the same is true in the opposite, ‘Stokes-flow’ limit of 
${Re} = \mathrm {const}$, 
${Re}_\parallel \rightarrow 0$. We will show that it is: at fixed 
${Re}$, the dynamo does cease to operate at sufficiently small 
${Re}_\parallel$ (although it may be reactivated by decreasing 
${Re}$ to values comparable to 
${Re}_\parallel$).
In the isotropic-MHD case, the correlation time of a white-noise-forced flow vanishes in the limit 
${Re}\rightarrow {0}$ (provided the forcing amplitude is scaled with viscosity) and the dynamo is well described by the Kazantsev–Kraichnan model. Provided there is sufficient scale separation between the forcing and resistive scales, the dynamo will continue to operate during this stage regardless of 
${Re}$ (Schekochihin et al. Reference Schekochihin, Cowley, Hammett, Maron and McWilliams2002a). To investigate whether a similar result holds for the Braginskii system, we perform a set of low-
${Re}$ isotropic-MHD and low-
${Re}_\parallel$ unlimited Braginskii-MHD simulations at reduced resolution (
$N_\mathrm {cell} = 112^3$). For each run, the forcing amplitude 
$\varepsilon$ is adjusted so that 
$u_\mathrm {rms}\sim 1$ in steady state. Details of these runs are given in the last block of table 1.
The evolution of the magnetic energy for simulations in this Stokes-flow regime is plotted in figure 20. As predicted by the Kazantsev–Kraichnan model, the behaviour of the isotropic MHD simulations changes little, provided 
$u_\mathrm {rms}$ is kept constant across all runs. They all have magnetic-energy growth and reach saturation. The unlimited Braginskii-MHD dynamo, on the other hand, operates only beyond a certain critical value of 
${Re}_\parallel \gtrsim {1}$, with an abrupt cut-off in the growth rate for smaller 
${Re}_\parallel$. Well above this cut-off, the scaling of the growth rate still follows the expected Kolmogorov (Reference Kolmogorov1941) scaling 
${Re}_\parallel ^{1/2}$ (figure 20, bottom of panel b). These results, as well as those from the analytic model solved in § 4.3.2, indicate that the dynamo is only viable for moderate values of the ratio 
${Re}_\parallel /{Re}$: too small a value results in a dynamo with too much mixing (relative to stretching) and thus excessive resistive dissipation.
Figure 20. (a) Evolution of the magnetic energy for various 
${Re} < 1$ (isotropic MHD, solid lines) and 
${Re}_\parallel < 1$, 
${Re} \gg 1$ (unlimited Braginskii-MHD, dotted and dash-double-dotted lines) simulations; 
$\mu ^{-1}_{{\rm B}}$ is varied from 
$20$ down to 
$0.5$ for simulations of unlimited Braginskii MHD. For 
$\mu ^{-1}_{{\rm B}} \lesssim 4$ and 
$\mu ^{-1} = 1500$, the unlimited Braginskii-MHD simulations do not exhibit a dynamo. (b) Time-averaged growth rate 
$\gamma = \textrm {d} \ln B^2_\mathrm {rms}/\textrm {d} t$ of the magnetic energy during the kinematic stage as a function of parallel Reynolds number for unlimited Braginskii MHD in the Stokes regime on a linear (top) and logarithmic (bottom) scale.
The extinction of dynamo in both limits 
${Re}_\parallel \rightarrow 0$, 
${Re}= \mathrm {const}$ and 
${Re}_\parallel =\mathrm {const}$, 
${Re}\rightarrow \infty$ can be thought of as the result of a strong anisotropization of the underlying turbulence, leading to motions that are more two-dimensional than three-dimensional (see § 4.1). Interestingly, simulations that have 
${Re}_\parallel$, 
${Re} < 1$ yet 
${Re}_\parallel /{Re}\sim 1$ still experience a viable dynamo; compare the dotted and dot-dashed light blue lines in figure 20. This shows that it is the relative size of the stretching to mixing time scales, rather than their absolute values, that determine whether there is dynamo (cf. § 4.3). Having an isotropic viscosity comparable to an anisotropic viscosity allows the momentum to diffuse more isotropically, undermining the anisotropic viscosity's tendency to make the flow more two-dimensional. Indeed, the case with 
${Re} \sim {Re}_\parallel$ is more similar to the isotropic 
${Pm} \gtrsim 1$ dynamo, where the contributions of both stretching and mixing come from the same range of scales.
4.4.1 Mixing through a critically balanced pressure-anisotropy cascade
Despite not supporting dynamo action for 
${Re}_\parallel \lesssim {1}$, the unlimited Braginskii-MHD system still features a cascade of turbulent energy to small scales. The kinetic- and magnetic-energy spectra of various unlimited Braginskii-MHD simulations in the Stokes-flow regime are displayed in figure 21(a). As 
${Re}_\parallel$ is decreased, the kinetic-energy spectrum beyond the parallel viscous scale approaches a power law with the exponent 
$-5/3$, a slope seen most clearly in the 
$\mu _{{\rm B}} = 0.5$ run (which is a failed dynamo; see the light-blue dot-dot-dashed line in figure 20). For that run, the scaling of the wavevector anisotropy (bottom-right panel) is measured (via (2.11)) to be 
$\ell _\parallel \propto \ell _\perp ^{2/3}$. These scaling exponents, 
$-5/3$ and 
$2/3$, match those predicted by the Goldreich & Sridhar (Reference Goldreich and Sridhar1995) theory for MHD turbulence in the presence of a strong guide field, where the linear (Alfvén-wave) and nonlinear (turbulence) time scales are taken to be comparable (‘critical balance’). Here, though, there is no guide field, and the Alfvén-wave time scale is much too long at these low magnetic-field strengths to matter dynamically. So why then do these scalings arise? If the observed scalings are representative of a critically balanced state, just what time scales is it a critical balance between?
Figure 21. (a) kinetic (left) and magnetic (right) energy spectra for various unlimited Braginskii-MHD simulations with small and large parallel Reynolds number 
${Re}_\parallel$. (b) perpendicular and parallel structure functions (left) and spectral anisotropy scalings (right) for the same runs (see (2.11)).
To answer these questions, we compute the shell-filtered kinetic-energy transfer functions (see (2.6) and (2.7)) for the unlimited Braginskii-MHD system in the Stokes-flow regime. These are shown in figure 22.Footnote 18 In the failed Braginskii-MHD dynamo that nevertheless exhibits a robust turbulent cascade (left panels), we find from the r.m.s. transfer functions (bottom row) that in the inertial range there is a scale-by-scale balance between the hydrodynamic nonlinearity (green line) and the (unlimited) Braginskii viscous stress (red line). Thus, we can picture a cascade of MHD turbulence consisting of Alfvén-like waves propagating along the local magnetic-field direction, with ‘critical balance’ occurring between the linear frequency associated with the pressure anisotropy, 
${\sim }k_\parallel \sqrt {{\rm \Delta} p}$, and the nonlinear frequency associated with the eddy-turnover time, 
${\sim }k_\perp u_\perp$. Effectively, the unchecked viscous stress supplants the (negligible) Maxwell stress in setting the characteristic linear frequency of the plasma.Footnote 19 Adopting the usual assumption of a conservative cascade with 
$k_\perp u^3_\perp \sim \textrm {const}$, this balance gives 
$k_\parallel \sqrt {{\rm \Delta} p} \sim k_\perp u_\perp \sim k_\perp ^{2/3}$. Noting that the pressure anisotropy in this simulation is dominated by its mean value (and thus is largely scale independent), the scaling 
$k_\parallel \sim k^{2/3}_\perp$ results (as is indeed seen in figure 21 for this run). By contrast, the transfer functions in the complementary run with stronger isotropic viscosity (run Usg
$^\ast$, which is a dynamo – see the light-blue dotted line in figure 20) are qualitatively similar to those displayed in the lower panels (
$\mu ^{-1}_{{\rm B}} = 20$, 
$\mu ^{-1} = 96$) of figure 10.
Figure 22. Root-mean-square shell-filtered kinetic-energy transfer functions 
$\mathcal {T}_k$ (see (2.7)) for the unlimited Braginskii-MHD system in the Stokes-flow regime employing 
$\mu _{{\rm B}}^{-1} = 0.5$ and 
$\mu ^{-1} = 1500$ (
$\mu ^{-1} = 4$) in the left (right) frame. The runs chosen here are one that corresponds to magnetic-field decay (USg, left) and one that is a viable dynamo (USg
$^*$, right).
5 Summary and discussion
In this work, we have studied the fluctuation dynamo in a weakly collisional, magnetized plasma using incompressible Braginskii-MHD simulations and semi-analytic modelling. The pressure anisotropy, and thus the parallel viscous stress, was either hard-wall limited to lie within the firehose- and mirror-instability thresholds or allowed to venture freely beyond those thresholds. While the latter option has traditionally been considered unphysical due to the plethora of evidence – both theoretical and observational – that kinetic-instability thresholds are well respected in collisionless, magnetized plasmas, its study is important for (at least) four reasons. First, the unlimited case teaches us about the fundamental dynamical properties of the pressure anisotropy stress, which has received relatively little attention in the past (particularly in the case of dynamo). Second, it offers an additional point of comparison to dynamo behaviour in isotropic MHD and limited Braginskii MHD; in particular, our finding of its similarity to the saturated state of the MHD dynamo affords a better understanding of that more traditional case. Third, many aspects of its evolution are remarkably similar to those found in the more ab initio hybrid-kinetic simulations of plasma dynamo presented in Paper I; this is fortunate, as the unlimited Braginskii-MHD simulations provide a better controlled and more economical test bed with which to diagnose the field and flow statistics in this regime. Finally, we have argued that a significant period in the dynamo amplification of the intracluster magnetic field occurs at a time when the plasma 
$\beta$ is too large for kinetic instabilities to regulate the pressure anisotropy efficiently enough to pin it tightly near its 
${\sim }1/\beta$ stability boundaries. During this phase, a constant collision frequency that partially restrains the pressure anisotropy (as in the unlimited simulations presented above) may in fact be the more realistic ‘closure’. Indeed, much of the evolution of the collisionless dynamo found in Paper I using a hybrid-kinetic approach occurred during such a phase, with a suppressed parallel rate of strain, an anisotropization of the flow velocity, an imperfect regulation of the pressure anisotropy, and a Kolmogorov-like cascade of perpendicular kinetic energy – all of which are manifest in our unlimited Braginskii-MHD runs.
Our main conclusions are as follows:
(i) The chaotic flow driven by large-scale forcing produces highly intermittent and structured magnetic fields, which are organized into folds and grow exponentially until the Lorentz tension force is strong enough to back-react dynamically on the velocity field. This folded structure, a hallmark of the
${Pm}\gtrsim {1}$ fluctuation dynamo, persists into the saturated state. These results hold regardless of whether the plasma is described by isotropic MHD or Braginskii MHD with either limited or unlimited pressure anisotropy, so long as the ratio of the isotropic and Braginskii viscosities (or the ratio of the parallel and perpendicular Reynolds numbers; see (v) below) is not too small.(ii) Hard-wall limiters on the parallel viscosity, intended to mimic the rapid regulation of pressure anisotropy by kinetic Larmor-scale instabilities otherwise not properly captured in Braginskii MHD, effectively reduce the Braginskii-MHD dynamo to its
${Re}\gg {1}$, ${Pm}\sim {1}$ MHD counterpart. With the exception of some minor differences, such as a slight suppression of ${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$ in firehose/mirror-stable regions (figure 23), the two are nearly indistinguishable. This conclusion is broadly consistent with the findings of Santos-Lima et al. (Reference Santos-Lima, de Gouveia Dal Pino, Kowal, Falceta-Gonçalves, Lazarian and Nakwacki2014), who used the CGL equations with limiters. Regions of the plasma lying near or beyond the mirror-instability threshold are subject to a magnetic tension that is effectively enhanced by a factor of $3/2$ and thus are qualitatively unaffected by the positive, limited viscous stress. The regions of the plasma lying near or beyond the firehose instability threshold, at which the effective tension from the Maxwell and viscous stresses is zero, appear to be unimportant to the dynamics of the dynamo, most likely because these regions have a smaller volume-filling factor and contain comparatively lower magnetic energy.(iii) When the dynamical feedback of parallel-viscous dissipation on the field-stretching motions is allowed (the unlimited Braginskii-MHD model), the dynamo takes on a different character. Not only is it slower, but many characteristics of the flow and magnetic field change very little from the kinematic stage to the saturated state. Further, most of these characteristics bear a striking resemblance to those found in the saturated, rather than kinematic, state of the
${Re}\gg {1}$, ${Pm}\gtrsim {1}$ isotropic-MHD dynamo. These include: the magnetic-energy spectrum (figure 9), the characteristic scales of the folded magnetic-field geometry (figure 12), the p.d.f. of the magnetic-field-line curvature (figure 13), and the p.d.f. of the alignment angles between the magnetic field and the rate-of-strain tensor eigenvectors (figures 16 and 17). In addition, a scale-dependent anisotropy of the velocity field was found in the saturated state of the limited Braginskii-MHD dynamo and throughout the entire evolution of the unlimited Braginskii-MHD dynamo (specifically, $\ell _\parallel \propto \ell ^{3/4}_\perp$; see figures 5 and 11); a similar anisotropy was found in the saturated state of the isotropic-MHD dynamo.(iv) Motivated by this resemblance between the saturated MHD state and ‘kinematic’ unlimited Braginskii-MHD dynamos and by the structural similarity of the magnetic-tension and Braginskii-viscous stresses, we have constructed a theory for the unlimited Braginskii-MHD dynamo analogous to that developed by Schekochihin et al. (Reference Schekochihin, Cowley, Taylor, Hammett, Maron and McWilliams2004a) to model the saturated sate of the
${Pm}\gg {1}$ MHD small-scale dynamo (§ 4.3). This theory introduces a field-biased rate-of-strain tensor into the Kazantsev–Kraichnan model of the fluctuation dynamo that captures the partial two-dimensionalization of the velocity-gradient statistics with respect to the local magnetic-field direction caused by the anisotropic viscosity. (In Schekochihin et al. (Reference Schekochihin, Cowley, Taylor, Hammett, Maron and McWilliams2004a), this partial two-dimensionalization is instead caused by the dynamically important Lorentz force in the saturated state.) This model predicts magnetic-energy spectra in remarkably good quantitative agreement with those found in our simulations (arguably, more quantitative an agreement than we had the right to expect).(v) Another prediction of our modified Kazantsev–Kraichnan model is that enhanced small-scale mixing and local two-dimensionalization of the flow shut down the Braginskii-MHD dynamo as the isotropic viscosity
$\mu \rightarrow 0$ (at fixed Braginskii viscosity). We have confirmed this behaviour numerically in the complementary limit of $\mu _{{\rm B}}\rightarrow \infty$ at constant $\mu$ – the ‘Braginskii Stokes-flow’ regime – where we see that the unlimited Braginskii-MHD dynamo fails at sufficiently small parallel Reynolds numbers (at fixed $u_\textrm {rms}$). Isotropic viscosity comparable to the Braginskii viscosity saves the dynamo in this situation by diffusively bleeding momentum into the field-parallel direction; indeed, we find that the Stokes-flow dynamo works in isotropic MHD regardless of the value of ${Re}$ so long as there is sufficient scale separation between the forcing and resistive scales. Combined with conclusion (iv) above, this implies that the Braginskii-MHD dynamo is only viable for moderate values of the ratio ${Re}_\parallel /{Re}$: too small a value results in too much field-line mixing by the perpendicular flows and gradients relative to field-line stretching and thus to excessive resistive dissipation of the magnetic field. This principle is quantified in our modified Kazantsev–Kraichnan model by the quantities $\sigma _\parallel$ and $\sigma _\perp$ (see § 4.3 and, in particular, (4.8)).
Figure 23. Probability distribution function of the parallel rate of strain (left) and the energy-weighted parallel rate of strain (right) for unlimited Braginskii-MHD (purple line; run U1) and limited Braginskii-MHD employing limiters on the firehose and mirror side (green line; run L1) and on the mirror side only (red line; run L1m). Also included for reference are three simulations of isotropic MHD (blue/orange/ yellow lines; runs MHD1/2/4). The p.d.f.s are shown after six correlation times of the random forcing, corresponding to the diffusion-free regime where the resistive scale has not been reached (top) and after fifty correlation times, corresponding to the end of the exponential-growth phase and onset of saturation (bottom). Inset: p.d.f. of the parallel rate of strain for a limited Braginskii-MHD simulation, conditioned on the stability criteria given by (1.6) at 
$t/{t_{\textrm {corr},f}}=50$.
Going forward from here, three avenues for progress are open. First, there are striking similarities in the morphology and statistics of the flows and fields and the evolution of the pressure-anisotropy p.d.f. between the hybrid-kinetic runs reported by Paper I and the unlimited-Braginskii runs presented herein. This suggests that some version of the unlimited Braginskii-MHD equations might be used to describe the ‘magnetized, kinetic regime’ of the plasma dynamo (i.e. one in which the effective collision frequency needed to constrain the pressure anisotropy within the firehose/mirror stability thresholds 
$\nu _\textrm {eff} \sim \beta |{\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}| \gtrsim \varOmega _{\rm i}$). Secondly, it should be assessed how accurate are the pressure-anisotropy limiters for describing the dynamo evolution in the ‘magnetized, fluid regime’ (
$\nu _\textrm {eff}\lesssim \varOmega _{\rm i}$), in which the pressure anisotropy can (in principle) be well regulated by the instabilities. If this regime is realized and the scalings conjectured in § 1.1.2 hold, the interpolation between these two regimes implies a parallel Reynolds number increasing with magnetic-field strength and thus the possibility of an explosive dynamo. This idea will be the subject of a separate publication.
Finally, we have not been able to assess definitively whether Braginskii viscosity breaks the tendency for the small-scale dynamo to saturate with the majority of its magnetic energy residing near resistive scales. On the one hand, it is interesting that the magnetic spectrum found in the unlimited Braginskii-MHD dynamo shows little tendency to concentrate power on resistive scales, and even less tendency to evolve in time in going from the kinematic stage to the saturation (see figure 9). On the other hand, because computational expense limits the maximum achievable scale separation between the viscous and resistive scales in our simulations, we cannot yet establish whether the peak of the spectrum is independent of 
${Rm}$, or whether the conjecture by Malyshkin & Kulsrud (Reference Malyshkin and Kulsrud2002) – that the interchange motions that are undamped by Braginskii viscosity unwrap the folded magnetic fields and thus promote their inverse cascade to larger scales – can be realized. Future numerical work, both fluid and kinetic, should maximize scale separation with the goal of definitively evaluating the ability of the fluctuation dynamo to generate saturated magnetic fields with large-scale coherence in weakly collisional plasmas such as the intracluster medium.
Acknowledgements
Support for D.A.S. and M.W.K. was provided by US DOE Contract DE-AC02-09-CH11466 and an Alfred P. Sloan Research Fellowship in Physics. Support for J.S. was provided by Rutherford Discovery Fellowship RDF-U001804 and Marsden Fund grant UOO1727, which are managed through the Royal Society Te Apārangi. A.A.S. was supported in part by the UK STFC Consolidated grant ST/N000919/1 and by UK EPSRC grants EP/M022331/1 and EP/R034737/1. The authors thank F. Rincon, P. Kempski, E. Quataert and V. Skoutnev for valuable conversations; the anonymous referees for the time and care they spent in reading our manuscript; and S. Cowley, who with J. Maron initiated a study of this problem in 2001 (Maron & Cowley Reference Maron and Cowley2001), for his insight and encouragement. A.A.S. thanks M. Strumik, whose (unpublished) investigations of the dynamical effects of the pressure-anisotropy stress informed some of our early thinking on this problem. The completion of this work was facilitated by the generous hospitality of the Wolfgang Pauli Institute in Vienna and of the Kavli Institute for Theoretical Physics in Santa Barbara (supported by the National Science Foundation under Grant No. NSF PHY-1748958).
Editor Steve Tobias thanks the referees for their advice in evaluating this article.
Declaration of interests
The authors report no conflict of interest.
Appendix A. Linear stability analysis of the unlimited Braginskii-MHD equations
In this appendix, we present the linear theory for the unlimited, incompressible, Braginskii-MHD system of equations (2.1). We demonstrate that, while these equations exhibit a parallel firehose instability, they do not correctly capture the mirror instability.
We assume an unforced, homogeneous, stationary equilibrium state with a uniform pressure anisotropy 
${\rm \Delta} p \ne 0$ and 
${\boldsymbol {B}}_0 = B_0\boldsymbol {\hat {z}}$, subject to small-amplitude perturbations (subscripted with a ‘1’) whose amplitudes are of the form 
$\exp (\gamma t + \textrm {i}{\boldsymbol {k}}\,{\boldsymbol {\cdot }}\,{\boldsymbol {x}})$. While this equilibrium is formally incompatible with the Braginskii closure 
${\rm \Delta} p = 3\mu _{{\rm B}} {\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$, we consider perturbations with sufficiently high frequencies (
$\omega \gg | {\boldsymbol \nabla } {\boldsymbol {u}}|$) and small scales (
$k \gg k_\mu$) that 
${\rm \Delta} p$ may be considered constant in space and time (a similar analysis using a 
${\rm \Delta} p$ obtained from the drift kinetic equation was carried out by Schekochihin et al. (Reference Schekochihin, Cowley, Kulsrud, Hammett and Sharma2005)).
Without loss of generality, let 
${\boldsymbol {k}} = k_\perp \boldsymbol {\hat {x}} + k_\parallel \boldsymbol {\hat {z}}$. The equations, written to first order in the perturbation amplitudes and assuming Laplacian diffusion, are
\begin{gather} ( \gamma + \mu k^2 ) {\boldsymbol{u}}_1 = - \textrm{i}{\boldsymbol{k}} p_1 + \textrm{i} k_\parallel B_0^2 {\boldsymbol{b}}_1 + \textrm{i} k_\parallel {\rm \Delta} p \left({\boldsymbol{b}}_1 - 2b_{1z}\boldsymbol{\hat{z}}\right) - \mu_{{\rm B}} k^2_\parallel u_{1z}\boldsymbol{\hat{z}}, \end{gather}
\begin{gather}{\boldsymbol{b}}_1 = \frac{\textrm{i} k_\parallel }{\gamma + k^2 \eta} {\boldsymbol{u}}_1, \end{gather}
where 
${\boldsymbol {b}}_1 \doteq {\boldsymbol {B}}_1/B_0$. If one were to consider other types of diffusion, 
$\mu$ and 
$\eta$ may be simply redefined (e.g. 
$\mu \rightarrow \mu _{{\rm H}}k^2$ for hyper-diffusion). Incompressibility (
${\boldsymbol {k}}\,{\boldsymbol {\cdot }}\,{\boldsymbol {u}}_1$ = 0) is used to determine 
$p_1$:
\begin{equation} p_1 = ( \textrm{i} k_\parallel \mu_{{\rm B}} u_{1z} - 2b_{1z}{\rm \Delta} p)\frac{k_\parallel^2}{k^2} , \end{equation}
where we have used the solenoidality constraint on the perturbed magnetic field, viz., 
${\boldsymbol {k}}\,{\boldsymbol {\cdot }}\,{\boldsymbol {b}}_1=0$. Substituting (A 3) into (A 1) leads to
\begin{equation} (\gamma + \mu k^2) {\boldsymbol{u}}_1 = \textrm{i} k_\parallel B_0^2 {\boldsymbol{b}}_1 + \textrm{i} k_\parallel {\rm \Delta} p\left[{\boldsymbol{b}}_1 - 2b_{1z}\boldsymbol{\hat{z}}\,{\boldsymbol{\cdot}}\,(\boldsymbol{\mathsf {I}}-{\boldsymbol{\hat{k}\hat{k}}})\right] - \mu_{{\rm B}} k_\parallel^2 u_{1z}\boldsymbol{\hat{z}} \,{\boldsymbol{\cdot}}\,(\boldsymbol{\mathsf {I}}-{\boldsymbol{\hat{k}\hat{k}}}). \end{equation}Combining (A 2) and (A 4) yields the matrix equation
\begin{align} &\left( \gamma + \eta k^2 \right) \left[ (\gamma + \mu k^2 )\boldsymbol{\mathsf {I}} + \mu_{{\rm B}}k_\parallel^2 (\boldsymbol{\mathsf {I}}-{\boldsymbol{\hat{k}\hat{k}}}) \,{\boldsymbol{\cdot}}\,\boldsymbol{\hat{z}}\boldsymbol{\hat{z}} \right]\boldsymbol{\cdot} {\boldsymbol{b}}_1 \nonumber\\ &\quad = - k_\parallel^2 B_0^2 {\boldsymbol{b}}_1 - k_\parallel^2 {\rm \Delta} p \left[\boldsymbol{\mathsf {I}} - 2(\boldsymbol{\mathsf {I}}-{\boldsymbol{\hat{k}\hat{k}}})\,{\boldsymbol{\cdot}}\, \boldsymbol{\hat{z}}\boldsymbol{\hat{z}}\right]\boldsymbol{\cdot} {\boldsymbol{b}}_1. \end{align}
Solenoidality requires that 
$b_{1z} = -(k_\perp /k_\parallel ) b_{1x}$, splitting the dispersion relation into two separate branches,
\begin{gather} \left(\gamma + \eta k^2 \right) \left(\gamma + \mu k^2\right) = - k_\parallel^2 \left(B_0^2 + {\rm \Delta} p \right) , \end{gather}
\begin{gather}\left(\gamma + \eta k^2 \right) \left(\gamma + \mu k^2 + \mu_{{\rm B}} \frac{k_\parallel^2k_\perp^2}{k^2} \right) = - k_\parallel^2 \left(B_0^2 + {\rm \Delta} p \frac{k_\parallel^2 -k_ \perp^2}{k^2} \right). \end{gather}
The first branch, corresponding to (A 6), involves magnetic-field perturbations that are perpendicular to the background magnetic field, with 
${\boldsymbol {\hat {b}}}_1 \doteq {\boldsymbol {b}}_1/|{\boldsymbol {b}}_1| = (1,0,0)$ or 
${\boldsymbol {\hat {b}}}_1 = (0,1,0)$ (i.e. 
$\delta {\boldsymbol {B}}_\perp$). The second branch, corresponding to (A 7), involves both transverse and longitudinal (i.e. 
$\delta B_\parallel$) magnetic-field fluctuations, so that 
${\boldsymbol {\hat {b}}}_1 = (-k_\parallel /k_\perp , 0 ,1)$.
The solution of (A 6) for the first branch is
\begin{equation} \gamma = -D_+ \pm \sqrt{ D_-^2 - k_\parallel^2(B_0^2+{\rm \Delta} p)}, \end{equation}
where 
$D_\pm \doteq (\mu \pm \eta ) k^2/2$. The dissipationless limit 
$\mu = \eta =0$ returns the firehose stability threshold 
${\rm \Delta} p > -B_0^2$ (see (1.6)). Instability occurs once the destabilizing pressure anisotropy nullifies the otherwise restoring magnetic-tension force that is responsible for the shear-Alfvén wave.
Equation (A 7) for the second branch is most clearly analysed in the limit of 
$\mu _{{\rm B}} k_\parallel ^2 \gg \gamma$ and 
$\mu = \eta = 0$, which gives
\begin{equation} \gamma \approx \frac{1}{\mu_{{\rm B}}} \frac{k^2}{k_\perp^2} \left( - B_0^2 + {\rm \Delta} p\frac{k_\perp^2-k_\parallel^2}{k^2} \right). \end{equation}
The stability criterion is 
${\rm \Delta} p(1-2k_\parallel ^2/k^2) < B_0^2$, beyond which the mode exhibits a positive but scale-independent growth rate in the limit 
$k_\parallel / k \ll 1$, in which the magnetic-field perturbation is dominated by 
$\delta B_\parallel$. This illustrates the inability of the Braginskii-MHD system to capture correctly the mirror instability, which in a kinetic calculation exhibits a growth rate that increases as 
$|k_\parallel |$ (at least until finite ion-Larmor-radius effects intercede and reduce growth) and has a stability threshold given by 
${\rm \Delta} p < B_0^2/2$ (Hasegawa Reference Hasegawa1969). Instead, the Braginskii mirror instability grows at a rate proportional to 
${\rm \Delta} p/\mu _{{\rm B}}$, i.e. comparable to the rate of strain of the field-stretching motions that are responsible for driving the pressure anisotropy in the first place (which makes its calculation above formally invalid, since constancy of the background rate of strain was assumed in it).
Appendix B. Magneto-immutability in unlimited Braginskii-MHD dynamo
The minimization of 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$ seen in our unlimited Braginskii-MHD simulations is reminiscent of what was found in the recent study of guide-field Braginskii-MHD turbulence by Squire et al. (Reference Squire, Schekochihin, Quataert and Kunz2019), in which the flow self-organizes to minimize changes in magnetic-field strength and the consequent production of pressure anisotropy – an effect named ‘magneto-immutability’ by those authors. Given the discussion in the main text concerning the anisotropization of the velocity field by the unlimited Braginskii viscosity, it is worth asking whether there is any physical relationship between the results in § 4 and those presented by Squire et al. (Reference Squire, Schekochihin, Quataert and Kunz2019). We also ask whether certain minor differences between the limited-Braginskii-MHD dynamo and its isotropic-MHD counterpart in their saturated states may be attributable to magneto-immutability effects.
B.1 Statistics of ${{\hat {b}}}{{\hat {b}}}\, {:}\, {\nabla } {u}$ and ${{BB}}:{\nabla } {u}$
To address these questions, we follow Squire et al. (Reference Squire, Schekochihin, Quataert and Kunz2019) in calculating the p.d.f.s of both the parallel rate of strain 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$ and the energy-weighted parallel rate of strain 
$\boldsymbol {BB}\!\boldsymbol {:}\!\boldsymbol {\nabla u}/B^2_\mathrm {rms}$. The former quantity gives insight into how magneto-immutability affects the statistics of the parallel rate of strain, while the latter can be used to see how magneto-immutability affects the growth of the magnetic field on average. The results are displayed in figure 23, whose upper panels show these p.d.f.s after six correlation times of the random forcing (corresponding to diffusion-free stage; see § 3.1) and the bottom panels show the same p.d.f.s after 
$50$ correlation times (corresponding to the end of the exponential-growth stage and the onset of saturation).
At the start of the kinematic stage, the p.d.f. from the limited Braginskii-MHD run is almost identical to that from the 
${Pm} = 1$, 
${Re} \gg 1$ MHD run (MHD4). By contrast, the unlimited Braginskii-MHD dynamo resembles more closely the 
${Pm}\gg {1}$ MHD dynamo (MHD1), suggesting a strong tendency towards a magneto-immutable state in which 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$ is dynamically regulated (although the former's p.d.f. is slightly broader, because Braginskii viscosity is overall less efficient at damping motions than is an isotropic viscosity of similar magnitude.) As the limited Braginskii-MHD dynamo evolves (green line in figure 23), the p.d.f. of the parallel rate of strain develops a large notch centred around 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u} = 0$, indicating a relative preference in Braginskii MHD for motions with 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u} \approx 0$. This is consistent with the finding that even systems with hard-wall pressure-anisotropy limiters show a tendency towards magneto-immutability (Squire et al. Reference Squire, Schekochihin, Quataert and Kunz2019). The notch in the p.d.f., clearly seen at 
$t/{t_{\textrm {corr},f}} = 50$ in the green and red curves, can be attributed to the unlimited Braginskii viscosity acting on the mirror- and firehose-stable regions of the plasma. This becomes clear if one conditions the p.d.f. to examine regions that lie either within or outside of the stability region as defined by (1.6); this conditional p.d.f. is displayed in the inset at the centre of figure 23. As portions of plasma enter the region of stability (1.6), unlimited viscous forces quickly damp the parallel rate of strain, condensing the p.d.f. in that region near zero and forming a sharp peak that contrasts with the otherwise wide p.d.f. in the unstable regions where the Braginskii viscosity is greatly reduced. This has the effect of rendering the p.d.f. of 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$ in the limited Braginskii simulation near saturation much thinner than its MHD counterpart (cf. green and blue lines in the bottom panels).
The simulations of Braginskii-MHD dynamo performed here, in conjunction with the p.d.f.s displayed in figure 23, can also be used to validate certain conjectures on the behaviour of the dynamo.
First, it was envisioned by Melville et al. (Reference Melville, Schekochihin and Kunz2016) that microscale instabilities could produce an environment more favourable to motions that amplify the magnetic field than to those that damp it, which would result in a dynamo potentially faster than its isotropic-MHD counterpart. This proposition is tested here in the extreme limit where hard-wall pressure-anisotropy limiters are employed only on the mirror side (run L1m), resulting in a situation where only motions that decrease the magnetic-field strength are strongly damped. One might expect that this case would lead to the fastest possible fluctuation dynamo for a given set of diffusion coefficients. Indeed, the left panels of figure 23 show that negative rates of strain are suppressed compared to the fully limited simulation (compare the red and green lines), suggesting the potential for faster dynamo growth in the case without a firehose limiter. However, these negative rates of strain only amount to modest changes to low-probability regions in the p.d.f. of the energy-weighted rate of strain, whose average drives magnetic-energy growth. The evolution of the box-averaged magnetic energy in this run is not substantially different than in those employing limiters on both regions of instability (see figure 2).
Secondly, it was argued by Malyshkin & Kulsrud (Reference Malyshkin and Kulsrud2002) that, for an isotropically tangled magnetic field on sub-viscous scales, the anisotropic viscous stress on the average resembles an isotropic viscous stress with effective viscosity 
$\mu _{{\rm B}}/5$. To test this hypothesis, the p.d.f.s of 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$ and 
$\boldsymbol {BB}\!\boldsymbol {:}\!\boldsymbol {\nabla u}/B^2_\mathrm {rms}$ from an isotropic-MHD simulation (run MHD2, 
$\mu =1/100$) are included in figure 23 as the orange lines. These p.d.f.s resemble closely those found in the unlimited Braginskii-MHD case with 
$\mu _{{\rm B}}=1/20$ (purple line), as argued. However, the magnetic-energy growth rate is greatly reduced in the unlimited Braginskii system (compare the blue dashed and purple solid lines in the left side of figure 8). As explained in § 3.3.2, this is because sub-parallel-viscous mixing motions allowed by the Braginskii viscosity promote resistive dissipation of the magnetic field, a fact not considered in the work of Malyshkin & Kulsrud (Reference Malyshkin and Kulsrud2002). Thus, while an effective isotropic viscosity of 
$1/5$ its anisotropic counterpart may indeed give an idea of the overall levels of 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$, its use does not accurately capture all the features of the Braginskii-MHD dynamo.
B.2 Nullification of ${{\hat {b}}}{{\hat {b}}}\, {:}\, {\nabla } {u}$ by small-scale modes
Further evidence of magneto-immutability can be seen in the 
$k$-space transfer of kinetic energy by the Braginskii viscous stress, which is shown in figure 24 for the kinematic stage of runs U1 (purple lines) and U4 (orange lines). While some of the energy extracted from the large-scale motions is damped away – note that 
$\mathcal {T}_k[{\boldsymbol \nabla } \,{\boldsymbol {\cdot }}\, ({\boldsymbol {\hat {b}}} {\boldsymbol {\hat {b}}} {\rm \Delta} p)] < 0$ at large scales – a small portion of it is transferred to small-scale fluctuations.Footnote 20 These small-scale fluctuations attempt to nullify the net positive parallel rate of strain by introducing negative rates of strain, as highlighted in the inset of figure 24. In both runs, these motions attempt to counteract the growth of magnetic energy driven by the large scales. In a system where the collisional relaxation of pressure anisotropy is governed by a constant 
$\mu _{{\rm B}}$ (as in our unlimited runs), the only means of regulating the pressure anisotropy and thereby avoiding strong damping of the large-scale motions is to re-organize the fields and flows to control 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$ dynamically.Footnote 21 In kinetic systems, the pressure anisotropy, and thus the parallel viscous stress, can be reduced instead by anomalously increasing the collision frequency through the pitch-angle scattering of particles off firehose/mirror fluctuations. Indeed, such an enhanced collisionality has been directly measured in hybrid-kinetic simulations of the firehose and mirror instabilities (Kunz et al. Reference Kunz, Schekochihin and Stone2014; Melville et al. Reference Melville, Schekochihin and Kunz2016), and is what underpins the very idea of hall-wall limiters, in which 
$\mu _{{\rm B}}$ is effectively modified to maintain a kinetically stable plasma. No such collisional regulation can occur in our unlimited Braginskii-MHD runs, and so the motions responsible for driving 
$p_\perp \ne p_\parallel$ must adjust. Interestingly, the high-
$Re$ isotropic MHD simulation in the saturated state also exhibits small-scale structures that attempt to nullify 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$ (see figure 4 and § 3.2), a feature not previously noted in the literature. Apparently, the unlimited Braginskii-MHD dynamo takes similar steps as the saturated MHD dynamo to reduce the overall stretching rate of the magnetic field.
Figure 24. Time-averaged 
$k$-shell-filtered kinetic-energy transfer function of the pressure anisotropy 
$\mathcal {T}_k[{\boldsymbol \nabla } \,{\boldsymbol {\cdot }}\, ({\boldsymbol {\hat {b}}} {\boldsymbol {\hat {b}}} {\rm \Delta} p)]$ (see (2.6)) and the shell-filtered energy-weighted rate of strain 
$\langle \boldsymbol {BB\!:\!\nabla u}^{[k]} \rangle /B_\mathrm {rms}^2$ (see (2.8)) for runs U1 (purple) and U4 (orange). Bins have size 
$\sqrt {2}k$, where 
$k$ is the centre of the bin. These quantities have been averaged over 60 correlation times during the kinematic stage of the unlimited Braginskii-MHD dynamo. Inset: enlargement of the region enclosed by the dotted rectangle. Note that regions of kinetic-energy injection coincide with regions of negative 
${\boldsymbol {BB\!:\!{\boldsymbol \nabla } u}}$; small-scale fluctuations attempt to nullify the net-positive parallel rate of strain and the consequent growth of magnetic energy driven by the large scales.
B.3 Regulation of pressure anisotropy
Signatures of this adjustment of 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}$ can be seen in figure 25, which displays the p.d.f. of the parallel rate of strain with respect to the reciprocal of the magnetic energy 
$B^{-2}$ (which serves as a proxy for the plasma 
$\beta$ seen in similar histograms of pressure anisotropy; see, e.g. Bale et al. Reference Bale, Kasper, Howes, Quataert, Salem and Sundkvist2009 and Reference St-Onge and KunzPaper I). Dash-dotted (dashed) lines trace the mirror (firehose) instability thresholds given by (1.6), with the stable region lying to the left of these lines. In the limited case, the bulk of the plasma stays beyond the mirror-instability threshold, also migrating toward larger magnetic energies as the dynamo progresses. Because the viscous stress is limited, its ability to damp out the motions driving the p.d.f. in the direction of mirror instability is attenuated. As a result, the portions of the plasma lying beyond the stability boundaries are subject to a greatly reduced viscous stress that is comparable to the (small) Maxwell stress. However, as the magnetic energy grows, the limited viscous stress is able to regain its dynamical importance. More and more plasma exhibits pressure anisotropy in the stable region (as also seen in figure 23) and evolves in a way that appears to respect the mirror boundary. No such behaviour is seen in the unlimited run, although the bulk of the plasma does seem to sit on the mirror threshold 
${\rm \Delta} p \approx B^2_0/2$ once saturation is reached.Footnote 22 Instead, there is a dynamical feedback on the motions with 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u} > 0$ that limits their ability to amplify the magnetic field and drive 
$p_\perp > p_\parallel$.
Figure 25. Probability distribution function of the parallel rate of strain with respect to 
$B^{-2}$ for the hard-wall limited run L1 (a) and unlimited Braginskii-MHD run U1 (b) at successive times. Dashed-dotted (dotted) lines trace the mirror (firehose) instability thresholds given by (1.6).
Figure 25 is to be compared to the analogous figure for pressure anisotropy versus plasma 
$\beta$ in the magnetized kinetic regime computed using a collisionless hybrid-kinetic code (figure 2 of Reference St-Onge and KunzPaper I). The behaviour throughout much of the kinetic calculation bears a stronger resemblance to that found in the unlimited Braginskii run than in the limited run. As explained in § 1.1.3, the effective collision frequencies required at such large plasma 
$\beta$ to pin the pressure anisotropy instantaneously at its stability thresholds are difficult to realize in a simulation with limited scale separation. Indeed, the simulations of Reference St-Onge and KunzPaper I have a box-averaged anomalous collision frequency that is much smaller than 
$\beta |{\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}\, {\boldsymbol :}\, {\boldsymbol \nabla } {\boldsymbol u}|$ throughout the kinematic state of evolution, which puts further in doubt whether hard-wall limiters can truly serve as a realistic model for the action of mirror and firehose instabilities in kinetic plasmas.
Appendix C. Derivation of anisotropic Kazantsev–Kraichnan model
For completeness, we provide here the full derivation of (4.7), the anisotropic Kazantsev–Kraichnan model presented and utilized in § 4.3. This derivation was omitted in Schekochihin et al. (Reference Schekochihin, Cowley, Taylor, Hammett, Maron and McWilliams2004a) for lack of space and never came to be published anywhere, but we feel that it is important to spell out all the steps and assumptions that went into it in order for the reader to be able to judge the level of plausibility of our arguments in § 4.3.
We start with the evolution equations
\begin{gather} \partial_t \tilde{B} = \tilde{{b}}^i\tilde{{b}}^m \partial_m \tilde{u}^i \tilde{B} - \eta \tilde{k}^2 \tilde{B}, \end{gather}
\begin{gather}\partial_t \tilde{k}_m = - \partial_m \tilde{u}^i \tilde{k}_{i}, \end{gather}
\begin{gather}\partial_t \tilde{{b}}^i = \tilde{{b}}^m \partial_m \tilde{u}^i - \tilde{{b}}^l \tilde{{b}}^m \partial_m \tilde{u}^l \tilde{{b}}^i, \end{gather}
where 
${b}^i$ denotes the 
$i$th component of the magnetic-field unit vector 
${\boldsymbol {\hat {b}}}$ and tildes are used in this section to denote random variables (cf. (4.4)). We assume that 
$\tilde {{\boldsymbol {u}}}$ is approximately linear in space,Footnote 23 white in time, and anisotropic with respect to the local magnetic-field direction:
\begin{equation} \tilde{u}^i(t,{\boldsymbol{x}}) = \tilde{\sigma}^i_m (t) x^m, \end{equation}
where the rate-of-strain tensor 
$\tilde {\sigma }^i_m \doteq \partial _m \tilde {u}^i$ has the two-time correlation function
\begin{equation} \overline{\tilde{\sigma}^i_m(t)\tilde{\sigma}^j_n(t')} = \tilde{\varGamma}^{ij}_{mn}(t)\delta(t-t') , \end{equation}where an overline denotes the ensemble average and
\begin{align} \tilde{\varGamma}^{ij}_{mn}(t) &= \kappa_2 \left[\delta^{ij}\delta_{mn} + a(\delta^i_m \delta^j_n + \delta^i_n \delta^j_m) + \chi_1 \delta^{ij}\tilde{{b}}^m\tilde{{b}}^n + \chi_2 \tilde{{b}}^i\tilde{{b}}^j \delta_{mn} + \chi_3 \tilde{{b}}^i\tilde{{b}}^j \tilde{{b}}^m \tilde{{b}}^n \right.\nonumber\\ &\left.\quad +\,\chi_4(\delta^i_m \tilde{{b}}^n \tilde{{b}}^j + \delta^i_n \tilde{{b}}^m \tilde{{b}}^j + \tilde{{b}}^i\tilde{{b}}^m \delta^j_n + \tilde{{b}}^i\tilde{{b}}^n \delta^j_m)\right] \end{align}
is the general fourth-rank tensor that is anisotropic with respect to the magnetic-field direction and symmetric under interchange of 
$i$, 
$m$ and 
$j$, 
$n$. Here, 
$\kappa _2$ is the second-order coefficient in the Taylor expansion of the field-anisotropic velocity correlation tensor defined by (4.1)
\begin{align} \tilde{\kappa}^{ij}({\boldsymbol{y}}) &= \kappa_0 \delta^{ij} - \frac{1}{2} \kappa_2 \left[ y^2 \delta^{ij} + 2a y^i y^j + \chi_1 \delta^{ij} (y_m\tilde{{b}}^m)^2 + \chi_2 \tilde{{b}}^i \tilde{{b}}^j y^2 + \chi_3 \tilde{{b}}^i \tilde{{b}}^j (y_m\tilde{{b}}^m)^2\right. \nonumber\\ &\quad \left.+\, 2\chi_4 (y_m\tilde{{b}}^m) \left( \tilde{{b}}^j y^i + \tilde{{b}}^i y^j \right) \right] + \dots , \end{align}
and the constants 
$\chi _i$ and 
$a$ parametrize the rate of strain. Two of them can be fixed by assuming an incompressible flow: 
$\tilde {\varGamma }^{ij}_{in}=0$ and 
$\tilde {\varGamma }^{ij}_{mj} = 0$, so
\begin{equation} a = -\frac{1+\chi_4}{1+d} , \quad \chi_1+\chi_2+\chi_3=-(d+2)\chi_4, \end{equation}
where 
$d$ is the dimensionality of the system. The isotropic case is recovered when 
$\chi _i = 0$ for all 
$i$. The tensors 
$\tilde {\varGamma }^{ij}_{mn}$ and 
$\tilde {\kappa }^{ij}({\boldsymbol {k}})$ (see (4.2)) are related by
\begin{equation} \tilde{\varGamma}^{ij}_{mn}(t) = \int \frac{\textrm{d}^d {\boldsymbol{k}}}{(2{\rm \pi})^d} \, k_mk_n \tilde{\kappa}^{ij}({\boldsymbol{k}}). \end{equation}
Formally, 
$\tilde {\varGamma }^{ij}_{mn}(t)$ depends on the random quantity 
$\tilde {{b}}_i$, and so the averaging performed in (C 3) is understood to be a time average taken over times longer than the velocity correlation time (which is zero), but shorter than the time over which 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}$ changes significantly. This is a reasonable simplification if 
${\boldsymbol {\hat {b}}}$ aligns itself along the direction of the maximal Lyapunov exponent of the flow, which occurs exponentially fast in time (Goldhirsch, Sulem & Orszag Reference Goldhirsch, Sulem and Orszag1987). Therefore, in the (Lagrangian) frame of the flow, the vector 
${\boldsymbol {\hat {b}}}$ changes slowly in time.
We aim to obtain an evolution equation for the magnetic-energy spectrum 
$M(k)$. This can be done by first deriving an evolution equation for the joint probability density function of the magnetic field 
${\boldsymbol {B}}$ and its wavenumber 
${\boldsymbol {k}}$,
\begin{equation} \mathcal{P}(t;B,{\boldsymbol{k}},{\boldsymbol{\hat{b}}}) = \overline{\tilde{\mathcal{P}}} \doteq \overline{\delta(B-\tilde{B}(t))\delta(k_m - \tilde{k}_m(t))\delta({b}^i - \tilde{{b}}^i(t))}, \end{equation}and then taking the appropriate moments; to wit,
\begin{equation} M(t,k) = \int \textrm{d} \varOmega_{\boldsymbol{k}}\, k^2 \int \textrm{d}^3 {\boldsymbol{\hat{b}}} \int\textrm{d} B\, B^2 \mathcal{P}(t; B,{\boldsymbol{k}},{\boldsymbol{\hat{b}}}), \end{equation}
where 
$\textrm {d} \varOmega _{\boldsymbol {k}}$ is an element of solid angle in wavenumber space. Taking the time derivative of (C 8) and using (C.1), we get
\begin{align} \frac{\partial \mathcal{P}}{\partial t} &= \overline{\frac{\partial \tilde{\mathcal{P}}}{\partial t}} = \overline{\left[-\frac{\partial \tilde{B}(t)}{\partial t} \frac{\partial }{\partial B} - \frac{\partial \tilde{k}_m(t)}{\partial t} \frac{\partial }{\partial k_m} - \frac{\partial \tilde{{b}}^i(t)}{\partial t} \frac{\partial }{\partial {b}^i} \right]\tilde{\mathcal{P}}} \nonumber\\ &= -\overline{\left[ \tilde{{b}}^i(t) \tilde{{b}}^m(t) \tilde{\sigma}^i_m(t) \tilde{B}(t) \frac{\partial }{\partial B} - \eta \tilde{k}^2(t)\tilde{B}(t) \frac{\partial }{\partial B} - \tilde{\sigma}^i_m (t)\tilde{k}_i(t) \frac{\partial }{\partial k_m}\right.} \nonumber\\ & \quad\overline{\left.\mbox{} + \tilde{\sigma}^i_m(t) \tilde{{b}}^m(t) \frac{\partial }{\partial {b}^i} - \tilde{\sigma}^l_m(t) \tilde{{b}}^l(t) \tilde{{b}}^m(t) \tilde{{b}}^i(t) \frac{\partial }{\partial {b}^i}\right]\tilde{\mathcal{P}}} \nonumber\\ &= - \hat{L}_i^m \overline{\tilde{\sigma}^i_m(t) \tilde{\mathcal{P}}} + \eta k^2 \frac{\partial }{\partial B}B \mathcal{P}, \end{align}where
\begin{equation} \hat{L}^m_i \doteq \frac{\partial }{\partial B} B {b}^i {b}^m - \frac{\partial }{\partial k_m} k_i + \frac{\partial }{\partial {b}^i} b^m - \frac{\partial }{\partial {b}^l} {b}^l {b}^i {b}^m. \end{equation}
To arrive at the final line of (C 10), the identity 
$a\delta (a-b) = b\delta (a-b)$ was used. Note that everything in the square brackets in the final line is non-random. To perform the remaining ensemble average, we make use of the Furutsu–Novikov formula (Furutsu Reference Furutsu1963; Novikov Reference Novikov1965), which generalizes Gaussian splitting to functions,
\begin{equation} \overline{\tilde{\sigma}_m^i(t) \tilde{\mathcal{P}}(t)} = \int_0^t \textrm{d} t'\, \overline{\tilde{\sigma}^i_m(t) \tilde{\sigma}^j_n(t')} \, \overline{\frac{\delta \tilde{\mathcal{P}}(t) }{\delta \tilde{\sigma}^j_n(t')}} . \end{equation}
The functional derivative with respect to 
$\tilde {\sigma }^j_n(t')$ can be calculated by formally integrating 
$\partial \tilde {\mathcal {P}}(t)/\partial t$ with respect to time,
The last term in first line is dropped because it disappears when 
$t=t'$, owing to causality (i.e. 
$\tilde {\mathcal {P}}$ cannot depend on future values of 
$\tilde {\sigma }$). Using this result in (C 12) alongside (C 3), we find
\begin{equation} \overline{\tilde{\sigma}_m^i(t) \tilde{\mathcal{P}}(t)} = -\frac{1}{2} \hat{L}^n_j \tilde{\varGamma}_{mn}^{ij}(t) \overline{\tilde{\mathcal{P}}(t)} \approx-\frac{1}{2} \hat{L}^n_j \overline{\tilde{\varGamma}_{mn}^{ij}(t) \tilde{\mathcal{P}}(t)} = -\frac{1}{2} \hat{L}^n_j \varGamma_{mn}^{ij} \overline{\tilde{\mathcal{P}}(t)} , \end{equation}
where 
$\varGamma ^{ij}_{mn}$ is defined similarly to (C 4) but with the substitution 
$\tilde {{b}} \rightarrow {b}$. The second (approximate) equality relies on the assumption of a slowly varying 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}$, as was used to obtain (C 3). Recall that, in the striped volume-filling region of a magnetic fold, only the sign of 
${\boldsymbol {\hat {b}}}$ changes in space, and so the dyad 
${\boldsymbol {\hat {b}}}{\boldsymbol {\hat {b}}}$, a quadratic quantity, changes on a much longer time scale than either the magnetic energy or the turbulent velocity. Once 
$\tilde {\varGamma }_{mn}^{ij}$ is brought within the average, it becomes 
$\varGamma _{mn}^{ij}$ via the delta functions in the definition of 
$\tilde {\mathcal {P}}$. The result of these manipulations is a closed equation for the joint probability density function
\begin{equation} \frac{\partial \mathcal{P}}{\partial t} = \frac{1}{2} \hat{L}^m_i \hat{L}^n_j {\varGamma}^{ij}_{mn} \mathcal{P} + \eta k^2 \frac{\partial }{\partial B} B\mathcal{P} . \end{equation} Equation (C 15) can be greatly simplified by noting that 
$\mathcal {P}(B,{\boldsymbol {k}},{\boldsymbol {\hat {b}}})$ must have the following factorization:
\begin{equation} \mathcal{P}(B,{\boldsymbol{k}},{\boldsymbol{\hat{b}}}) = \delta(|{\boldsymbol{\hat{b}}}|^2 -1) \delta({\boldsymbol{\hat{b}}}\,{\boldsymbol{\cdot}}\, {\boldsymbol{k}}) \hat{P}(B,k) . \end{equation}
The two delta functions result, respectively, from 
${\boldsymbol {\hat {b}}}$ being a unit vector and from the solenoidality constraint 
${\boldsymbol {\hat {b}}}\,{\boldsymbol {\cdot }}\,{\boldsymbol {k}} = 0$. The remaining factor in (C 16), 
$\hat {P}(B,k)$, is a result of the statistics being homogeneous and the relative alignment of 
${\boldsymbol {\hat {b}}}$ and 
${\boldsymbol {k}}$ being fixed.Footnote 24 In order to express (C 15) in terms of 
$\hat {P}$, we use the chain rule to write
\begin{equation} \hat{L}_j^n = -k_j \frac{\partial}{\partial k_n} + \left({b}^n \frac{\partial}{\partial {b}^j} -{b}^j {b}^n {b}^q \frac{\partial}{\partial {b}^q} \right) + {b}^j {b}^n \left(\frac{\partial}{\partial B}B - d - 2\right), \end{equation}
where 
$d$ is the dimensionality of the system, and then calculate the combination 
$\hat {L}^n_j \varGamma ^{ij}_{mn}\mathcal {P}(B,{\boldsymbol {k}},{\boldsymbol {\hat {b}}})$. For any function 
$f$,
\begin{align} \hat{L}^n_j \delta({\boldsymbol{\hat{b}}}\,{\boldsymbol{\cdot}}\, {\boldsymbol{k}})f &=\delta({\boldsymbol{\hat{b}}}\,{\boldsymbol{\cdot}}\, {\boldsymbol{k}})\hat{L}^n_j f - k_j {b}^n \delta'({\boldsymbol{\hat{b}}}\,{\boldsymbol{\cdot}}\,{\boldsymbol{k}}) f+ ({b}^n k_j - {b}^j {b}^n {\boldsymbol{\hat{b}}} \,{\boldsymbol{\cdot}}\, {\boldsymbol{k}})\delta'({\boldsymbol{\hat{b}}}\,{\boldsymbol{\cdot}}\,{\boldsymbol{k}})f \nonumber\\ &= \delta({\boldsymbol{\hat{b}}}\,{\boldsymbol{\cdot}}\,{\boldsymbol{k}}) \hat{L}_j^n f- {b}^j {b}^n ({\boldsymbol{\hat{b}}}\,{\boldsymbol{\cdot}}\,{\boldsymbol{k}})\delta'({\boldsymbol{\hat{b}}}\,{\boldsymbol{\cdot}}\, {\boldsymbol{k}})f \nonumber\\ &= \delta({\boldsymbol{\hat{b}}} \,{\boldsymbol{\cdot}}\, {\boldsymbol{k}})(\hat{L}^n_j + {b}^n {b}^j)f, \end{align}
where we have used 
$x \delta '(x) = -\delta (x)$ to obtain the final equality. Similarly,
\begin{equation} \hat{L}^n_j \delta(|{\boldsymbol{\hat{b}}}|^2 -1)f = \delta(|{\boldsymbol{\hat{b}}}|^2-1)(\hat{L}^n_j + 2 {b}^n {b}^j)f.\end{equation}Combining (C 18) and (C 19) leads to
\begin{align} \hat{L}^n_j\varGamma^{ij}_{mn}\mathcal{P} &= \delta({\boldsymbol{\hat{b}}}\,{\boldsymbol{\cdot}}\,{\boldsymbol{k}})\delta(|{\boldsymbol{\hat{b}}}|^2-1)(\hat{L}^n_j + 3{b}^n{b}^j)\varGamma^{ij}_{mn}\hat{P}(B,k) \nonumber\\ &= \delta({\boldsymbol{\hat{b}}}\,{\boldsymbol{\cdot}}\,{\boldsymbol{k}})\delta(|{\boldsymbol{\hat{b}}}|^2-1)\left[-\frac{k_jk_n}{k^2} \varGamma_{mn}^{ij}k \frac{\partial }{\partial k} + {b}^j {b}^n \varGamma^{ij}_{mn}\left(\frac{\partial }{\partial B} B - d +1\right) \right.\nonumber\\ &\quad\left. +\, {b}^n \frac{\partial \varGamma^{ij}_{mn}}{\partial {b}^j} - {b}^j{b}^n{b}^q \frac{\partial \varGamma^{ij}_{mn}}{\partial {b}^q} \right]\hat{P}(B,k). \end{align}
Taking into account the solenoidality constraint 
${\boldsymbol {\hat {b}}}\,{\boldsymbol {\cdot }}\, {\boldsymbol {k}} = 0$ imposed by the prefactor in (C 20), the following combinations in (C 20) may be calculated:
\begin{align} \frac{k_jk_n}{k^2}\varGamma^{ij}_{mn} &=\kappa_2 \left[a\delta^i_m + (a+1) \frac{k_ik_m}{k^2} + \chi_4 {b}^i {b}^m \right], \end{align}
\begin{align} {b}^j {b}^n \varGamma^{ij}_{mn} &= \kappa_2 \left[(a+\chi_4)\delta^i_m + (1+a+\chi_1+\chi_2+\chi_3 +3\chi_4){b}^i{b}^m \right],\\ {b}^n \frac{\partial}{\partial {b}^j}\varGamma^{ij}_{mn} &= \kappa_2 \left[\delta^i_m (\chi_1 + (d+2)\chi_4)\right. \nonumber\\ &\quad\left. +\, {b}^i{b}^m(\chi_1+ (d+1)\chi_2 + (d+3)\chi_3 + (d+4)\chi_4) \right], \end{align}
\begin{align} {b}^j {b}^n {b}^q \frac{\partial}{\partial {b}^q} \varGamma^{ij}_{mn} &= \kappa_2 \left[2\chi_4 \delta^i_m + 2(\chi_1+\chi_2 + 2\chi_3 + 3\chi_4){b}^i{b}^m \right]. \end{align}Using these formulae in (C 20) gives
\begin{align} \hat{L}^n_j\varGamma^{ij}_{mn}\mathcal{P} &= \kappa_2 \delta({\boldsymbol{\hat{b}}}\,{\boldsymbol{\cdot}}\,{\boldsymbol{k}})\delta(|{\boldsymbol{\hat{b}}}|^2-1)\left\{ -\left[a\delta^i_m + (a+1)\frac{k_ik_m}{k^2} + \chi_4 {b}^i {b}^m\right] k \frac{\partial}{\partial k} \right.\nonumber\\ &\quad +\left[(a+\chi_4)\delta^i_m + (1+a+\chi_1+\chi_2+\chi_3+3\chi_4){b}^i{b}^m \right]\frac{\partial}{\partial B}B \nonumber\\ &\quad +\left[-a(d-1) + \chi_1 + \chi_4\right]\delta^i_m\nonumber\\ &\quad \vphantom{\left[a\delta^i_m + (a+1)\frac{k_ik_m}{k^2} + \chi_4 {b}^i {b}^m\right]}+ \left.\left[-(1+a)(d-1) -d \chi_1 - (2d-1)\chi_4\right] {b}^i {b}^m \right\}\hat{P}(B,k). \end{align}
Further applying the operator 
$\hat {L}_i^m$ and expending much effort along the same lines gives us the expression for the first term in (C 15)
\begin{align} &\hat{L}_i^m\hat{L}^n_j\varGamma^{ij}_{mn}\mathcal{P}\end{align}
\begin{align} \quad & = \left\{(2a+1)k^2 \frac{\partial^2}{\partial k^2} + (1+2a+\chi_1 + \chi_2 + \chi_3 + 4\chi_4) \frac{\partial}{\partial B}B\frac{\partial}{\partial B}B\right.\nonumber\\ &\left.\qquad -2(a+\chi_4) \frac{\partial}{\partial B}B k \frac{\partial}{\partial k} + \left[d+(3d-1)a-\chi_1 - \chi_4\right]k \frac{\partial}{\partial k}\right.\notag\\ &\qquad \left.- (1+3a+\chi_1 + 3\chi_4)(d-1)\frac{\partial}{\partial B}B + \left[(d-1)a-\chi_1-\chi_4\right](d-1) \right\}\hat{P}(B,k). \end{align}Normalizability of the p.d.f. requires that
\begin{align} 1 &= \int \textrm{d}^d{\boldsymbol{\hat{b}}} \, \delta(|{\boldsymbol{\hat{b}}}|^2 -1) \int \textrm{d}^d {\boldsymbol{k}}\, \delta({\boldsymbol{\hat{b}}}\,{\boldsymbol{\cdot}}\, {\boldsymbol{k}}) \int \textrm{d} B \, \hat{P}(B,k) \nonumber\\ &= \frac{S_{d-1}S_{d-2}}{2}\int_0^\infty \textrm{d} k \,k^{d-2} \int \textrm{d} B \, \hat{P}(B,k), \end{align}
where 
$S_n$ is the surface area of unit 
$n$-sphere (e.g. 
$S_1 = 2{\rm \pi}$, 
$S_2 = 4{\rm \pi}$). Taking this normalization into consideration, we define 
$P(B,k) \doteq S_{d-1}S_{d-2} k^{d-2} \hat {P}(B,k)/2$ and use
\begin{align} k^{d-2}k\frac{\partial}{\partial k} \frac{1}{k^{d-2}}P & = \left[k\frac{\partial}{\partial k} - (d-2)\right]P, \end{align}
\begin{align}k^{d-2}k^2 \frac{\partial^2}{\partial k^2} \frac{1}{k^{d-2}} P & = \left[k^2 \frac{\partial^2}{\partial k^2} - 2(d-2) k \frac{\partial}{\partial k} + (d-2)(d-1)\right]P \end{align}to turn (C 15) into
\begin{align} \frac{\partial P }{\partial t} &= \frac{1}{2}\kappa_2 \left\{(1+2a)\frac{\partial}{\partial k} k^2 \frac{\partial}{\partial k} + ( 1+2a + \chi_1 + \chi_2 + \chi_3 + 4 \chi_4 ) \frac{\partial}{\partial B}B \frac{\partial}{\partial B}B \right. \nonumber\\ &\quad - 2(a+\chi_4) \frac{\partial}{\partial B}B \frac{\partial}{\partial k}k - \left[(d-2) + (d-3)a + \chi_1 + \chi_4\right] \frac{\partial}{\partial k}k \nonumber\\ &\quad\left. -\, (1+a+\chi_1 +\chi_4 )(d-1)\frac{\partial}{\partial B}B\right\} P + \eta k^2 \frac{\partial }{\partial B} B P . \end{align}We now enforce incompressibility. Substituting (C 6a,b) into (C 27) leads to the final form of the evolution equation of the joint p.d.f.,
\begin{align} \frac{\partial P}{\partial t} &= \frac{\kappa_2}{2(d+1)} \left\{(d-1-2\chi_4)\frac{\partial}{\partial k}k^2 \frac{\partial}{\partial k} + (d-1)(1-d\chi_4)\frac{\partial}{\partial B} B \frac{\partial}{\partial B} B \right. \nonumber\\ &\quad +2(1+d\chi_4) \frac{\partial}{\partial B} B \frac{\partial}{\partial k} k - \left[(d-1)^2 +4 \chi_4 + (d+1)\chi_1\right]\frac{\partial}{\partial k}k \nonumber\\ &\quad\left.-\,(d-1)\left[d(1+\chi_4) + (d+1)\chi_1\right] \frac{\partial}{\partial B}B \right\} P + \eta k^2 \frac{\partial }{\partial B} B P. \end{align}
The magnetic-energy spectrum is 
$M(k) = (1/2) \int _0^\infty \textrm {d} B\, B^2 P(B,k)$. Therefore, taking the 
$B^2$ moment of (C 28) leads in three dimensions (
$d=3$) to
\begin{align} \frac{\partial M}{\partial t} = \frac{\kappa_2}{4} \left[ (1-\chi_4) \frac{\partial }{\partial t} k^2 \frac{\partial }{\partial k} - 2(2+\chi_1 -2 \chi_4) \frac{\partial }{\partial k} k + 2(5+4\chi_1-3\chi_4 )\right] M - 2\eta k^2 M. \end{align}
Equation (4.7) for the time evolution of the magnetic spectrum written as a function of 
$\gamma _\perp$, 
$\sigma _\perp$, and 
$\sigma _\parallel$ follows from (C 29) after noting that
\begin{align} \gamma_\perp &\doteq \int \frac{\textrm{d}^3{\boldsymbol{k}}}{(2{\rm \pi})^3} \, k_\perp^2 \kappa_\perp = \frac{3-\chi_4}{2}\kappa_2, \end{align}
\begin{align}\sigma_\perp &\doteq \frac{1}{\gamma_\perp}\int \frac{\textrm{d}^3 \,{\boldsymbol{k}}}{(2{\rm \pi})^3} k_\parallel^2 \kappa_\perp = \frac{2(1+\chi_1)}{3-\chi_4}, \end{align}
\begin{align}\sigma_\parallel &\doteq \frac{1}{\gamma_\perp}\int \frac{\textrm{d}^3{\boldsymbol{k}}}{(2{\rm \pi})^3} \,k_\parallel^2 \kappa_\parallel = \frac{1-3\chi_4}{2(3-\chi_4)}, \end{align}
