The turbulent dynamo

Abstract

The generation of a magnetic field in an electrically conducting fluid generally involves the complicated nonlinear interaction of flow turbulence, rotation and field. This dynamo process is of great importance in geophysics, planetary science and astrophysics, since magnetic fields are known to play a key role in the dynamics of these systems. This paper gives an introduction to dynamo theory for the fluid dynamicist. It proceeds by laying the groundwork, introducing the equations and techniques that are at the heart of dynamo theory, before presenting some simple dynamo solutions. The problems currently exercising dynamo theorists are then introduced, along with the attempts to make progress. The paper concludes with the argument that progress in dynamo theory will be made in the future by utilising and advancing some of the current breakthroughs in neutral fluid turbulence such as those in transition, self-sustaining processes, turbulence/mean-flow interaction, statistical and data-driven methods and maintenance and loss of balance.

1. Introduction

1.1. Dynamo theory for the fluid dynamicist

It is really just a matter of perspective. To the fluid dynamicist, dynamo theory may appear as a rather esoteric and niche branch of fluid mechanics – in dynamo theory much attention has focused on seeking solutions to the induction equation rather than those for the Navier–Stokes equations. Conversely, to a practitioner dynamo theory is a field with myriad subtleties; in a severe interpretation the Navier–Stokes equations and the whole of neutral fluid mechanics may be regarded as forming a useful invariant subspace of the dynamo problem, with – it has to be said – non-trivial dynamics. In this perspective, I shall attempt to present the important and interesting developments in dynamo theory from the point of view of a fluid dynamicist, pointing out common themes. I shall focus on explaining how recent developments in fluid mechanics can contribute to future breakthroughs for magnetohydrodynamic dynamos.

This is not a review in the typical sense. Although I shall present the key results and features of dynamo theory, I shall not be exhaustive by any means. This perspective is focused on those areas of dynamo theory that I believe are both accessible and of interest to fluid dynamicists, drawing analogies with other areas of fluids where necessary. It is also concentrated on those areas that I believe offer the greatest scope for imminent breakthroughs. These are not necessarily those areas of research that lead to the most accurate modelling of any given astrophysical object; although it is undoubtedly the case (as described in the next section) that the generation of magnetic fields in these objects nearly always forms the motivation for dynamo investigations. Further details of dynamo theory are contained in myriad reviews and monographs, for example those of Moffatt (1978), Krause & Raedler (1980), Brandenburg & Subramanian (2005), Jones (2008) and the recently published monograph of Moffatt & Dormy (2019) and review of Rincon (2019).

We begin, however, by giving motivation for the study of dynamos – much of which arises from observations of cosmical magnetic fields, including those of planets, stars, galaxies and disks.

1.2. Motivation

1.2.1. The geomagnetic field

The Earth's magnetic field is presently predominantly dipolar with a mean surface strength of approximately $40\ \mathrm {\mu }$T. Currently the dipole axis is offset by approximately $10^\circ$ from the rotation axis. Paleomagnetic records indicate that the magnetic field has persisted for greater than three billion years and has always had a significant dipole component (see e.g. Dormy, Valet & Courtillot 2000; Aubert, Tarduno & Johnson 2010). Figure 1 shows how the record is punctuated by episodic reversals of the polarity of the dipole component; whilst the magnetic field reverses (in a time of the order of $10^4$ years) the field energy decreases and the field becomes smaller scale and more multipolar. The figure clearly shows the time scale between reversals of the field is exceptionally long, with some events (superchrons) lasting over ten million years. The Earth's magnetic field can currently be measured at the surface to spherical harmonic degree 13, with higher harmonics being screened by remnant magnetism of the crust. In addition to long term variation of the field, it is also subject to ‘secular variation’, which is seen in current and historical observations (Jackson, Jonkers & Walker 2000). Here, the spatial structure and strength of the field varies with time scales ranging from years to centuries, with significant features including the ‘westward drift’ of magnetic flux patches and changes in the length of day. These changes arise because of the interaction of magnetic fields with flows in the electrically conducting fluid outer core of the Earth.

Figure 1. The geomagnetic record showing the history of magnetic field reversals (source Wikipedia). The black/white bars relate to magnetic fields of opposite polarities.

1.2.2. Solar system planets

Most, if not all, solar system planets currently possess or have possessed dynamo-generated magnetic fields. Even smaller satellites such as Ganymede and the Moon show signs of current or historic dynamo action. It seems as though magnetic field generation is possible whether the planet is terrestrial, a gas giant, or an ice giant. It is therefore expected, and there is some observational evidence to support the theory, that many exoplanets should also be capable of dynamo action (Shkolnik et al. 2008).

As is the case for the Earth, dynamo action usually takes place in the interior of planets but the magnetic fields are measured as a potential field having diffused through poorly conducting regions. These potential fields are usually decomposed into spherical harmonics with amplitudes given by the so-called Gauss coefficients (see e.g. Schubert & Soderlund 2011). Briefly, the gas giants Jupiter and Saturn possess strong dipole dominated magnetic fields. Jupiter's field, confirmed by the Pioneer 10 flyby, has a largely axial dipole and, with a mean surface strength of $550\ \mathrm {\mu }$T, has the strongest field in the solar system. Before the recent Juno mission the observations of the magnetic field had a poor resolution (up to spherical harmonic degree three), although recent flybys by the Juno satellite are beginning to establish that the magnetic field has much more structure at smaller scales and a distinct hemispheric asymmetry (see e.g. Connerney 2018). It is worth mentioning at this point that the Juno mission will probably give us the closest direct observation of a naturally occurring dynamo generated magnetic field; it will be fascinating to follow the progress of the mission. Saturn's magnetic field, first revealed by Pioneer 11, has been measured to spherical harmonic degree three and has a mean surface strength of $30\ \mathrm {\mu }$T. The axis of the dipole is remarkably aligned with the rotation axis (with an offset of less than $1^\circ$) meaning that the results are consistent with the field being axisymmetric. As we shall see in § 2.4.3, it is not possible for such a field configuration to be generated by dynamo action (Cowling 1933); this prompts the widely held belief that Saturn's magnetic field exists solely to have annoyed Cowling.

The terrestrial planets – possessing iron-alloy central cores, silicate mantles and rocky crusts – also largely exhibit dynamo generated magnetic fields. Mercury's magnetic field, the weakest in the solar system at $0.3\ \mathrm {\mu }$T is dipole dominated, with the dipole largely aligned with the rotation axis; similarly the magnetic field of the Jovian satellite Ganymede is an axially aligned dipole of about $1\ \mathrm {\mu }$T. Meanwhile the Moon and Mars do not have magnetic fields that are currently maintained by dynamo action, but there is evidence of extinct dynamos in both cases (see Schubert & Soderlund (2011) and the references therein). It is unclear in both cases when, or why, the dynamo switched off, although theories for the cessation of such dynamos usually involve the cooling of the body leading to the unsustainability of thermal convection or the inner core of the body growing to such a size to make dynamo action inefficient. Venus has no measurable intrinsic magnetic field. It may be that the absence of an inner core in Venus means that dynamo action is not possible or that the core is not convective at all. Finally the ice giants Neptune and Uranus have magnetic fields as discovered by the Voyager 2 flybys. These planets possess fundamentally different multipolar magnetic fields, with mean surface field strengths of $30\ \mathrm {\mu }$T, and with the dipole component exhibiting a significant tilt from the rotation axis.

1.2.3. The Sun and stars

Arguably the most remarkable example of natural dynamo action is the solar activity cycle (Hathaway 2015; Brun & Browning 2017). Solar activity has been observed for many centuries, with both ancient Chinese and Athenian observations of sunspots being recorded. Sunspots have been systematically observed since the early seventeenth century, when Galileo turned the newly invented telescope towards our nearest star. It is now clear that solar magnetic activity exhibits an astonishingly systematic spatio-temporal behaviour. Sunspots appear in flux belts confined between the equator and latitudes of about $\pm 30^\circ$. They show cyclic activity with a period of approximately eleven years; as the cycle progresses the latitude at which activity is found moves in a wave from mid-latitude to the equator, before dying out – as shown in the solar butterfly diagram of figure 2. Note that the activity in this diagram is basically symmetric about the equator; this sunspot activity has been linked to magnetic fields through the Zeeman splitting of spectral lines. The Sun has an azimuthally averaged field of a few hundred $\mathrm {\mu }$T and the large-scale radial magnetic field changes sign across the equator and changes parity every eleven years, giving a 22-year magnetic cycle. The sunspot field reverses at sunspot minimum, whilst the coronal field reverses out of phase at sunspot maximum. The sun also possesses small-scale magnetic fields, which are only weakly coupled with the solar cycle.

Figure 2. The butterfly diagram of the solar cycle. This shows the positions of the spots for each rotation of the Sun since 1874. Magnetic activity first appears at mid-latitudes, widens and then moves towards the equator as each cycle progresses (image courtesy of D. Hathaway).

The long-term dynamics of solar activity can be established from both direct and indirect observations. Telescopic observations demonstrate the modulation of the basic eleven year cycle and also reveal the presence of a period in the seventeenth century, known as the Maunder Minimum, when sunspots were almost completely absent. Interestingly, as the Sun emerged from this minimum, sunspots were to be found solely in the southern hemisphere (Ribes & Nesme-Ribes 1993). Indirect observations of terrestrial isotopes (see Usoskin 2017), whose production rate is anti-correlated with magnetic activity, show that the solar cycle persisted through the Maunder Minimum (Beer, Tobias & Weiss 1998) and that minima such as the Maunder Minimum are recurrent events, that appear in clusters – so called ‘super-modulation’ (Weiss & Tobias 2016; Beer, Tobias & Weiss 2018).

Our understanding of the origin of the solar magnetic field owes much to observations of magnetic field generation in other stars. Stellar observations can obviously be used to calibrate dynamo theories, but will also give a clue as to past and future behaviour. The results from such observations are varied and bewildering to the non-expert (and sometimes to the expert). For our purposes it is sufficient to state the main results. For solar-like main-sequence stars, broadly speaking magnetic activity increases with rotation rate ($\varOmega$). There is a systematic increase in magnetic field strength with decreasing Rossby number of the star (here defined as $Ro= u_{rms}/2 \varOmega \ell$ where $u_{rms}$ is a typical velocity, $\ell$ is a characteristic length scale given by mixing length theory), until activity levels off past a threshold in Rossby number (typically $Ro \approx 0.1\text {--}0.3$). Interestingly, this threshold appears to be independent of the mass of the star. It is also the case that many of these stars exhibit activity cycles similar to the solar cycle. The period of the cycle is also a function of Rossby number, with the cycle period decreasing with Rossby number (i.e. faster rotators have shorter period). Finally for these stars the morphology of the generated field also appears to depend on rotation rate, with faster rotators being more dominated by their zonal fields (See et al. 2016). Other stars show a wide variety of behaviour, depending on their mass, age and spectral type (for example whether they exhibit convection in a core or in an envelope plays a major role). Young stars (e.g. fully convective T-Tauri stars) have strong fields (average field strengths of around $0.2$ T). Suffice to say that stellar magnetic fields are ubiquitous and the properties of the generated field depend on the nature of the turbulence and the rotation rate of the star. The interested reader is directed to the excellent review by Brun & Browning (2017) for more details.

1.2.4. Galaxies and galaxy clusters

The properties of magnetic field in galaxies (including our own galaxy) have been extensively measured, with a variety of observational techniques – see Brandenburg & Subramanian (2005), Beck (2015) and the reviews referenced therein. Our own galaxy has an estimated local field strength of about $5 \times 10^{-10}$ T, which is a typical amplitude for fields in galaxies. The magnetic field in galaxies usually has an ordered and a tangled component; in spiral galaxies the large-scale order manifests itself over several kiloparsecs. Usually the random component is strongest in the spiral arms, whilst the regular field extends into the interarm regions.

Galaxy clusters are the largest bound systems in the universe and are found to have magnetic fields; for these the total magnetic fields are estimated to be of the order of $10^{-11}$–$10^{-9}$ T. For most astrophysical objects it is now widely accepted that dynamo action is the only mechanism that can explain their magnetic properties and the persistence of the field for many ohmic decay times. For the galactic dynamo there appears still to be some debate (as the magnetic decay times for such vast objects are enormous), although there are a number of very significant arguments pointing to the dynamo origin of such fields (see e.g. Brandenburg & Subramanian 2005; Kulsrud & Zweibel 2008).

1.2.5. Accretion discs

Accretion discs are gaseous discs of material spinning around a central object (for example young stars, white dwarfs, neutron stars and black holes). Although direct observations of dynamo-generated magnetic fields in such objects is difficult (although see Donati et al. 2005), proxies such as emission or the magnetism of meteorites formed in the disc around the young Sun give indirect evidence. Direct observations may place upper limits on the strength of any generated magnetic fields. However, it is believed that magnetic fields play an important role in the collimation of jets and the sustainment of accretion (as we shall see later). The dynamo in such discs is very interesting theoretically – being essentially nonlinear (Tobias, Cattaneo & Brummell 2011a) – with many of the characteristics of the transition to dynamo action being similar to those for the transition to turbulence in wall-bounded flows (see e.g. Waleffe 1997; Barkley 2016). We shall discuss in details these essentially nonlinear dynamos in § 7.

2. Fundamentals

‘So many dynamos’

A well-known Palindrome

In this section we introduce the fundamentals of dynamo theory, identifying important analogies for fluid dynamicists. More details can be found for example in the excellent expositions and reviews referred to below.

2.1. Equations and boundary conditions

2.1.1. The induction and momentum equation

Dynamo theory typically (minimally) involves the construction of solutions to the pair of coupled partial differential equations for the evolution of the velocity field, ${\boldsymbol u}$ and the magnetic field ${\boldsymbol B}$ in an electrically conducting fluid under the magnetohydrodynamic (MHD) formalism. A complete derivation of the equations and a discussion of their applicability in any given physical situation is beyond the scope of this perspective; the interested reader is directed to Moffatt & Dormy (2019) and Jones (2008) for details of the derivation. Briefly, under the MHD approximation, which combines the non-relativistic Maxwell equations of electrodynamics with Ohm's law for a moving conductor, these take the form

(2.1)\begin{gather} \frac{\partial {\boldsymbol u}}{\partial t} + {\boldsymbol u} \boldsymbol{\cdot} \boldsymbol{\nabla} {\boldsymbol u} = - \frac{1}{\rho}\boldsymbol{\nabla} p + \frac{1}{\rho}\, {\boldsymbol j} \times {\boldsymbol B}+\nu \nabla^2 {\boldsymbol u} + \frac{\boldsymbol F}{\rho}, \end{gather}

(2.2)\begin{gather}\frac{\partial {\boldsymbol B}}{\partial t} = \boldsymbol{\nabla} \times ( {\boldsymbol u} \times {\boldsymbol B} ) -\boldsymbol{\nabla} \times (\eta \boldsymbol{\nabla} \times {\boldsymbol B} ), \end{gather}

(2.3)\begin{gather}\boldsymbol{\nabla} \boldsymbol{\cdot} {\boldsymbol B} = 0. \end{gather}

Here, the standard fluid parameters are $\rho$, which is the density of the fluid, $p$ which is the pressure and $\nu$ is the kinematic viscosity of the fluid, whilst ${\boldsymbol F}$ represents all the body forces (such as buoyancy or mechanical driving) acting to drive the fluid motion. The additional force term in the Navier–Stokes equations is termed the Lorentz force ${{\boldsymbol F}}_{Lor} = {\boldsymbol j} \times {\boldsymbol B}$ and arises owing to the interaction of the magnetic field with the current ${\boldsymbol j}$ flowing through the conductor. The current is itself given by the pre-Maxwell version of Ampère's law, i.e.

(2.4)\begin{equation} {\boldsymbol j} = \frac{1}{\mu} \boldsymbol{\nabla} \times {\boldsymbol B},\end{equation}

and hence

(2.5)\begin{equation} {{\boldsymbol F}}_{Lor} = \frac{1}{\mu \rho} (\boldsymbol{\nabla} \times {\boldsymbol B}) \times {\boldsymbol B}.\end{equation}

Here, and henceforth in this paper, $\mu$ is the permeability. The evolution equation for the magnetic field given by (2.2) is termed the induction equation. Here, $\eta = (\mu \sigma )^{-1}$ is the magnetic diffusivity which is a property of the conducting fluid, with electrical conductivity $\sigma$. Magnetic diffusivity is large for poor conductors. In circumstances in which the magnetic diffusivity is constant, (2.2) simplifies (utilising (2.3)) to

(2.6)\begin{equation} \frac{\partial {\boldsymbol B}}{\partial t} = \boldsymbol{\nabla} \times ( {\boldsymbol u} \times {\boldsymbol B} ) + \eta \nabla^2 {\boldsymbol B}.\end{equation}

2.1.2. Magnetic boundary conditions

The induction and momentum equations are usually (although not always) solved in a finite domain subject to the imposition of boundary conditions on the fluid velocity and magnetic field. The boundary conditions for the fluid velocity are straightforward and standard; usually one considers impenetrable boundaries with either stress-free or no-slip conditions on the tangential component of the velocity.

The magnetic boundary conditions are more problematic and usually involve some degree of simplification. Typically this involves considering the fluid in a domain $V$ with an exterior bounding surface $S$, external to which is either an insulator or a perfect conductor. If the region outside of the evolution of the dynamo fluid is perfectly conducting, surface charges, $\rho _S$, and surface currents, ${\boldsymbol j}_S$, are allowed on the boundary. If we define $[.]$ as the jump across the surface $S$ and ${\boldsymbol n}$ as an outward pointing normal to that surface, then integrating the pre-Maxwell equations across the surface gives

(2.7a–d)\begin{equation} [ {\boldsymbol n} \boldsymbol{\cdot} {\boldsymbol E} ] = \frac{\rho_S}{\epsilon}, \quad [ {\boldsymbol n} \boldsymbol{\cdot} {\boldsymbol B} ] = 0, \quad [ {\boldsymbol n} \times {\boldsymbol B} ] = \mu {\boldsymbol j}_S , \quad [ {\boldsymbol n} \times {\boldsymbol E} ] = 0, \end{equation}

where $\epsilon$ is the permittivity. If there are no surface currents or charges (i.e. for non-perfectly conducting boundaries) then ${\boldsymbol B}$ is continuous, provided $\mu$ is constant. If the region external to the solution domain does not allow currents (i.e. if this region is an insulator) then this continuity of ${\boldsymbol B}$ defines the problem. If the outside region allows currents then the normal derivative of ${\boldsymbol n} \boldsymbol {\cdot } {\boldsymbol B}$ is also continuous, although the normal derivatives of the tangential components of ${\boldsymbol B}$ are not necessarily continuous (see the detailed discussion in Jones 2008).

However, if the outside of the domain is a static perfect conductor it is normal to assume that there is no trapped magnetic field there and so (for no normal flow conditions)

(2.8a,b)\begin{equation} {\boldsymbol n} \boldsymbol{\cdot} {\boldsymbol B} =0, \quad {\boldsymbol n} \times {\boldsymbol j} =0. \end{equation}

This gives

(2.9)\begin{equation} B_z = \frac{\partial B_x}{\partial z} = \frac{\partial B_y}{\partial z} = 0, \end{equation}

at a Cartesian boundary $z=\textrm {const.}$, and

(2.10)\begin{equation} B_r = \frac{\partial (rB_{\theta})}{\partial r} = \frac{\partial(rB_{\phi})}{\partial r} = 0, \end{equation}

at a spherical boundary, $r=\textrm {const.}$

2.2. What is a dynamo?

Simply put, a self-exciting hydromagnetic dynamo is a self-consistent solution of the coupled Navier–Stokes and induction equation for which the magnetic energy,

(2.11)\begin{equation} M(t) = \int_V \frac{B^2}{2 \mu}\, \textrm{d} V, \end{equation}

remains finite as $t \rightarrow \infty$. Here $V$ is the volume over which the dynamo equations are solved, which could in principle be finite and bounded by a surface $S$ or (in an idealisation) taken to be all space. If the volume is finite it is traditional to take the region external to the domain as being an insulator (${\boldsymbol j}_{ext.}=0$) so that the magnetic field is maintained entirely by the current distribution within the domain and so ${\boldsymbol B} = { {O}}(|{\boldsymbol x}|^{-3})$ as $|{\boldsymbol x}| \rightarrow \infty$ (see e.g. Moffatt 1978).

2.3. Energetics and conservation laws

As in hydrodynamics, much can be understood by deriving global conservation laws, valid in the absence of driving and dissipation. For inviscid hydrodynamic flows in the absence of body forces (such as gravity) and boundary forces the total kinetic energy and kinetic helicity defined as

(2.12a,b)\begin{equation} E_k = \frac{1}{2}\int_V \rho {\boldsymbol u} \boldsymbol{\cdot} {\boldsymbol u} \, \textrm{d} V , \quad H_k = \int_V \rho {\boldsymbol u} \boldsymbol{\cdot} {\boldsymbol \omega} \, \textrm{d} V, \end{equation}

are conserved (Moreau 1961; Moffatt 1969).

In the presence of a magnetic field these quantities are no longer conserved. Transfer of energy may take place between kinetic and magnetic energies via the action of induction and the Lorentz force respectively. The magnetic and kinetic energy equations take the form, after a little vector calculus and ignoring surface terms,

(2.13)\begin{equation} \frac{\textrm{d}}{\textrm{d} t}\int_V \frac{{\boldsymbol B}^2}{2 \mu} \, \textrm{d} V = -\int_V {\boldsymbol u}\boldsymbol{\cdot} (\kern1.8pt {\boldsymbol j} \times {\boldsymbol B}) \, \textrm{d} V - \eta \mu \int_V {{\boldsymbol j}^2} \, \textrm{d} V,\end{equation}

and for an incompressible flow

(2.14)\begin{equation} \frac{\textrm{d}}{\textrm{d} t}\int_V \frac{\rho {\boldsymbol u}^2}{2} \, \textrm{d} V = +\int_V {\boldsymbol u}\boldsymbol{\cdot} (\kern1.5pt {\boldsymbol j} \times {\boldsymbol B}) \, \textrm{d} V - \int_V 2 \nu \rho {{\boldsymbol S}^2} \, \textrm{d} V , \end{equation}

where the dissipative terms have now been included, and $S_{ij} = \frac {1}{2}({\partial u_i}/{\partial x_j}+{\partial u_j}/{\partial x_i})$ is the rate of strain tensor for an incompressible fluid.

Adding (2.13) and (2.14) for an ideal (inviscid and perfectly conducting) fluid immediately reveals that the total (kinetic plus magnetic) energy is conserved, with magnetic energy only being created at the expense of kinetic energy. The ${\boldsymbol u}\boldsymbol {\cdot } (\kern1.5pt {\boldsymbol j} \times {\boldsymbol B})$ term is therefore responsible for this transfer of energy and arises from inductive effects in the induction equation and the Lorentz force in the momentum equation.

In addition to the total energy there are two other quadratic invariants of the ideal system, namely the magnetic helicity and the cross-helicity. The cross-helicity, given by

(2.15)\begin{equation} H^c = \int_V {\boldsymbol u} \boldsymbol{\cdot} {\boldsymbol B} \, \textrm{d} V , \end{equation}

is clearly not sign definite; neither is cross-helicity dissipation. Cross-helicity may be either amplified or damped locally and so less attention has focussed on the implications of its conservation in the non-dissipative case (although see Biskamp 2003; Yokoi 2013).

However, the implications of the conservation of magnetic helicity and its role in the saturation of nonlinear dynamos has received much attention and for that reason we devote some time to it here.

2.3.1. Conservation of magnetic helicity

Because $\boldsymbol {\nabla } \boldsymbol {\cdot } {\boldsymbol B} = 0$, it is often useful to write ${\boldsymbol B} = \boldsymbol {\nabla } \times {\boldsymbol A}$. Here ${\boldsymbol A}$ is termed the vector potential for the magnetic field. Clearly, ${\boldsymbol A}$ is only defined up to a choice of gauge, so that transforming ${\boldsymbol A} \rightarrow {\boldsymbol A} + \boldsymbol {\nabla } \psi$ for any scalar $\psi$ leaves the magnetic field unchanged. It is also then convenient to ‘uncurl’ (2.2) to give

(2.16)\begin{equation} \frac{\partial {{\boldsymbol A}}}{\partial t} =({{\boldsymbol u}}\times {{\boldsymbol B}})-\eta \,\boldsymbol{\nabla} \times {{\boldsymbol B}} + \boldsymbol{\nabla} \phi, \end{equation}

where $\phi$ is related to the choice of gauge.

Various choices of gauge are possible; the Coulomb gauge has

(2.16a,b)\begin{equation} \boldsymbol{\nabla} \boldsymbol{\cdot} {\boldsymbol A} = 0; \quad \nabla^2 \phi = -\boldsymbol{\nabla} \boldsymbol{\cdot} ({{\boldsymbol u}}\times {{\boldsymbol B}}). \end{equation}

The winding gauge (suitable for calculations in Cartesian geometries) has

(2.17a,b)\begin{equation} \nabla_H \boldsymbol{\cdot} {\boldsymbol A} = 0; \quad \nabla_H^2 \phi^\prime = -\boldsymbol{\nabla}_H \boldsymbol{\cdot} ({{\boldsymbol u}}\times {{\boldsymbol B}}), \end{equation}

where $\nabla _H = (\partial _x,\partial _y,0)$ and $\phi ^\prime = \phi - \eta ({\partial A_z}/{\partial z})$ (Prior & Yeates 2014). A numerically convenient gauge involves setting

(2.17)\begin{equation} \phi = \frac{\partial \psi}{\partial t}, \end{equation}

as described in Brandenburg & Subramanian (2005)

Now, magnetic helicity is defined as

(2.18)\begin{equation} H = \int_V {\boldsymbol A} \boldsymbol{\cdot} {\boldsymbol B} \, \textrm{d} V ,\end{equation}

and its evolution can be easily shown to be given by

(2.19)\begin{equation} \frac{\textrm{d} H}{\textrm{d} t}=-2 \eta \mu \int_V {\boldsymbol j} \boldsymbol{\cdot} {\boldsymbol B} \, \textrm{d} V + F_s,\end{equation}

where $F_s$ represents the surface flux of magnetic helicity. Hence, for a perfectly conducting fluids, with no loss or gain of magnetic helicity through the boundaries (i.e. $F_s = 0$), magnetic helicity is conserved. This has implications for dynamo action, which requires the generation or destruction of magnetic helicity, as we shall see. Magnetic helicity is a measure of the topological complexity and linkage of field lines, and in the absence of diffusive processes (by which the field can reconnect) this complexity is maintained. Furthermore magnetic helicity appears to be a more robust invariant than say total energy in the presence of small diffusion. The dissipative term for total energy, $-\eta \mu \int _V \,{{\boldsymbol j}^2}\, \textrm {d} V$, may remain finite as $\eta \rightarrow 0$, because ${\boldsymbol j}^2$ gets large in this limit. However, the dissipation of magnetic helicity, $-2 \eta \mu \int _V {\boldsymbol j} \boldsymbol {\cdot } {\boldsymbol B} \, \textrm {d} V$, appears to tend to zero as $\eta \rightarrow 0$, so magnetic helicity is well conserved. Of course the situation changes if magnetic helicity is allowed to enter or leave the domain of interest; it may do so via either advective or diffusive processes, as we shall see in § 6.5.2.

In hydrodynamic turbulence the presence of quadratic invariants has consequences for the nature of the cascades. Briefly, the same reasoning applies to MHD. As argued above, for small dissipation energy decays faster than magnetic helicity and cross-helicity (Biskamp 2003). Therefore, in a turbulent state, energy cascades towards small scales (analogous to the energy cascade in three-dimensional (3-D) hydrodynamics). However, the magnetic helicity cascades toward large scales (analogous to the energy cascade in 2-D hydrodynamics). The inverse cascade of magnetic helicity may lead to the formation of large-scale magnetic fields – this fact may prove to be important for the generation of systematic fields by nonlinear dynamos.

2.4. The induction equation and kinematic dynamos: the basics

For fluid dynamicists an obvious useful analogy can be made between the induction (2.2) and the incompressible vorticity equation given by

(2.20)\begin{equation} \frac{\partial {\boldsymbol \omega}}{\partial t} = \boldsymbol{\nabla} \times ( {\boldsymbol u} \times {\boldsymbol \omega} ) +\nu \nabla^2 {\boldsymbol \omega}. \end{equation}

For experts on vorticity dynamics, this analogy can lead to significant insight into the dynamics of the magnetic field. Magnetic flux tubes may in certain circumstances have a similar dynamics to that of vortex tubes; this analogy becomes a formal correspondence when the magnetic field and vorticity are weak. However, it is important not to push the analogy too far. Equation (2.20) is a nonlinear evolution equation for the vorticity (since the advecting velocity is itself related to the vorticity), whereas (2.6) is formally linear in the magnetic field, with the system only becoming nonlinear when coupled to the momentum equation via the Lorentz force. As a rule of thumb, if the magnetic field is weak compared with the velocity it behaves analogously to the vorticity; if it is of a similar strength it behaves more like the velocity.

Important limits of the induction equation are the so-called diffusive and perfectly conducting limits. If ${\boldsymbol u} = 0$, (2.6) reduces to the vector diffusion equation,

(2.21)\begin{equation} \frac{\partial {\boldsymbol B}}{\partial t} = \eta \nabla^2 {\boldsymbol B}, \end{equation}

so that for no fluid motion the field must diffuse away (assuming there are no fields at infinity). Hence motion is needed to maintain magnetic field. The time scale for diffusion of field with a typical length scale $\ell _B$ is given by $\tau _D=\ell _B^2 /\eta$. To give some idea of some typical diffusive time scales, we note that $\tau _d \sim 2$ seconds for $\ell _B = 1$ m and $\eta =0.04$ m$^2$ s$^{-1}$, which may be appropriate for a liquid sodium experiment. Whereas for a magnetic field on the scale of the Earth's core and a relevant diffusivity gives $\tau _d \sim 10^4$ years; the Sun has a diffusion time for the large-scale magnetic field of $10^9$ years, comparable with its age. For galaxies the diffusion time of magnetic field is significantly longer than the age of the universe!

In the absence of diffusion (the so-called perfectly conducting limit where $\sigma \rightarrow \infty$ and $\eta =0$), (2.6) becomes

(2.22)\begin{equation} \frac{\partial {\boldsymbol B}}{\partial t} = \boldsymbol{\nabla} \times ( {\boldsymbol u} \times {\boldsymbol B} ). \end{equation}

Sometimes this is called the frozen flux limit; the magnetic flux $\int _{\mathcal {S}} {\boldsymbol B} \boldsymbol {\cdot } \textrm {d} {\boldsymbol S}$ through the surface $\mathcal {S}$ bounded by any closed curve $\mathcal {C}$ moving with the fluid, remains fixed (Alfvén's theorem). Hence we can think of magnetic field as being frozen into the fluid, in a similar manner to vortex lines in an inviscid fluid, Alfvén's theorem is the magnetic counterpart of Kelvin's circulation theorem.

2.4.1. Importance of the magnetic Reynolds number $Rm$

Clearly from the above discussion, magnetic field will decay away unless the advective term in the induction equation is large enough to overcome diffusive effects. The relative importance of the two terms $\boldsymbol {\nabla } \times ({\boldsymbol u} \times {\boldsymbol B})$ and $\eta \nabla ^2 {\boldsymbol B}$ can be established by non-dimensionalising. We choose a typical length scale $L$ which is the size of the object or region under consideration and a typical fluid velocity $U$. On introducing scaled variables $t = ({L}/{U}) \tilde {t}$, ${\boldsymbol x} = L {\tilde {\boldsymbol x}}$, ${\boldsymbol u} = U \tilde {\boldsymbol u}$, and dropping tildes the induction equation becomes

(2.23)\begin{equation} \frac{\partial {\boldsymbol B}}{\partial t} = {\boldsymbol \nabla} \times ({\boldsymbol u} \times {\boldsymbol B}) + Rm^{-1} \nabla^2 {\boldsymbol B},\end{equation}

where $Rm = {U L}/{\eta }$ is the dimensionless magnetic Reynolds number.

In general, large $Rm$ means induction dominates over diffusion, whilst small $Rm$ means diffusion wins out over induction, and as we shall see minimum values of $Rm$ for dynamo action (so-called dynamo bounds) can sometimes be found.

It is worth noting at this point, however, that $Rm$ should only be used as a guide to determine the relative importance of advection and diffusion. In defining $Rm$ it has explicitly been assumed that $L$ is a typical length scale for both the magnetic field and the velocity (i.e. that $\ell _B = \ell _U = L$). This makes sense if $L$ is the size of the astrophysical object (i.e. the largest length scale available). Of course this may not be the case, with the potential for $\ell _B$ to be very different from $\ell _U$. For example if $\ell _B \gg \ell _U$ then the relative importance of advection and diffusion is given by $Rm_T = U \ell _B^2 / \eta \ell _U \sim \omega \ell _B^2/\eta$, where $\omega$ is a typical vorticity amplitude. This may be large even if $\ell _U$ and $U$ are small. Such a basic misunderstanding of the limitations of the information encoded in the magnetic Reynolds number may lead to the incorrect dismissal of certain small-scale flows as the possible origin of large-scale fields. Large-scale fields, as they are weakly diffusive, require very little induction for their maintenance.

2.4.2. A useful technique: axisymmetric field decomposition

It is clear from the form of the induction equation that a non-axisymmetric flow immediately creates a non-axisymmetric field, however, the converse is not true. If both the flow and field are axisymmetric then one can decompose the flow and field by setting (in cylindrical polars $(s,\phi ,z)$)

(2.24)\begin{gather} {\boldsymbol u}=s \varOmega \hat{{\boldsymbol \phi}} + {\boldsymbol u}_P = s \varOmega \hat{{\boldsymbol \phi}} + \boldsymbol{\nabla} \times \left(\frac{\psi}{s}\right) \hat{{\boldsymbol \phi}}, \end{gather}

(2.25)\begin{gather}{\boldsymbol B} = B \hat{{\boldsymbol \phi}} + {\boldsymbol B}_P =B \hat{{\boldsymbol \phi}} + \boldsymbol{\nabla} \times (A \hat{{\boldsymbol \phi}}). \end{gather}

Here, $s = r \sin \theta$ (where $r$ and $\theta$ relate to spherical polars $(r,\theta ,\phi )$). This defines the differential rotation $\varOmega$ the streamfunction $\psi$, the zonal field $B$ (sometimes termed toroidal field in an axisymmetric setting) and the scalar potential $A$. The induction equation then simplifies to

(2.26)\begin{gather} \frac{\partial A}{\partial t} + \frac{1}{s} ({\boldsymbol u}_P \boldsymbol{\cdot} \boldsymbol{\nabla} ) (sA) = \eta \left(\nabla^2 - \frac{1}{s^2}\right) A, \end{gather}

(2.27)\begin{gather}\frac{\partial B}{\partial t} + s ({\boldsymbol u}_P \boldsymbol{\cdot} \boldsymbol{\nabla} ) \left(\frac{B}{s}\right) = \eta \left(\nabla^2 - \frac{1}{s^2}\right) B + s{\boldsymbol B}_P \boldsymbol{\cdot} \boldsymbol{\nabla} \varOmega. \end{gather}

The form of these equations will be exploited in the next section to prove so-called anti-dynamo theorems. Note that, although both the $A$ and $B$ equations have advective and diffusive terms, the field stretching term, $s{\boldsymbol B}_P \boldsymbol {\cdot } \boldsymbol {\nabla } \varOmega$, only appears in (2.27) if gradients in angular velocity are present. (2.26) shows that there is no corresponding source term in the equation for $A$.

2.4.3. Anti-dynamo theorems

It is fair to say that the psychology of the dynamo practitioner has been strongly shaped by the early results of dynamo theory – results that showed the complexity and difficulty of achieving dynamo-generated fields. Key are the so-called anti-dynamo theorems that show the impossibility of dynamo action for large classes of magnetic fields and velocity fields with certain symmetries. We shall not reproduce all of the demonstrations and proofs here, since they are readily available from many sources (for example Dormy & Soward 2007; Jones 2008; Moffatt & Dormy 2019), but the importance of these for the subsequent direction of the development of the field cannot be overstated.

Cowling's theorem: an axisymmetric field vanishing at infinity cannot be maintained by dynamo action (Cowling 1933).

Proof. As noted above, a non-axisymmetric velocity field immediately generates non-axisymmetric magnetic field and so it is necessary to consider only axisymmetric flows and fields. The evolution of the field is therefore given by (2.26)–(2.27). Multiplying (2.26) by $s^2 A$, using the divergence theorem, and integrating over all space gives

(2.28)\begin{equation} \frac{\textrm{d}}{\textrm{d} t} \int \frac{1}{2} s^2 A^2 \, \textrm{d} V = - \eta \int | \boldsymbol{\nabla} (sA) |^2 \, \textrm{d} V , \end{equation}

assuming that surface terms vanish at infinity. This equation clearly shows that $sA$ decays to zero as $t \rightarrow \infty$; note $sA$ can not decay to a constant since this would imply $A\rightarrow \infty$ as $s\rightarrow 0$. Once $A$ has decayed there is no source term in the toroidal field equation (i.e. 2.27). Similar arguments then ensure the subsequent decay of $B/s$ and therefore the impossibility of the creation of an axisymmetric magnetic field by dynamo action. It is important to note that it is the absence of source terms in (2.26) that causes the problems for dynamo action. Relaxation of the constraint of axisymmetry allows for the re-inclusion of such a source term.

The effect of Cowling's theorem ruling out such simple symmetric solutions made the search for any dynamo solutions (which were necessarily three dimensional!) seem a formidable task. Indeed, for a long time it was not clear that any such dynamo solutions existed. The situation is encapsulated in the following story taken from Krause (1993) (which attributes the source as Paul Roberts) ‘Walter Elsasser and Einstein were friends in Germany before they both emigrated to the US in the 1930s. Several years after Elsasser had settled there (in the late 1930s in fact), he became interested in the origin of the geomagnetic field. Einstein paid him a visit, and (as people do) asked ‘What are you working on these days?’. Elsasser told him, and Einstein invited him to explain dynamo theory to him. Elsasser set-up the problem and then told Einstein about Cowling's theorem. Einstein's response was, ‘If such simple solutions are impossible, self-excited fluid dynamos cannot exist’. For once, the great man's craving for simplicity seems to have misled him’.

Other antidynamo theorems: there are many extensions of Cowling's antidynamo theorem to other geometries and to slightly different set-ups. For example, in Cartesian coordinates $(x,y,z)$, no magnetic field that vanishes at infinity and is independent of $z$ can be generated by dynamo action; the proof proceeds along similar lines to that given above (see e.g. Jones 2008). Note again that this is a restriction on the form of a dynamo-generated magnetic field.

Anti-dynamo theorems placing constraints on the form of the fluid velocities that can lead to dynamo action have been proven by Bullard & Gellman (1954), Backus (1958) and in Cartesian coordinates by Zel'dovich (1957). These results essentially show that a velocity field must have all three components in order to be capable of acting as a dynamo (see the long discussion and derivation in Moffatt 1978).

In conclusion, both the field and the flow must be sufficiently complicated for dynamo action to occur. A minimal requirement is that the field must be three-dimensional and the flow must not be purely poloidal (i.e. cannot be planar).

2.4.4. Bounds on dynamo action

Even for flows and fields that are not ruled out as dynamos on symmetry grounds, there are bounds that constrain dynamo action in finite domains; a crude measure of the expected efficiency of a dynamo is given by the magnetic Reynolds number $Rm$, which gives the ratio of advection to diffusion at a particular scale in the flow (usually taken to be the integral or system scale). The importance of $Rm$ can be formalised in bounds derived by Backus (1958) and Childress (1969), both of which use the evolution equation for the magnetic energy $E_M$, given by (see (2.13))

(2.29)\begin{equation} \mu \frac{\partial E_M}{\partial t} = - \eta \int | \boldsymbol{\nabla} \times {\boldsymbol B}|^2 \, \textrm{d} V + \int (\boldsymbol{\nabla} \times {\boldsymbol B}) \boldsymbol{\cdot} ({\boldsymbol u} \times {\boldsymbol B}) \, \textrm{d} V . \end{equation}

For field generation in a sphere of radius $a$, matching to a decaying potential dynamo action requires

(2.30)\begin{equation} Rm = \frac{a^2 e_{max}}{\eta} \ge {\rm \pi}^2, \end{equation}

where $e_{max}$ is the maximum of the rate of strain tensor (Backus 1958). Note here that $Rm$ is defined in terms of the maximum strain (and not a typical velocity amplitude), which seems natural. A slightly different bound from (2.29) (Childress 1969) requires

(2.31)\begin{equation} {Rm = \frac{a u_{max}}{\eta} \ge {\rm \pi}}. \end{equation}

Other bounds on dynamo action are also possible. However, it can be shown that a steady or periodic dynamo can exist in a bounded conductor with an arbitrarily small value of the kinetic energy. Hence there is no lower bound on dynamo action when $Rm$ is defined using the root mean square (rather than maximum) velocity, without placing limitations on the rate of strain (Proctor 2015).

Motivated by the desire to construct experimental dynamos in the laboratory, there is currently an effort to optimise the efficiency of dynamo flows given certain constraints. This involves taking the machinery developed for understanding the transition to turbulence, such as adjoint optimisation methods and using them to maximise the efficiency of dynamo velocities (Willis 2012).

2.5. Kinematic dynamos: some simple flows that work

The kinematic problem considers the induction equation in isolation, for a prescribed velocity field ${\boldsymbol u}$. It is then natural to consider flows that are either steady, periodic in time or statistically steady. Because the induction equation is linear in the magnetic field, solutions take the form of magnetic fields that grow or decay exponentially on average and the task for the dynamo theorist is to determine which flows lead to growing solutions (i.e. which can act as dynamos) and then perhaps examine the form of the solutions as a function of $Rm$.

2.5.1. Kinematic dynamos

Having convinced ourselves that overly symmetric fields and flows are not good for dynamo action and that sufficient stretching is required, it is time to discuss some flows that do actually work as dynamos. These are not presented in chronological order, as we shall start with the simpler case of flows in an infinite domain, before moving onto flows in spheres and spherical shells.

The Ponomarenko flow (Ponomarenko 1973):

Here, we consider the simplest possible flow that leads to dynamo action. It takes the form of a localised discontinuous ‘screw’ or vortex flow that in cylindrical coordinates $(s, \phi , z)$ is given by

(2.32)\begin{equation} {\boldsymbol u} = \begin{cases} s \varOmega {\hat{{\boldsymbol \phi}}} + U \hat{\boldsymbol z}, \quad s < a,\\ 0, \quad s > a, \end{cases} \end{equation}

where $\varOmega$ and $U$ are constants. Note that this flow is not planar, owing to the presence of the throughflow $U \hat {\boldsymbol z}$, and has kinetic helicity, $H = {\boldsymbol u} \boldsymbol {\cdot } {\boldsymbol \omega } = 2 U \varOmega .$ We shall see the importance of kinetic helicity for large-scale dynamos in § 5. Strong shear naturally occurs at (the physically unrealistic) discontinuity at $s=a$. In principle there are two independent parameters defining the flow $U$ and $\varOmega$. These can be re-expressed in terms of an overall amplitude of the flow given by $U_{amp} = (U^2 + a^2 \varOmega ^2)^{1/2}$ and the pitch angle $\chi = U / a \varOmega$.

The trick for elegant solution for the kinematic dynamo modes in such a configuration is to note first that the flow is steady and, together with the linearity of the induction equation, this implies that magnetic field will either grow or decay exponentially in time, with potentially a complex growth-rate $\lambda = \sigma + {\rm i} \omega$. Secondly, the flow is independent of $z$ and $\phi$ and therefore monochromatic magnetic fields can be sought in those directions. It is therefore advantageous to seek solutions of the form ${\boldsymbol B} = {\boldsymbol b}(s) \exp ( \lambda t + {\rm i} m \phi + {\rm i}kz ),$ where $m$ and $k$ are the azimuthal and vertical wavenumbers, which must be non-zero to avoid the anti-dynamo theorems discussed earlier.

This model is extremely illuminating as it can be solved pseudo-analytically (see Jones 2008) for details; marginally stable solutions can be found by setting $\mbox {Re}({\lambda }) = 0$ (for a given $Rm = a U_{amp}/\eta$, $\chi$, $ka$ and $m$). The Ponomarenko flow is indeed a dynamo! It can be thought of as a prototype dynamo that is a model of vortical plume. Dynamo action sets in at $Rm_{crit}=17.72$, for $ka_{crit} = -0.3875$, $m=1$, $a^2 \omega / \eta = -0.41$ and $\chi = 1.31$, so the optimal pitch is $ {O}(1)$. The magnetic field is strongest near $s=a$, where it is generated by (the unphysical) shear. Although it is important to examine how dynamo action onsets, it is also, as we shall see, of great interest to determine the behaviour at large $Rm$. Asymptotic solutions of the Ponomarenko dynamo (in the form of Bessel functions) show that (Gilbert 1988) the fastest growing modes are given by

(2.35a,b)\begin{equation} |m| = (6(1 + \chi^{-2}))^{-3/4} \left( \frac{a^2 \varOmega}{2 \eta} \right)^{1/2},\quad \sigma = 6^{-3/2} \varOmega (1 + \chi^{-2})^{-1/2}. \end{equation}

In the presence of viscosity, velocities with discontinuities are not realistic and so the Ponomarenko dynamo has been extended to the continuous case where ${\boldsymbol u} = s \varOmega (s) \hat {{\boldsymbol \phi }} + U(s) \hat {\boldsymbol z}$, where $\varOmega (s)$ and $U(s)$ are smooth functions (Gilbert 1988).

The Roberts flow (Roberts 1972a):

Perhaps the most illuminating, kinematic dynamo flow is the so-called G.O. Roberts flow. This flow, like the Ponomarenko flow, is specially crafted to circumvent both the Cartesian version of Cowling's anti-dynamo theorem and the planar velocity anti-dynamo theorem, whilst still remaining fairly tractable.

The G.O. Roberts flow is the special case of the well-studied ABC (Arnol'd, Beltrami and Childress) flow in Cartesian coordinates $(x,y,z)$ in an infinite domain given by

(2.36)\begin{equation} {\boldsymbol u} = (C \sin z + B \cos y, A \sin x + C \cos z, B \sin y + A \cos x). \end{equation}

For the Roberts flow $A=B=1$, $C=0$, so that ${\boldsymbol u}(x,y) = (\cos y, \sin x, \sin y + \cos x).$ This flow is two-dimensional (in the sense that it only depends on two coordinates), but it has all three components (thus not falling foul of the planar flow anti-dynamo theorem).

Before discussing the dynamo properties of the Roberts flow, it is useful to describe some of its basic hydrodynamic properties. The Roberts flow is integrable in the sense that it can be written in terms of a single steady streamfunction $\psi (x,y)$, so that

(2.37)\begin{equation} {\boldsymbol u} = \left( \frac{\partial \psi}{\partial y}, -\frac{\partial \psi}{\partial x}, \psi \right), \end{equation}

where $\psi = \sin y + \cos x$. The flow takes the form of an array of helical cells with throughflow, see figure 3(a). It has a typical horizontal spatial scale of $2 {\rm \pi}$ and an infinite vertical scale. Moreover the winding sense of each helix in the array is the same. This ensures that the normalised relative kinetic helicity, defined to be

(2.38)\begin{equation} \mathcal{H}_{rel} = \int \frac{| {\boldsymbol u}\boldsymbol{\cdot} {\boldsymbol \omega}| }{| {\boldsymbol u}\boldsymbol{\cdot} {\boldsymbol u} |^{1/2}| {\boldsymbol \omega}\boldsymbol{\cdot} {\boldsymbol \omega}|^{1/2}}\,\textrm{d} V , \end{equation}

is unity (a so-called maximally helical flow); the importance of helicity for the large-scale dynamo properties of the flow will be discussed at some length in § 5.3.1.

Figure 3. (a) Contours of the streamfunction $\psi$ for the G.O. Roberts flow. Positive (and zero) contours are solid and negative contours are dashed. (b) Growth rate $\sigma$ as a function of wavenumber, $k_z$ for various $Rm$ (after Roberts 1972a). (c,d) Scaled magnetic energy in the plane $z=0$ for two different $Rm=16$ and $512$. As $Rm$ is increased the field is expelled into magnetic boundary layers of width $ {O}(Rm^{-1/2})$. Note only the domain between $0$ and $2 {\rm \pi}$ is shown. The magnetic energy is scaled between $0$ and $1$. Figure courtesy of A. Clarke.

Roberts utilised the same considerations as Ponomarenko (although a year previous) to search for magnetic field solutions to the induction equation of the form

(2.39)\begin{equation} {\boldsymbol B} = {\boldsymbol b}(x,y) \exp ( \lambda t + {\rm i}k_z z), \end{equation}

where ${\boldsymbol b}(x,y)$ is periodic in $x$ and $y$ (including the possibility that ${\boldsymbol b}(x,y)$ has a mean part). The solution of the 2-D problem requires spectral methods yielding a matrix eigenvalue problem for the growth rate $\sigma = \mbox {Re}(\lambda )$ as a function of the wavenumber $k_z$ and $Rm$. For each value of $Rm$ there is an optimal value of $k_z$ that maximises the growth rate, as shown in figure 3(b).

The form of the magnetic field for this flow is very illuminating. Figure 3(c,d) shows the magnetic energy in the plane $z=0$. For this steady flow the field is generated by the flow between the stagnation points. As $Rm$ is increased the field is expelled into magnetic boundary layers of width $ {O}(Rm^{-1/2})$ enhancing diffusion (although the length of the filamentary field structures is determined by the geometry of the flow). This is therefore an example of a small-scale dynamo – the field generated is dominated by structure at a scale smaller than a typical scale in the velocity; such a dynamo relies on the stretching overcoming the dissipative effects of diffusion.

In a tour de force paper, Soward (1987) analysed the behaviour of the Roberts dynamo at high $Rm$ utilising asymptotic methods. He showed that the maximal growth rate for this dynamo $\sigma \rightarrow 0$ as $Rm \rightarrow \infty$ (although very slowly – indeed

(2.40)\begin{equation} \sigma \sim \frac{\log \log Rm}{\log Rm}, \end{equation}

for large $Rm$). In the parlance of dynamo theory this makes the Roberts flow a slow dynamo – we shall discuss slow, fast and quick dynamos later.

Finally for this section we stress again that the Roberts flow is a small-scale dynamo that generates field via stretching. Although the flow is helical, this is not its defining characteristic here. Indeed, Roberts also considered a flow with no net helicity, viz. ${\boldsymbol u}=(\sin 2y, \sin 2x, \sin (x+y))$. Although this is a less efficient dynamo than the helical Roberts flow, it is nonetheless a (slow) dynamo. As noted by H.K. Moffatt ‘Helicity is not essential for dynamo action, but it helps’.

2.5.2. Kinematic dynamos in a spherical domain

The examples above are instructive (and similar types of flows will be utilised later to illustrate other dynamo properties). A natural question to pose, however, is whether similar types of flow can lead to dynamos in a bounded domain – the most natural examples of which are spheres or spherical shells.

In spherical geometry one may decompose an incompressible velocity field in terms of two scalar fields (Bullard & Gellman 1954), i.e.

(2.41)\begin{equation} {\boldsymbol u} = \boldsymbol{\nabla} \times ( T(r,\theta,\phi,t)\,\hat{\boldsymbol r} ) + \boldsymbol{\nabla} \times \boldsymbol{\nabla} \times ( S(r,\theta,\phi,t)\, \hat{\boldsymbol r} ).\end{equation}

Here, $T$ and $S$ are the toroidal and poloidal components. Alternatively we write

(2.42)\begin{equation} {\boldsymbol u} = \sum_{l,m} {\boldsymbol t}_l^m + {\boldsymbol s}_l^m, \end{equation}

where ${\boldsymbol t}_l^m$ and ${\boldsymbol s}_l^m$ are given by

(2.43)\begin{gather} {\boldsymbol t}_l^m = \boldsymbol{\nabla} \times ( t_l^m (r,t) Y_l^m (\theta, \phi)\,\hat{\boldsymbol r}), \end{gather}

(2.44)\begin{gather}{\boldsymbol s}_l^m = \boldsymbol{\nabla} \times \boldsymbol{\nabla} \times (s_l^m (r,t) Y_l^m (\theta, \phi)\, \hat{\boldsymbol r}) \end{gather}

where $-l \le m \le l$, and $Y_l^m$ is the spherical harmonic.

The flow is defined by choosing the values of $l$ and $m$ and the corresponding radial functional form for the scalar fields $t_l^m(r,t)$ and $s_l^m(r,t)$. In a landmark paper, Bullard and Gellman, nearly twenty years before the flows discussed in the last section, chose ${\boldsymbol u} = \epsilon {\boldsymbol t}_1^0 + {\boldsymbol s}_2^2$ and set $t_1^0 = r^2 (1-r)$, and $s_2^2 = r^3 (1-r)^2$. Having made a similar expansion to (2.39) for the magnetic field, they used spectral interaction rules to determine the growth rate of the field. They reported dynamo action for this flow for high enough $Rm$, but unfortunately the dynamo growth was spurious (as shown by subsequent higher resolution calculations). We shall return to this unfortunate property of dynamo calculations in § 2.5.3. Although the reported dynamo action was incorrect, the work of Bullard and Gellman revitalised the field after the depressing wilderness years of anti-dynamo theorems, giving hope that self-excited dynamos were indeed possible. Generalisations and extensions of the Bullard and Gellman type dynamos have been studied by Kumar & Roberts (1975) and Dudley & James (1989).

If the non-axisymmetric components of such flows are small ($ {O}(\epsilon )$) then it is possible to perform an asymptotic expansion with the axisymmetric components of flow and field dominating over the non-axisymmetric parts. This is the nearly axisymmetric dynamo of Braginskii (1975), which is an example of a self-consistent mean field model (see § 5).

2.5.3. A note of caution

The evolution equation for the magnetic field as described by the induction equation takes the form of a competition between magnetic field stretching (advection) and diffusion. Often these are both exponential processes; whether magnetic field grows or decays is then determined by the small difference between the efficiency of these processes. It is therefore of vital importance that both of these processes are accurately represented by the numerical scheme utilised to solve the induction equation. Failure to do so will inevitably lead to the incorrect determination of the dynamo properties of the flow. This was first demonstrated by the spurious solutions to the Bullard-Gellman dynamo (Bullard & Gellman 1954). Owing to computational limitations, the resolution chosen for the spectral scheme was not sufficient to resolve the dissipative structures (current sheets) in the magnetic field. Hence the dissipation was underestimated and dynamo action was claimed when no such sustained magnetic field generation was possible. It is always this way round. Insufficient resolution in a dynamo calculation will lead to the system appearing to be a dynamo when in fact it is not. This is troubling to computational dynamo theorists, although not as troubling to some as it should be.

It should be clear from the above discussion that any misrepresentation of the magnetic diffusion in a dynamo calculation is to be avoided at all costs. This includes not only misrepresentations that arise owing to lack of resolution, but also those that arise owing to the numerical scheme. It is often the case in fluids calculations that sub-grid processes are modelled (say by the inclusion of a hyperdiffusion, a numerical fix or by nominally solving the diffusionless problem with a stable scheme and using numerical errors to take the form of the dissipative processes (see e.g. Miesch 2015). In hydrodynamics, where there may be many competing terms in the nonlinear evolution equation for the velocity and dissipation may be small, such schemes may not be too damaging – but for dynamo calculations extreme care must be taken in interpreting the results from such schemes. Let me stress again that if one is examining the competition between two processes in the linear induction equation, one of which is known to be incorrectly represented, there is a strong chance that the results are incorrect.

3. So what is the problem then?

The previous sections demonstrated that dynamo action is possible in simple steady flows in the kinematic regime, i.e. for a prescribed velocity field, in both infinite and finite domains. Despite the restrictions placed on dynamo fields and the velocities that generate them, growing solutions for the magnetic field are possible.

The rest of this Perspective will focus on the current issues that are troubling dynamo theorists. In this section we shall introduce each issue, indicating its importance for our understanding, before describing in subsequent sections the past and current attempts at a resolution of each issue and concluding in § 10 with possible future lines of research (many of which are based on current techniques utilised in hydrodynamics).

3.1. Turbulence – high and low $Pm$

The flows considered above are defined at a single spatial scale and are steady. Of course naturally occurring flows in geophysics and astrophysics are neither. In general, owing to the vast length scale of such flows, the Reynolds numbers of the dynamo flows are enormous and the flows are extremely turbulent – hence the title of this Perspective. Moreover similar arguments may pertain to the magnetic Reynolds number $Rm$. The ratio of these two non-dimensional numbers is given by the magnetic Prandtl number $Pm = \nu /\eta$, which is a property of the fluid/plasma. In geophysics and astrophysics $Pm$ is usually either very large or very small; in numerical experiments it is usually $ {O}(1)$. The combination of turbulence, with its large range of spatial and temporal scales to be resolved and the naturally occurring extreme values of $Pm$ presents a formidable problem to the theorist and the numericist.

When $Pm$ is large ($Rm \gg Re$) the magnetic field can be generated on scales much smaller than the viscous dissipation scale. How dynamos (at least kinematically) behave in this regime is the preserve of fast dynamo theory, which is discussed in § 4.1. The situation is reversed (and more complicated) when $Pm$ is small ($Re \gg Rm$); in this case the magnetic field dissipates in the inertial range of the turbulence; with substantial implications for the dynamo – this case is discussed in §§ 4.2.1 and 4.2.3.

3.2. Organisation

The dynamos described so far tend to be small-scale dynamos in the sense that they generate field kinematically on a scale smaller than the typical velocity scale (sometimes on the diffusive scale, which gets very small as $Rm$ increases). However, observed geophysical and astrophysical magnetic fields display organisation on the scale of the astrophysical object and are sometimes called large-scale dynamos. A natural question, therefore, is ‘What is the origin of this organisation and how does the mechanism leading to organised flows compete with that leading the production of small-scale fields?’ Can this competition be understood within a kinematic framework, or are nonlinear effects from the momentum equation required? Is it the case that the organisation arises as a consequence of turbulent interactions or despite them? These questions are discussed in § 5, where the concept of mean-field electrodynamics will be introduced and critiqued.

3.3. Saturation

From the point of view of a fluid dynamicist, the focus of dynamo theorists on the induction equation is somewhat mysterious, although, to be fair, the difficulties in finding dynamo solutions described above go some way to explaining this preoccupation. However, in the present era of massive computations, it makes no sense to make the kinematic assumption and take the velocity as prescribed. What is required is to solve the coupled Navier–Stokes and induction equations. Questions that can be addressed within this framework include: How does a turbulent small-scale dynamo saturate? Do organised fields saturate in a similar manner? Are there dynamos that operate through instabilities of a magnetic field driving a flow – if so, can these lead to subcritical dynamo action? These questions are addressed in § 6

3.4. The role of rotation – rapid or otherwise

As we shall see in § 5, the dynamo generation of organised magnetic field requires a breaking of reflexional symmetry of the system (Moffatt 1978). In geophysical and astrophysical systems this naturally occurs via the influence of rotation (often in combination with stratification). Rotation is responsible for providing correlations that lead to the generation of a net electromotive force, in much the same way as it can drive large-scale flows via the introduction of correlations and non-trivial Reynolds stresses (see e.g. Vallis 2006). For dynamos, rotation plays a key role in determining the solutions of the induction equation. Moreover, once the field is generated it becomes dynamic in the momentum equation. In the case of rapid rotation, such as is the case in Earth's core and in rapidly rotating stars, this is a particularly interesting and delicate issue. Magnetic fields may break leading-order geostrophic balances (leading to so-called magnetostrophic balance) or, even if the primary balance is still geostrophic, the fields may play a crucial role in the prognostic equation for unbalanced motions, as discussed in § 8. Furthermore, the presence of strong rotation can lead to dynamo-generated magnetic field acting as a conduit to turbulence via nonlinear processes. The magnetic field can act so as to relax the constraints engendered by strong rotation, leading to more efficient (convective) driving of turbulence and strongly subcritical behaviour. All these issues are discussed in § 8.

4. Small-scale magnetic field generation

In this section we describe current research into small-scale dynamos. One problem that appears to have been solved to the satisfaction of dynamo theorists is the fast dynamo problem.

4.1. One-scale velocity fields and fast dynamo theory

Simply put, fast dynamo theory is concerned with the kinematic generation of magnetic field (on any scale) at high $Rm$ (such as pertains in virtually all astrophysical objects). Consider a velocity field defined at a single scale $\ell _0$ with a characteristic velocity $u_0$. Then, defining $Rm = u_0 \ell _0 /\eta$ in the usual way, the fast dynamo problem is concerned with the behaviour of the growth rate $\sigma (Rm)$ of the dynamo at high $Rm$. There are two possibilities, either $\sigma \rightarrow 0$ as $Rm \rightarrow \infty$ in which case the dynamo is described as ‘slow’. An alternative is that $\sigma \rightarrow \textrm {const.} > 0$ as $Rm \rightarrow \infty$ – a so-called fast dynamo. The two possible options for the growth-rate curve are illustrated in figure 4.

Figure 4. Growth rate $\sigma$ as a function of $Rm$ for a slow dynamo (black curve) and a fast dynamo (red curve).

Moffatt & Proctor (1985) demonstrated that the eigenmodes associated with fast dynamo action may exist, providing that they have a scale of variation $ {O}(Rm^{-1/2})$ as $Rm \rightarrow \infty$, nearly everywhere in the fluid domain. Much of our understanding of the behaviour of fast dynamos arises from the field of dynamical systems and mixing; progress has been made by examining the simpler problem where the flow is modelled as a discontinuous (in time) map. Indeed there are strong parallels between the two problems. I will not go into details here, but the interested reader should consult the excellent Childress & Gilbert (1995).

A central result arising from dynamical systems approaches to fast dynamo action is that which bounds the asymptotic growth rate of the dynamo by the topological entropy of the flow (Klapper & Young 1995); see also Finn & Ott (1988). This is important as it immediately rules out the possibility of fast dynamo action for integrable flows, such as the steady $2 \frac {1}{2}$-D flows discussed earlier. Chaotic particle paths are required for a flow to be a fast dynamo. Chaos may be introduced in two obvious ways, either by making the flow fully three-dimensional or by introducing time dependence to the $2 \frac {1}{2}$-D flows (whilst for simplicity still keeping the flows at a single spatial scale). Moving to three dimensions makes computations at high $Rm$ extremely challenging, owing to the severe constraints imposed by the requirement to resolve structures on the scale $Rm^{-1/2}$ in three dimensions. Nonetheless, progress has recently been made in investigating the dynamo properties of the ABC flow given earlier as

(4.1)\begin{equation} {\boldsymbol u}_{ABC} = (C \sin z + B \cos y, A \sin x + C \cos z, B \sin y + A \cos x). \end{equation}

This flow is chaotic, as shown by the Poincaré sections of the particle paths in Dombre et al. (1986) and the finite-time Lyapunov exponents shown in figure 5(a). Calculations of the growth rate as a function of $Rm$ have periodically been made, since the earliest calculations (see e.g. Arnold & Korkina 1983; Galloway & Frisch 1984) extending to higher $Rm$ as computational power has increased. Figure 5(b) taken from Bouya & Dormy (2015) shows that with current computing facilities $Rm \sim 10^5$ is possible for this flow. The figure shows that even at this value of $Rm$ the dynamo is not in its asymptotic regime.

Figure 5. Finite-time Lyapunov exponents for the $ABC=1$ flow (after Brummell, Cattaneo & Tobias 2001). The ABC $=1$ flow is unusual in having rather large integrable regions and small chaotic regions; Poincaré sections for this flow can be found in Dombre et al. (1986). (b) Growth rate as a function of $Rm$ for the ABC dynamo (after Bouya & Dormy 2015).

A much more promising way to investigate fast dynamo action is to introduce chaos via time dependence in a $2 \frac {1}{2}$-D flow. This of course has the benefit of allowing computations to proceed in two dimensions and so facilitate the investigation of the high $Rm$ regime. This approach was pioneered by Otani (1988, 1993) and Galloway & Proctor (1992), who constructed similar flows. Here, I give details dynamo action in the Galloway–Proctor circularly polarised (GPCP) flow, which takes the form

(4.2)\begin{equation} {\boldsymbol u} = \boldsymbol{\nabla} \times (\psi (x,y,t) \hat{\boldsymbol z}) + \psi (x,y,t)\hat{\boldsymbol z}, \end{equation}

where

(4.3)\begin{equation} \psi = \sin (\,y + \sin t) + \cos (x + \cos t). \end{equation}

This is based on the Roberts flow, with the time dependence being introduced via a rotation of the cellular pattern around a circle. This introduces a significant region of chaos (although regions of integrability still remain) as shown in the finite-time Lyapunov exponents of figure 6).

Figure 6. (a) Lyapunov exponents as a function of starting position in the plane $z=0$ for the GPCP flow. (b) Growth rate as a function of $Rm$ for fixed $k_z = 0.57$ for the GPCP flow. Courtesy of A. Clarke (after Galloway & Proctor 1992).

The magnetic field for the GPCP flow is again generated on the small diffusive length scale $\ell _B \sim Rm^{-1/2}$ as $Rm \rightarrow \infty$. The growth rate (which is a function of vertical wavenumber $k_z$ and $Rm$) has the form shown in figure 6(b) for fixed $k_z$. Interestingly, for large enough $Rm$ the preferred wavenumber becomes independent of $Rm$, as does the optimised growth rate. Notice also that this dynamo reaches its asymptotic growth rate for moderate $Rm$ and so is an example of a quick dynamo (Tobias & Cattaneo 2008a), discussed below. Although it is impossible to prove numerically that these time-dependent flows are fast dynamos, all the evidence certainly points in this direction. It is now widely accepted that sufficiently chaotic flows at a single scale will act as fast dynamos. Of course it is possible to introduce time dependence to 3-D single-scale cellular flows such as the ABC flow, which has the tendency to increase the regions of chaos and hence their dynamo efficiency. Such dynamos have been investigated in both the linear and nonlinear regimes (Brummell et al. 2001).

One aspect of fast dynamo theory that is not widely appreciated is that it is only applicable for high $Pm$ fluids, as noted by Tobias & Cattaneo (2008a). The paradigm of holding the flow fixed and increasing $Rm$ is equivalent to increasing the $Pm$ of the fluid. Clearly increasing $Rm$ whilst holding $Pm$ fixed requires an equivalent increase in $Re$, which will lead to a change in the form of the flow for any realistic forcing mechanism. In particular, increasing $Re$ almost inevitably leads to turbulence with a wide range of spatial and temporal scales.

4.2. Multi-scale velocity fields

In this section we examine kinematic dynamos where the underlying flow has a spectrum of spatial scales. As discussed above, one can think of two cases; at high $Pm$ the magnetic field dissipates at scales much smaller than the smallest eddy and one can rely somewhat on fast dynamo theory based on considering the smallest eddy alone (since the smallest eddy is the one with the fastest turnover time). At low $Pm$ the magnetic field dissipates in the inertial range of the turbulence and life becomes more complicated. We shall start by considering the simplified (some might say over-simplified) case of random velocity fields.

4.2.1. Random dynamos - the Kraichnan–Kazantsev formulation

The discussion developed in this section summarises that of Tobias, Cattaneo & Boldyrev (2012). Turbulent velocity fields have a coherent and random component, both of which may contribute to the dynamo properties. The simplest model of (kinematic) dynamo action driven by a purely random flow on a range of scales is the so-called Kazantsev model (Kazantsev 1968). It is an example of a solvable model for the statistics of the magnetic field based on those for a prescribed velocity, and as such may be viewed as an early example of direct statistical simulation (see § 10).

The prescribed velocity takes the form of a Gaussian, $\delta$-correlated (in time) velocity field, with zero mean and a covariance given by $\langle u_i({\boldsymbol x+r}, t)u_j({\boldsymbol x}, \tau ) \rangle =\kappa _{ij}({\boldsymbol x},{\boldsymbol r})\delta (t-\tau )$. Geophysical and astrophysical flows that lead to dynamo action do have means, are in general anisotropic and inhomogeneous (and this should be a feature of any description of dynamo action). However, analytic progress can be made for the dynamo problem by assuming that the underlying flow is isotropic and homogeneous, in which case the velocity correlation function has the form

(4.4)\begin{equation} \kappa_{ij}({\boldsymbol r})=\kappa_N(r)\left(\delta_{ij}-\frac{r_ir_j}{r^2} \right)+\kappa_L(r)\frac{r_ir_j}{r^2},\end{equation}

where $r = |{\boldsymbol r}|$. In addition, for incompressible velocity fields we have that $\kappa _N=\kappa _L+(r\kappa _L^\prime )/2$, and so the velocity statistics can be characterised by the single scalar function $\kappa _L(r)$. Progress is made by defining a corresponding expression for the magnetic covariance $\langle B_i({\boldsymbol x+r}, t)B_j({\boldsymbol x}, t) \rangle =H_{ij}({\boldsymbol x},{\boldsymbol r}, t),$ where for similar reasons to above

(4.5)\begin{equation} H_{ij}=H_N(r,t)\left(\delta_{ij}-\frac{r_ir_j}{r^2} \right)+H_L(r,t)\frac{r_ir_j}{r^2}.\end{equation}

Similarly, $\boldsymbol {\nabla } \boldsymbol {\cdot } {\boldsymbol B}=0$ gives $H_N=H_L+(rH_L^\prime )/2$, and so the correlator is completely determined by $H_L(r,t)$. The evolution equation for $H_L$, in terms of the input function $\kappa _L(r)$, follows from deriving the equation for the magnetic correlator

(4.6)\begin{equation} {\partial_t H}_L=\kappa H_L^{\prime \prime} + \left(\frac{4}{r}\kappa+\kappa' \right)H_L^{\prime} +\left( \kappa''+\frac{4}{r}\kappa' \right)H_L, \end{equation}

where $\kappa (r)=2\eta +\kappa _L(0)-\kappa _L(r)$ is the ‘renormalised’ velocity correlation function. Remarkably, changing variable via $H_L=\psi (r,t) r^{-2}\kappa (r)^{-1/2}$ leads to a related equation that formally coincides with the Schrödinger equation in imaginary time, i.e.

(4.7)\begin{equation} \partial_t \psi=\kappa(r)\psi''-V(r)\psi.\end{equation}

Here, $\psi$ describes the wave function of a quantum particle with variable mass, $m(r)=1/[2\kappa (r)]$, in a 1-D potential ($r>0$) given by

(4.8)\begin{equation} V(r)=\frac{2}{r^2}\kappa(r)-\frac{1}{2}\kappa''(r)-\frac{2}{r}\kappa'(r)- \frac{(\kappa'(r))^2}{4\kappa(r)}. \end{equation}

Equation (4.7) has been investigated in detail for various choices of prescribed $\kappa (r)$ (see Ruzmaikin, Sokoloff & Zeldovich 1990; Chertkov et al. 1999; Arponen & Horvai 2007). For a more thorough review see Tobias et al. (2012) and Rincon (2019), although the main results are summarised below.

A homogeneous, isotropic turbulent flow with a wide range of scales, and a well-defined inertial range can be characterised by the second-order structure function $\varDelta _2(r) = \langle |({\boldsymbol u}({\boldsymbol x},t) - {\boldsymbol u}({\boldsymbol x}+{\boldsymbol r},t))\boldsymbol{\cdot}{\boldsymbol e}_r|^2 \rangle$, where ${\boldsymbol e}_r={\boldsymbol r}/r$. The inertial and dissipative ranges are then described by the scaling exponents of $\varDelta _2(r)$ with $\varDelta _2(r) \sim r^{2 \alpha }$, where $\alpha$ is termed the roughness exponent of the flow. In the dissipative sub-range we expect $\alpha =1$ as the velocity is a smooth function of position, whilst for the inertial range the velocity is rough and $\alpha < 1$ – for Kolmogorov turbulence $\alpha = 1/3$. If the slope of the energy spectrum in the inertial range is given by $E_k \sim k^{-p}$ then $p$ is related to the roughness exponent by $p = 2 \alpha +1$.

We shall briefly describe dynamo action in the $Pm \gg 1$ case, where the the velocity is smooth on the dissipative scale of the field, and the $Pm \ll 1$ case where the velocity is rough there. For the smooth case, exponentially growing solutions of (4.6) can be found with the magnetic energy spectrum $E_M$ peaked at the magnetic dissipation scale, just as for the single-scale flows considered earlier. The spectrum for the magnetic energy in the range $1/l_\eta < k < 1/l_\nu$ has $E_M \sim k^{3/2}$, independent of the spectral index for the velocity in the inertial range (Kulsrud & Anderson 1992). This regime for a smooth velocity is also referred to as the Batchelor regime, as it corresponds to that of Batchelor (1959) for passive-scalar advection.

The low $Pm$ case is more complicated, and we shall not go into much detail here. Briefly, when $Pm \ll 1$, the magnetic field dissipates in the inertial range, where $\kappa (r) \sim r^{1+\alpha }$ (see e.g. Tobias et al. 2012). Hence in the inertial range the Schrödinger (4.7) has an effective potential with the following properties. At small scales the effective potential is regularised by magnetic diffusion; growing dynamo solutions correspond to bound states for the wave function $\psi$. These are guaranteed to exist when $0 < \alpha < 1$. Hence, dynamo action is always possible even in the case of a rough velocity at low $Pm$ (Boldyrev & Cattaneo 2004). However, it is important to note that the effective potential always remains $\propto 1/r^2$ in the inertial range; its depth decreases as $\alpha \rightarrow 0$. Hence, it is harder to drive dynamo action the rougher the velocity – this has serious consequences for our ability to generate magnetic fields in liquid metal experiments, as discussed in § 9.

Equation (4.7) can be solved asymptotically and solutions demonstrate that the growth time of the dynamo is of the order of the turnover time of the eddies at the diffusive scale ($\ell _\eta$). Moreover, it shows (Boldyrev & Cattaneo 2004) that there is a dramatic increase in the effort (computational or experimental) as the velocity becomes rougher (1+$\alpha$ gets smaller). Hence the effort required to drive a dynamo in a rough velocity also increases.

For these random flows we may now describe the behaviour of the critical magnetic Reynolds number $Rm_c$ as one moves from large $Pm$ to small $Pm$. At high $Pm$ the effort necessary to drive a dynamo is modest. As $Pm$ decreases through unity the viscous scale becomes smaller than the resistive scale and the dynamo begins to operate in the inertial range. There is then an increase in the effort required to sustain dynamo action. However, once the dynamo is in the inertial range, further decreases in $Pm$ do not make any difference to the effort required (as measured by $Rm_c$). We also note here that the Kazantsev model relates to the addition of multiplicative noise in the induction equation. Mathematically this case was considered by, for example, Laval et al. (2006) who were interested in the effects of turbulence on the dynamo threshold. They showed that there are two regimes – an intermittent regime, with a power law distribution for $B$ and at higher $Rm$, a second threshold where the field saturates towards a stable non-zero value.

This picture is largely consistent with numerical models of dynamo action in random flows (Yousef, Brandenburg & Rdiger 2003; Schekochihin et al. 2004, 2005) as we shall see in the next section.

The Kazantsev model as proposed is extremely restrictive, although it can easily be extended to cases where the correlation time of the turbulence is finite (Vainshtein & Kichatinov 1986; Kleeorin, Rogachevskii & Sokoloff 2002), provided the growth time of the dynamo is long compared with this correlation time. Anisotropic versions of the Kazantsev model have also been constructed by Schekochihin et al. (2004b). These solvable models will also prove useful in understanding the generation of organised field (as discussed in § 5.3.4).

4.2.2. Numerical solutions of random dynamos – the low $Pm$ problem

Owing to the importance of understanding how dynamos onset in turbulent flows at low $Pm$ for experiments (see § 9), there has been much numerical effort in this direction. These numerical calculations are fraught with difficulty as extremely large calculations are required to capture the separation of spatial scales. Tobias et al. (2012) calculate that the resolution of a numerical model required to answer the question of the behaviour of the critical value of $Rm$ for dynamo action at low $Pm$ is beyond the reach of current computational resources (requiring a calculation of size $10\ 000^3$ points). However, numerical calculations are beginning to yield some indication of the role of low $Pm$.

In order to relate the Kazantsev models to numerical kinematic numerical simulations, where the Navier–Stokes equations are solved (with no Lorentz force) to provide the input to the induction equation, the Kazantsev models have been extended take into account departures from Gaussianity in the flow, with similar conclusions being drawn as for Gaussian flows. The results are summarised in figure 7, which shows $Rm_c$ as a function of $Re$ for a collection of such calculations (Schekochihin et al. 2007). At large and moderate $Pm$ these results are consistent with the predictions of the Kazantsev model described above, as shown in figure 7. As noted above, calculations rapidly becomes prohibitive at small $Pm$, and so this regime is not really accessible to direct numerical simulation (DNS), unless large-eddy simulations (LES) are utilised (Ponty et al. 2007). However, care must be taken here – in this regime the dynamo growth is controlled by the roughness exponent of the flow and so one would require an LES scheme with exceptional representation of the roughness.

Figure 7. Onset of dynamo action at moderate $Pm$, from Schekochihin et al. (2007). (a) Growth rate of magnetic energy as a function of $Pm$ for five values of $Rm$. (b) Growth/decay rates in the parameter space $(Re, Rm)$. Also shown are the interpolated stability curves $Rm_c$ as a function of $Re$ based on the Laplacian. Hyperviscous runs are shown separately.

Turbulence is, however, significantly more complicated than random flows that are $\delta$-correlated in time, limiting the applicability of the Kazantsev model. For these flows it should be the case that characteristics other than the spectral slope of the velocity (and hence the roughness exponent) control the dynamo growth. We discuss this case in the next section.

4.2.3. Flows with coherence

In this section we consider the case more relevant to geophysical and astrophysical flows, where the flow has two components, a random component as described above and a systematic component whose coherence time is long compared with its turnover time. Such coherent structures, such as long-lived vortices, are ubiquitous in flows where rotation and stratification are important. The structures also exist on a wide range of spatial scales and so it is natural to ask what determines the kinematic dynamo growth rate in a multi-scale flow with coherent structures.

The theory for such flows (which do not now have short correlation times) was developed by Tobias & Cattaneo (2008a,b). This was achieved in two stages, the first step involved demonstrating that for such flows, knowledge of the form of the spectrum is not enough to determine the dynamo properties. They considered a flow with both long-lived vortices and a random component with a well-defined spectrum. They took advantage of the G.O. Roberts trick for dynamos with $2 \frac {1}{2}$-D velocities, by synthesising the flows from solutions to the active scalar equations. They also generated a second (random) flow with the same spectrum as the first by randomising the phases of the spectral components. They found that the flow with coherent long-lived structures is a much more efficient dynamo than the flow that is purely random. Here, by efficient we mean that at the same $Rm$ the dynamo growth rate is higher for the coherent flow. The sustained systematic stretching from the long-lived coherent structures is pivotal in generating field efficiently. In particular, helical vortices are able to generate magnetic field structures – in this case the field generated takes the form of a helix reminiscent of that generated by a Ponomarenko dynamo.

It is now reasonable to assume that a multi-scale kinematic dynamo will be dominated by the coherent parts of the flow – the long-lived structures (rather than the temporally $\delta$-correlated random components) will determine the dynamo growth rate and the form of the field. So, if a dynamo consists of a superposition of such coherent eddies with a range of spatial scales and turnover times (all shorter than their correlation time), what factors determine the small-scale dynamo growth rate?

Progress can be made by assuming that each eddy acts as dynamo largely independently of the eddies at very different scales. This assumption was validated in a model problem by Cattaneo & Tobias (2005). A further step is to assume that each of the dynamo eddies are ‘quick dynamos.’ As defined by Tobias & Cattaneo (2008a), a ‘quick dynamo’ is one that reaches a neighbourhood of its maximal growth rate quickly as a function of $Rm$. Here, ‘quickly’ is a somewhat imprecise term – although a rule of thumb would be for the growth rate to come within $10\,\%$ of its maximum for $Rm \sim 10 Rm_c$; by this definition most dynamos are quick dynamos. If this is the case, then the dynamo is driven by the coherent eddy which has the shortest turnover time $\tau$ and has $Rm \geq { {O}}(1)$. Both the local $Rm$ and $\tau$ are functions of the spatial scale and so the location in the spectrum of the eddy responsible for dynamo action depends on the spectral slope as well as the magnetic Reynolds number on the integral scale.

The difference between dynamos that are dominated by random components and those that have long-lived coherent structures has recently been confirmed by Seshasayanan, Dallas & Alexakis (2017a). They considered the onset of dynamo action in a randomly forced flow subject to the effects of rotation. As the rotation is increased, the flow develops more spatial and temporal coherence with the coherent vortices eventually winning out over the random element of the flow. This has the effect of reducing the critical $Rm$ at which dynamo action occurs, even at low $Pm$ – the coherence engendered by the rotation turns a low $Pm$ dynamo into a high $Pm$ dynamo as predicted by Tobias & Cattaneo (2008a).

5. Organised magnetic field generation

It is often the case that astrophysical magnetic fields display some degree of organisation (or order) either spatially (displaying order on spatial scales large compared with the typical eddies in the turbulent flows) or temporally (displaying temporal coherence on time scales much longer than typical time scales in the turbulence) or sometimes both (for example the spatio-temporal behaviour of the solar cycle).

It is therefore natural to consider theories that describe the evolution of the ‘organised’ components of the magnetic field (and potentially the velocity field). This could be for fields that are organised either spatially or temporally. In order to make progress it is necessary to give some mathematical precision to the concept of an organised field. This brings in the concept of averaging and forces the theorist to turn to statistical theories that are designed to provide information about these average quantities. Such theories are ubiquitous in fluid mechanics (for example Kraichnan 1965; Frisch 1995; Bühler 2014) although, owing to the historical path of research in dynamo theory, the statistical theories in the two disciplines have tended to develop along different tracks. I shall argue strongly that future progress in dynamo theory can only be made by reconciling these approaches and returning to the paradigms preferred by fluid dynamicists. Nonetheless much can be learned from the dynamo approach, which I shall briefly review here.

5.1. The nature of averaging

We shall be concerned with the derivation and solution of equations for the average properties of our state variables (for example ${\boldsymbol B}$ and ${\boldsymbol u}$). We shall primarily be interested in decomposing state variables into mean (average) and fluctuating parts so that for example

(5.1)\begin{equation} {\boldsymbol B} = {\overline{\boldsymbol B}} + {\boldsymbol b}'.\end{equation}

Here, the overbar represents a linear averaging process that (in general) obeys the Reynolds averaging rules, i.e.

(5.2)\begin{equation} \overline{{\boldsymbol B_1} + {\boldsymbol B_2}} = \overline{{\boldsymbol B}_1} + \overline{{\boldsymbol B}_2}, \end{equation}

and

(5.3)\begin{equation} {\overline{\overline{\boldsymbol B}}} = \overline{\boldsymbol{B}}.\end{equation}

Hence, averaging (5.1) gives

(5.4)\begin{equation} \overline{\boldsymbol b'} = 0. \end{equation}

In terms of products it is also convenient, although not necessary, if the averaging procedure satisfies

(5.5)\begin{equation} {\overline{\overline{\boldsymbol B}{\boldsymbol b'}}}=0, \end{equation}

and

(5.6)\begin{equation} {\overline{{\overline{\boldsymbol B_1}}\,{\overline{\boldsymbol B_2}}}}={\overline{\boldsymbol B_1}}\,{\overline{\boldsymbol B_2}}.\end{equation}

There are many different forms of averaging that satisfy the above Reynolds averaging rules, and some useful forms that do not. We briefly comment on a few here.

5.1.1. Spatial and temporal averaging

Perhaps the most utilised form of averaging assumes that the fluctuations are characterised by a length scale $\ell _0$ (perhaps the scale of the energy containing eddies or fluctuating magnetic field). This is assumed small compared with the larger scale $L$ of the variation in the averaged quantities. For a system-scale dynamo $L$ will typically be ${ {O}}(R)$ where $R$ is the length scale of the domain. It is then possible to define an intermediate scale $a$ that satisfies $\ell _0 \ll a \ll L$. The spatial average can then be defined as (Moffatt 1978)

(5.7)\begin{equation} {\overline{\boldsymbol B}} \equiv \langle {\boldsymbol B} \rangle_a \equiv \frac{3}{4 {\rm \pi}a^3}\int_V {\boldsymbol B}({\boldsymbol x}+{\boldsymbol \xi},t){\rm d}^3 {\boldsymbol \xi}, \end{equation}

where $V$ is a sphere of radius $a$. A fairly drastic, although perhaps also very natural, example of spatial averaging is to average over one spatial coordinate. For example the butterfly diagram of figure 2 is constructed by averaging the surface sunspot data over longitude and plotting the averaged field as a function of latitude and time. For this type of averaging there is a natural separation of spatial scales, with the axisymmetric component of field formally being separated in azimuthal wavenumber space from the non-zero wavenumbers. This type of ‘zonal’ averaging is often utilised in planetary and stellar physics where zonal symmetry of the underlying turbulent statistics is a natural assumption (see e.g. Braginskii 1975). In local Cartesian numerical models, we shall see that it is sometimes natural to average over horizontal coordinates.

A similar separation and averaging procedure is available if the quantity to be averaged varies on two time scales. For example if the fluctuations have a characteristic short time scale $t_0$ (perhaps as defined by the correlation time) and the averaged quantities vary on a long time scale $T$ then one can average over an intermediate time scale $\tau$ with $t_0 \ll \tau \ll T$ by defining

(5.8)\begin{equation} {\overline{\boldsymbol B}} \equiv \langle {\boldsymbol B} \rangle_\tau \equiv \frac{1}{2 \tau}\int_{-\tau}^{\tau} {\boldsymbol B}({\boldsymbol x},t+\tau')\, \textrm{d} \tau'. \end{equation}

5.1.2. Ensemble averaging

All of the averaging processes described above rely on some separation of scales (either spatial or temporal) between the average quantities and the fluctuations (this occurs naturally in the case of coordinate averaging). In most turbulent systems, this separation is difficult to achieve, and there is energy at all scales. Whether this really is a difficulty for the theory is open to debate (see e.g. Brandenburg & Subramanian 2005). A more natural averaging procedure, from a statistical viewpoint, is to take an average over realisations of the turbulence – a so-called ensemble average. If the flow is ergodic then this method of averaging should (in combination with the underlying symmetries of the turbulence) yield the same answers as temporal or spatial averaging in the relevant asymptotic limits (Moffatt 1978).

5.1.3. Averaging by filtering and general quasilinear (GQL) averaging

A discussion of averaging procedures would not be complete without including filtering, either spectral or Gaussian (Germano 1992). These filters are chosen to smooth small-scale spatial features and are often representative of the action of observations at finite spatial resolution. It is worth noting that such procedures typically do not satisfy the Reynolds rules of averaging; for example Gaussian filtering does not satisfy (5.3) whilst spectral filtering typically does not satisfy (5.6). However, Gaussian filtering can be used in dynamo calculations (see e.g. Hollins et al. 2017) whilst spectral filtering can be made to satisfy the Reynolds averaging rules by restricting the interactions to remove triad decimation in pairs (Kraichnan 1985). This filtering forms the basis of the generalised quasilinear approximation introduced by Marston, Chini & Tobias (2016).

5.2. Kinematic considerations – mean-field electrodynamics

‘Each success only buys an admission ticket to a more difficult problem’.

Henry Kissinger

Armed with a suitable averaging procedure, it is possible to make significant progress in deriving equations for the evolution of averaged quantities. As before we shall proceed by solving for the evolution of the magnetic field via the averaged induction equation for a prescribed velocity field. This formulation, which is now known as kinematic mean-field electrodynamics has the potential for deep insight; however, I shall argue that it is only via deriving the averaged equations for the full coupled Navier–Stokes/induction system that future progress will be made.

We proceed by splitting both the velocity and magnetic field into averaged and fluctuating parts so that

(5.9a,b)\begin{equation} {\boldsymbol u} = {\overline{\boldsymbol U}} + {\boldsymbol u}', \quad {\boldsymbol B} = {\overline{\boldsymbol B}} + {\boldsymbol b}'. \end{equation}

Now, taking the average of the induction equation (using the Reynolds averaging rules) yields

(5.10)\begin{equation} \frac{\partial \overline{\boldsymbol B}}{\partial t} = \boldsymbol{\nabla} \times \overline{({\boldsymbol u} \times {\boldsymbol B})} +\eta \nabla^2 {\overline{\boldsymbol B}}.\end{equation}

Here, the non-trivial induction term is given (again using the averaging rules) by

(5.11)\begin{equation} {\overline{({\boldsymbol u} \times {\boldsymbol B})}} = \overline{\boldsymbol U}\times \overline{\boldsymbol B} + \overline{{\boldsymbol u'} \times {\boldsymbol b'}}. \end{equation}

So there are two contributions to average induction; the first comes from induction of the average flow and average field (a so-called mean/mean interaction), whilst the second,

(5.12)\begin{equation} {\boldsymbol{\mathcal{E}} = \overline{{\boldsymbol u'} \times {\boldsymbol b'}}}, \end{equation}

arises from the average interaction of the fluctuating velocity and fluctuating magnetic field. This term is commonly known as the turbulent electromotive force (EMF).

Determination of this contribution to the evolution equation for the mean magnetic field forms a large part of the theoretical body of work in mean-field dynamo theory; in rather the same way that modelling the Navier–Stokes Reynolds stress tensor is crucial for many theories for mean-flow/turbulence interactions in hydrodynamic turbulence.

Of course, in principle for a given $\overline{\boldsymbol U}$ and ${\boldsymbol u}'$ (remember we are in the kinematic approximation!), ${\boldsymbol {\mathcal {E}}}$ can be calculated exactly by solving the equation for the fluctuating magnetic field (which is derived by subtracting (5.10) from the full induction (2.2)) to give

(5.13)\begin{equation} \frac{\partial {\boldsymbol b}'}{\partial t}=\boldsymbol{\nabla} \times (\overline{{\boldsymbol U}} \times {\boldsymbol b}')+\boldsymbol{\nabla} \times ({\boldsymbol u}' \times {\overline{\boldsymbol B}})+\boldsymbol{\nabla} \times \boldsymbol{\mathcal{G}} + \eta \nabla^2 {\boldsymbol b}', \end{equation}

where

(5.14)\begin{equation} \boldsymbol{\mathcal{G}} = {\boldsymbol u}' \times {\boldsymbol b}' - \overline{{\boldsymbol u}' \times {\boldsymbol b}'}. \end{equation}

Notice that the fluctuation (5.13) has four terms on the right hand side. The first two of these involve the interaction of means with fluctuations whilst the fourth is a linear term in the fluctuations. The third term given in (5.14) involves the interaction of fluctuations with fluctuations (to yield fluctuations) and may prove to be problematic for theoretical progress. For this reason this has been termed the pain in the neck (PIN) term. In the – perhaps less expressive – hydrodynamic literature fluctuations are often identified with eddies and the corresponding term in the fluctuating Navier–Stokes equation is often termed the eddy/eddy $\rightarrow$ eddy nonlinearity.

5.3. Calculation of the EMF

At this point in the derivation, no approximations have been made and there is little controversy. Indeed, one could solve (5.10) coupled with (5.13) – equivalent to solving the full induction equation. One could argue about whether it makes sense to separate the scales into averages and fluctuations, but if these are treated in the same manner then the point is moot as to whether this is a debate worth having.

However, in order for the approach of mean-field electrodynamics to be useful, the solution of the full problem is to be avoided (and indeed is impossible with current computational resources). Progress can only be made by putting the averages and the fluctuations on a different theoretical footing – and so treating them differently. Solving for the averages is not computationally difficult, and so this can be performed efficiently. However, in order to perform these calculations we need to be able to calculate the turbulent EMF, which arises from the average effect of the fluctuation/fluctuation interaction. In general it would be helpful to have a theory that relates the EMF to the average quantities, either via a differential equation or a simple prescription.

There are two possible ways to proceed. The first is to derive an evolution equation for the (generalised) EMF (which can be thought of as a low-order statistic of the flow). This approach has much in common with some techniques used in hydrodynamic turbulence theory.

The second approach is to milk the formal linearity of the kinematic induction equation for all it is worth, and see where it gets us. The answer is a surprising distance. As noted by Moffatt (1978) ‘$\ldots$ there is now a satisfactory body of theory for the determination of $\boldsymbol {\mathcal {E}}$. The reason for this (comparative) degree of success can be attributed to the linearity in ${\boldsymbol B}$ of the induction equation. There is no counterpart of this linearity in the dynamics of turbulence’. It is important though to appreciate the limitations of this approach.

5.3.1. Deductions from linearity

Clearly, calculation of the EMF is possible from evaluation of ${\boldsymbol b}'$ (recalling that the velocity field is prescribed in the kinematic regime), so the solution of (5.13) is desirable. However, as noted previously, the solution of this equation together with that of the mean equation is equivalent to solving the full problem. It is, therefore, useful to see what can be learned from the structure of (5.13). This equation is formally a linear equation for ${\boldsymbol b}'$, with a forcing term given by $\boldsymbol {\nabla } \times ({\boldsymbol u}' \times {\overline{\boldsymbol B}})$. Solutions to this equation for ${\boldsymbol b}'$ are linearly (although not homogeneously) related to $\overline{\boldsymbol B}$; that is

(5.15)\begin{equation} {\boldsymbol b}' = {\boldsymbol b}'_{ss}+{\mathcal{L}}(\overline{\boldsymbol B}). \end{equation}

Here, ${\mathcal {L}}$ is a linear operator and ${\boldsymbol b}'_{ss}$ is that fluctuation field that exists in the absence of a mean field $\overline{\boldsymbol B}$. At high $Rm$ this ‘small-scale dynamo field’ is inescapable as we have discussed at length above. However, it is not clear that the presence of this field poses any difficulty for the calculation of the EMF. Clearly ${\boldsymbol b}'_{ss}$ is, by definition, that field that exists in the absence of any large-scale field $\overline{\boldsymbol B}$. It seems unlikely therefore that, except in extremely contrived situations, it can contribute to an EMF at high $Rm$; since its contribution to the EMF would inevitably drive a large-scale field, which is supposed absent. In the nonlinear regime, this argument clearly does not hold (as the small-scale dynamo subspace may saturate and become unstable to large-scale perturbations as we shall discuss later). Indeed some mechanisms for large-scale field generation rely on interactions arising nonlinearly from a saturated small-scale dynamo (Squire & Bhattacharjee 2016).

Hence, it seems that, although the fluctuation dynamo field in the kinematic regime is not linearly and homogeneously related to $\overline{\boldsymbol B}$, the EMF is. This is not to say that the presence of small-scale dynamo action is irrelevant; indeed, as we shall see in § 5.4.2, in the kinematic regime at high $Rm$ the growth rate of the dynamo (at both large and small scales) is completely controlled by the presence of a small-scale dynamo (Nigro et al. 2017).

The linear dependence of the EMF on the averaged field suggests an integral representation of the form

(5.16)\begin{equation} {\mathcal{E}}_i({\boldsymbol x},t) = \int \int K_{ij}({\boldsymbol x},t;{\boldsymbol \xi},\tau) \overline{{B}_j}({\boldsymbol x}+{\boldsymbol \xi},t+\tau) \, {\rm d}^3 {\boldsymbol \xi}\, {\rm d} \tau, \end{equation}

for some kernel $K_{ij}$.

If a separation of spatial scales (for example) also pertains then it is possible to expand $\overline{\boldsymbol B}$ in a Taylor series, i.e.

(5.17)\begin{equation} \overline{{B_i}} ({\boldsymbol x+{\boldsymbol \xi}}) = {\overline{B_i}} ({\boldsymbol x}) + \xi_j \frac{\partial {\overline{B_i}}}{\partial x_j}({\boldsymbol x}) + \frac{1}{2} \xi_j \xi_k \,\frac{\partial^2 {\overline{B_i}} }{\partial x_j \partial x_k}({\boldsymbol x}) + \cdots. \end{equation}

For scale separation, $|{{\boldsymbol \xi }}|$ is small compared with a typical scale of $\overline{\boldsymbol B}$, and higher-order terms may be neglected. As noted by Moffatt (1978) terms in the Taylor expansion involving time derivatives of the mean field can be written in terms of the mean field and spatial derivatives through back substitution into (5.10).

It is therefore natural, having made this approximation, to postulate a series expansion of the EMF in terms of spatial gradients of the average field of the form

(5.18)\begin{equation} \mathcal{E}_i = \alpha_{ij} {\overline{B_j}} + \beta_{ijk} \frac{\partial {\overline{B_j}}}{\partial x_k}+\cdots ,\end{equation}

where the coefficients $\alpha _{ij}$ and $\beta _{ijk}$ are known as turbulent transport coefficients.

In general, $\alpha _{ij}$ and $\beta _{ijk}$ need to be calculated from the turbulent flow statistics (usually this requires some approximations as we shall see). However, some general statements about their form can be made immediately from (5.18).

We note that $\boldsymbol {\mathcal {E}}$ is a polar vector, whereas $\overline{\boldsymbol B}$ is an axial vector. From this we can immediately deduce that $\alpha _{ij}$ is a pseudo-tensor, which can be decomposed into symmetric and anti-symmetric parts as

(5.19)\begin{equation} \alpha_{ij} = \alpha_{ij}^s - \epsilon_{ijk} \gamma_k. \end{equation}

Hence, the ‘$\alpha$-term’ in (5.18) takes the form $\alpha _{ij}^s{}{\overline{B_j}} + ({\boldsymbol \gamma } \times \overline{{\boldsymbol B}})_i$. The anti-symmetric part of the $\alpha$-tensor therefore contributes an extra mean velocity to the mean-field equations. The symmetric part of the $\alpha$-tensor obviously depends on the flow. In a system with preferred directions, given by say gravity ${\boldsymbol g}$ and rotation ${\boldsymbol \varOmega }$, it is given by Krause & Raedler (1980)

(5.20)\begin{equation} \alpha_{ij} = \alpha_1 \delta_{ij} \hat{\boldsymbol g} \boldsymbol{\cdot} \hat{{\boldsymbol \varOmega}} + \alpha_2 \hat{g}_i \hat{\varOmega}_j + \alpha_3 \hat{g}_j \hat{\varOmega}_i. \end{equation}

This form is useful; however, more enlightening is to consider the case where the turbulence is statistically isotropic (although still not reflectionally symmetric). In this case

(5.21)\begin{equation} \alpha_{ij} = \alpha \delta_{ij}, \end{equation}

with ${\boldsymbol \gamma }=0$. The constant $\alpha$ is a pseudo-scalar, which can only be non-zero if the turbulence lacks reflexional symmetry.

Similar considerations can lead to progress in understanding the second term in the expansion. For the simple case of isotropic turbulence, the $\beta _{ijk}$-tensor is also isotropic and so

(5.22)\begin{equation} \beta_{ijk} = \beta \epsilon_{ijk}. \end{equation}

Here, $\beta$ is a pure scalar (i.e. it may be non-zero even for turbulence with no broken symmetry).

I stress here that, in general, these tensors (even if deemed to be useful) will not be homogeneous or isotropic. The presence of rotation, stratification and mean flows tend to give a preferred direction to the turbulence and the tensors will not take a particularly simple or enlightening form. There is now a significant body of work ascribing importance to the generation of an EMF via particular parts of the tensors arising from different physical interactions. I shall argue later that what really matters are the low-order statistics of the MHD turbulence (in this case the EMF itself). Little progress can be made by ascribing and naming effects.

However, it is still worth following the theory through to its logical conclusion; it remains to determine the transport coefficients $\alpha _{ij}$ and $\beta _{ijk}$, ideally in terms of the prescribed flow. There have been many attempts at this that we shall summarise in § 5.3.3. Before doing this, we shall briefly examine why this theory has proven so popular by demonstrating solutions to the kinematic mean-field equations.

5.3.2. Solution of the mean-field (filtered) equations

With the assumption that the transport coefficients are isotropic, so that $\alpha _{ij} = \alpha \delta _{ij}$ and $\beta _{ijk} = \beta \epsilon _{ijk}$ the mean-field dynamo equations take the form

(5.23)\begin{equation} \frac{\partial \overline{\boldsymbol B}}{\partial t} = \boldsymbol{\nabla} \times (\overline{\boldsymbol U} \times \overline{\boldsymbol B}) +\boldsymbol{\nabla} \times (\alpha \overline{\boldsymbol B}) - \boldsymbol{\nabla} \times (\beta \boldsymbol{\nabla} \times \overline{\boldsymbol B}) +\eta \nabla^2 {\overline{\boldsymbol B}}. \end{equation}

If $\beta$ is constant, $\boldsymbol {\nabla } \times (\beta \boldsymbol {\nabla } \times \overline{\boldsymbol B}) = - \beta \nabla ^2 {\overline{\boldsymbol B}}$ so the equation becomes

(5.24)\begin{equation} \frac{\partial \overline{\boldsymbol B}}{\partial t} = \boldsymbol{\nabla} \times (\overline{\boldsymbol U} \times \overline{\boldsymbol B}) +\boldsymbol{\nabla} \times (\alpha \overline{\boldsymbol B}) +(\eta+\beta) \nabla^2 {\overline{\boldsymbol B}}. \end{equation}

We can now see the role of the transport coefficients clearly for this choice of isotropic turbulence. The $\boldsymbol {\nabla } \times (\alpha \overline{\boldsymbol B})$ term is an inductive term leading to the generation of magnetic field (we shall spell this out in the next section) whilst for this simple configuration the $\beta$ term acts like an enhanced diffusivity. We note that other components of the $\beta$ tensor can act so as to generate field for non-isotropic flows. Given that there are 27 of these coefficients, it is fair to say that the role of them all has yet to be established. It is not even clear that it is desirable so to do, as argued above.

The adoption of the mean-field dynamo equations circumvents Cowling's theorem! All the non-axisymmetric parts of the dynamo are encoded in the small scales, which are not solved for. Hence axisymmetric solutions for the mean fields are allowed (and informative). Combining (5.24) with the axisymmetric formalism equations given in (2.26)–(2.27) one derives the axisymmetric mean-field dynamo equations

(5.25)\begin{gather} \frac{\partial A}{\partial t} + \frac{1}{s} ({\boldsymbol u}_P \boldsymbol{\cdot} \boldsymbol{\nabla} ) (sA) = \alpha B + \eta_T \left(\nabla^2 - \frac{1}{s^2}\right) A, \end{gather}

(5.26)\begin{gather}\frac{\partial B}{\partial t} + s ({\boldsymbol u}_P \boldsymbol{\cdot} \boldsymbol{\nabla} ) \left(\frac{B}{s}\right) = \boldsymbol{\nabla} \times (\alpha{\boldsymbol B}_P) \boldsymbol{\cdot} \hat{{\boldsymbol \phi}} + \eta_T \left(\nabla^2 - \frac{1}{s^2}\right) B + s{{\boldsymbol B}_P} \boldsymbol{\cdot} \boldsymbol{\nabla} \varOmega, \end{gather}

where now it is understood that $A$ and $B$ represent the poloidal and toroidal parts of the mean magnetic field fields and $\eta _T = \eta + \beta$.

Cowling's antidynamo theorem therefore is defeated by the presence of a source term $\alpha B$ in the equation for the poloidal field. This source term, which is now called the $\alpha$-effect had been previously derived by Parker (1955). In the equation for the evolution of the azimuthal field both the shear (the $\omega$-effect) and the $\alpha$-effect can operate to generate this component of the field from the meridional field. Depending on the relative importance of the $\alpha$ and $\omega$ effects, mean-field dynamos may be termed $\alpha ^2$, $\alpha \omega$ or $\alpha ^2 \omega$ dynamos. There are two ways of generating azimuthal field $B$ from poloidal field ${\boldsymbol B}_P$: the $\alpha$-effect or the $\omega$-effect. If the $\alpha$-effect dominates, the dynamo is called an $\alpha ^2$-dynamo. If the $\omega$-effect dominates it is an $\alpha \omega$ dynamo. We can also have $\alpha ^2 \omega$ dynamos where both mechanisms operate equally.

The importance and utility of the mean-field ansatz can be seen by analysing the kinematic mean-field dynamo solutions in a 2-D local Cartesian geometry. We seek solutions independent of $y$, with no meridional flow and constant shear, $U'$, in the $\alpha \omega$ limit (i.e. we set $\boldsymbol {\nabla } \times (\alpha {\boldsymbol B}_P) \boldsymbol {\cdot } \hat {{\boldsymbol \phi }} = 0$), so that

(5.27a,b)\begin{gather} {\boldsymbol B} = \left(- \frac{\partial A}{ \partial z}, B, \frac{\partial A}{ \partial x} \right) , \quad {\boldsymbol u} = (0, U' z, 0), \end{gather}

(5.28a,b)\begin{gather}\frac{\partial A}{\partial t} = \alpha B + \eta_T \nabla^2 A, \quad \frac{\partial B}{\partial t} = J(A,U'z) + \eta_T \nabla^2 B, \end{gather}

where $J(\,f,g) = \,f_x g_z - g_x \,f_z$. We seek wave-like solutions for $A$ and $B$ proportional $\exp (\sigma t + i {\boldsymbol k} \boldsymbol {\cdot } {\boldsymbol x})$, with ${\boldsymbol k}= (k_x,0,k_z)$ to yield a dispersion relation of the form

(5.29)\begin{equation} (\sigma + \eta_T k^2)^2 = {\rm i} k_x \alpha U', \end{equation}

with $k^2 = k_x^2+k_z^2$, and so

(5.30)\begin{equation} \sigma = \frac{1+{\rm i}}{\sqrt{2}} (\alpha U' k_x)^{1/2} - \eta_T k^2. \end{equation}

In a finite Cartesian domain of size $d$ in the $z$-direction, one sets $k_z= {\rm \pi}/ d$ and large-scale dynamo action will occur if the dimensionless dynamo number $D = \alpha U' d^3/ \eta _T^2$ exceeds a threshold. The bifurcation to large-scale dynamo action occurs in a Hopf bifurcation to travelling waves. If $D > 0$ the waves travel in the negative $x$-direction, whilst if $D < 0$ they travel in the positive $x$-direction.

At this point it is immediately clear why mean-field theory has proved so popular. For an axisymmetric (indeed 1-D) model the theory can yield travelling wave solutions reminiscent of the waves of sunspot activity seen in the solar observations. The theory can easily be extended to bounded Cartesian domains – although care must be taken in carrying out the instability calculations owing to the difference between convective and absolute instabilities (Bassom & Gilbert 1997; Tobias, Proctor & Knobloch 1997), and to spherical domains (Roberts 1972b).

5.3.3. Calculations of kinematic transport coefficients: theory and computation

We described above how solution of the mean-field equations can lead to large-scale magnetic fields for simple choices of the transport coefficients. We also described in § 5.3.1 that basic properties of the transport coefficients ($\alpha _{ij}$ and $\beta _{ijk}$) can be determined via consideration of symmetries (and on making the assumptions of homogeneity).

However, it is important to calculate these transport coefficients for solvable and realistic models of turbulence (both in the kinematic and dynamic regime). We begin by describing the kinematic evaluation of these coefficients – this is where the theory gets controversial, as all the current analytical calculations of transport coefficients necessarily involve some form of approximation.

First-order smoothing:

The simplest ansatz in which to calculate the transport coefficients in the kinematic regime invokes the quasilinear approximation. This involves the discard of the fluctuation/fluctuation $\rightarrow$ fluctuation interactions in the equation for the fluctuating field; that is we set $\boldsymbol {\mathcal {G}}=0$ in (5.13) to yield

(5.31)\begin{equation} \frac{\partial {\boldsymbol b}'}{\partial t}=\boldsymbol{\nabla} \times (\overline{{\boldsymbol U}} \times {\boldsymbol b}')+\boldsymbol{\nabla} \times ({\boldsymbol u}' \times {\overline{\boldsymbol B}})+\eta \nabla^2 {\boldsymbol b}'. \end{equation}

The discard of the nonlinear term $\boldsymbol {\mathcal {G}}$ can be understood as assuming that the Kubo number (to be defined below) of the MHD turbulence is small (in the sense that the fluctuation/fluctuation interaction is small compared with the mean/fluctuation interaction). This may occur in many circumstances of interest, although the conditions for which the approximation is valid is difficult to determine a priori. Formally, however, this is certainly the case when the fluctuating field is small compared with the mean field i.e. $|{\boldsymbol b}'| \ll |{\overline{\boldsymbol B}}|$. It is also applicable when $Rm \ll 1$ in which case one may also neglect the $\partial {\boldsymbol b}'/\partial t$ term to leave

(5.32)\begin{equation} \eta \nabla^2 {\boldsymbol b}'=-\boldsymbol{\nabla} \times (\overline{{\boldsymbol U}} \times {\boldsymbol b}')-\boldsymbol{\nabla} \times ({\boldsymbol u}' \times {\overline{\boldsymbol B}}), \end{equation}

which can then easily be solved using standard linear transform techniques (see Moffatt (1978) for a detailed description) to yield ${\boldsymbol b}'$ and hence, as ${\boldsymbol u}'$ is known, $\boldsymbol {\mathcal {E}}$.

The nonlinear term $\boldsymbol {\mathcal {G}}$ can also be formally neglected a priori if the correlation time $\tau _c$ of the turbulence (which is given in the kinematic regime) is sufficiently short. In particular if the Kubo number $u_{rms} \tau _c/\ell \ll 1$ then (5.31) is justified. If in addition $Rm \gg 1$ then

(5.33)\begin{equation} \frac{\partial {\boldsymbol b}'}{\partial t}\approx\boldsymbol{\nabla} \times (\overline{{\boldsymbol U}} \times {\boldsymbol b}')+\boldsymbol{\nabla} \times ({\boldsymbol u}' \times {\overline{\boldsymbol B}}). \end{equation}

Again, this linear system may be solved by transform methods. In the simplified case with no local mean flow one may write

(5.34)\begin{equation} {\boldsymbol{\mathcal{E}}}(t) = \overline{{\boldsymbol u}(t) \times \int_0^t \boldsymbol{\nabla} \times ( {\boldsymbol u}(t') \times \overline{{\boldsymbol B}(t')} ) \,\textrm{d} t'}. \end{equation}

In a statistically steady state and for isotropic turbulence this simplifies to give

(5.35)\begin{equation} {\boldsymbol{\mathcal{E}}}(t) = \int_0^t (\hat{\alpha}(t-t') \overline{{\boldsymbol B}(t')} - \hat{\eta}_T(t-t') \overline{\boldsymbol{\nabla} \times {\boldsymbol B}(t')})\, \textrm{d} t'. \end{equation}

where $\hat {\alpha }(t-t') = -\frac {1}{3}\overline {{\boldsymbol u}'(t) \boldsymbol {\cdot } {\boldsymbol \omega }'(t')}$ and $\hat {\eta }_T(t-t')= \frac {1}{3}\overline {{\boldsymbol u}'(t) \boldsymbol {\cdot } {\boldsymbol u}'(t')}$.

For a mean field that is sufficiently slowly varying compared with the turbulence the integral kernels can be approximated by $\delta$-functions and so

(5.36)\begin{equation} {\boldsymbol{\mathcal{E}}}(t) = \alpha \overline{{\boldsymbol B}} - \hat{\eta}_T \overline{\boldsymbol{\nabla} \times {\boldsymbol B}}, \end{equation}

with

(5.37)\begin{gather} \alpha = - \tfrac{1}{3} \int_0^t \overline{{\boldsymbol u}'(t) \boldsymbol{\cdot} {\boldsymbol \omega}'(t')} \, \textrm{d} t' \approx - \tfrac{1}{3} \tau_c \overline{{\boldsymbol u}' \boldsymbol{\cdot} {\boldsymbol \omega}'}, \end{gather}

(5.38)\begin{gather}\eta_T = \tfrac{1}{3} \int_0^t \overline{{\boldsymbol u}'(t) \boldsymbol{\cdot} {\boldsymbol u}'(t')} \, \textrm{d} t' \approx \tfrac{1}{3} \tau_c \overline{{\boldsymbol u}' \boldsymbol{\cdot} {\boldsymbol u}'}. \end{gather}

First-order smoothing, which is applied within the kinematic framework, is designed so that the EMF may be related to the mean magnetic field through properties of the prescribed velocity field. The smoothing is such that the fluctuating magnetic field is ‘slaved’ to the velocity field and the mean magnetic field and so the EMF can be readily calculated. At high $Rm$, for correlation times that are not short, this slaving is not maintained and the magnetic fluctuations rapidly become decorrelated with the velocity field. This has two major effects. The first is that the EMF cannot be evaluated explicitly and one must resort to numerical procedures for its calculation (as described in the next section). The second, and more severe, effect is that the loss of correlation can lead to a diminution of the EMF relative to the fluctuations. This may have serious consequences for large-scale dynamo theory in the kinematic regime at high $Rm$ as described below – although, as we shall see, in section there will be some effects that can mitigate the loss of correlation and to some extent rescue the theory.

Numerical evaluation of the transport coefficients:

With the growth of computing power, it is now possible to evaluate the transport coefficients numerically, albeit at moderate $Rm$ for 3-D flows. This is, of course, possible both when the magnetic field is weak and does not act back on the flow and therefore the transport coefficients are in the kinematic regime, or in the fully nonlinear regime so that the effects of the Lorentz force on transport coefficients can be calculated – an important topic discussed in § 6.5.

Constructing a reliable method for the calculation of the transport coefficients is not trivial (especially when $Rm$ is not small) since it involves the calculation of averages and care must be taken to ensure that these averages are converged. Calculation of the $\alpha$-coefficient is conceptually straightforward. Consider (5.12) and (5.18). These are two separate expressions for the EMF. In the case where the averages are spatial, the flow and magnetic field can be decomposed into a spatial mean (in this case a uniform field) and spatial fluctuations. Now a uniform field has no spatial derivatives and so (5.18) simplifies to

(5.39)\begin{equation} \mathcal{E}_i = \alpha_{ij} {\overline{B_j}}. \end{equation}

Hence, it is possible to impose constant mean fields with different orientations and calculate $\mathcal {E}_i = \overline {{\boldsymbol u'} \times {\boldsymbol b'}}$ numerically to identify the different components of the $\alpha$-tensor. For a very turbulent flow at high(ish) $Rm$ it turns out that $\overline {{\boldsymbol u'} \times {\boldsymbol b'}}(t)$ in the kinematic regime may be a rapidly varying timeseries which takes a broad distribution of values with a large variance and a relatively small mean (see e.g. Cattaneo & Hughes (2006) who computed the alpha coefficient in a plane layer of rotating Boussinesq convection). This is an indication that the systematic behaviour in the kinematic regime may be in competition with (and perhaps dominated by) the disorganised or fluctuating dynamo behaviour (see § 5.4 for a more complete discussion of this); however, it is fair to say that the relative effectiveness of the fluctuating and organised magnetic fields may depend on both the form of the spectrum of the velocity and the value of $Rm$ as demonstrated by Cattaneo & Tobias (2014) – at very high $Rm$ the dynamo is controlled by smaller eddies from further down the spectrum with shorter correlation times, which means that averaging is a more effective procedure and the distribution of the EMF is narrowed. Thus, in addition to spatial averaging, temporal (or ensemble) averaging should also be utilised to ensure well converged statistics.

Computing the turbulent diffusivity $\beta _{ijk}$ is somewhat more problematic since that relates the EMF to gradients of the organised magnetic field. This causes problems as there is no simple method of imposing a magnetic field with a spatial gradient that is not prone to grow exponentially via dynamo action at high $Rm$ in the kinematic regime. Three methods for evaluating the turbulent diffusivity have been suggested; the ‘test field’ method, the ‘turbulent Ångström’ method and the ‘Lagrangian statistics’ method.

The test field method, introduced by Schrinner et al. (2005) involves the solution of an auxiliary equation and is a simple extension of the procedure outlined above for the calculation of the $\alpha$-effect. Here the aim is to calculate the EMF for a variety of applied test magnetic fields with different orientations and field gradients. If the test fields are chosen to be constants then the $\alpha$-tensor can be recovered in a similar manner to that described above. In principle the $\beta$-coefficients can be backed out by imposing a series of spatially varying test fields and solving for the fluctuations, and hence the EMF. The test field method has been criticised (Cattaneo & Hughes 2009) in that it requires a prohibitively large computational cost to evaluate the transport coefficients at high $Rm$. Nonetheless, at the moderate $Rm$ values that are currently achievable for current numerical computations convergence of results does appear achievable if ensemble averages are used in addition to spatial averages.

The turbulent Ångström method (Tobias & Cattaneo 2013a) is based on that used to measure the thermal conductivity of solids. The method consists of applying an oscillatory source of magnetic field (usually via an oscillating current) with a given frequency larger than the largest natural frequency of the turbulence. Here the turbulent diffusivity can be calculated from the response measured at the frequency of the imposed oscillating field. This method is attractive as the turbulent diffusivity can be calculated to any required precision, simply by analysing a long enough timeseries for the magnetic field, although this may be computationally expensive.

The final method, which involves calculating Lagrangian statistics, is formally valid in the infinite $Rm$ limit and is outlined in Moffatt (1978). Briefly, both the $\alpha$ and $\beta$ tensors may be evaluated numerically by calculating averages over Lagrangian trajectories. This method can (and has) been extended to the finite $Rm$ case by adding stochastic fluctuations to the Lagrangian trajectories (Drummond & Horgan 1986). This elegant method appears to merit further investigation.

5.3.4. Organised Kazantsev dynamos

We have described above how organised magnetic fields may be expected to emerge for flows that lack reflectional symmetry, and how, under certain assumptions, the evolution of the large-scale field may be determined in the kinematic regime by turbulent transport coefficients that encode the eddy/eddy $\rightarrow$ mean interactions. Even earlier, in § 4.2.1, we also discussed how the Kazantsev model for random (short correlation time) flows could yield insight into the onset of small-scale dynamos.

It turns out that the Kazantsev model can be extended to random flows that lack reflectional symmetry; in this case the velocity and magnetic correlators are defined using two functions – one, as before, related to the energy density (either kinetic or magnetic), the other to the helicity (either kinetic or magnetic) (Vainshtein & Kichatinov 1986). The two parts of the magnetic correlator are then described by a pair of coupled Schrödinger-like equations (Vainshtein & Kichatinov 1986; Berger & Rosner 1995). The analysis is now more complicated, yet it is possible to demonstrate that the system remains self-adjoint (Boldyrev, Cattaneo & Rosner 2005). These authors demonstrated that, for sufficient kinetic helicity, extended states are found that do not decay exponentially at infinity. These are the ‘large-scale dynamo solutions’, analogous to the mean-field solutions. It can, however, be shown that at large $Rm$ the largest growth rates remain associated with the localised bound states, so that the overall dynamo growth rate remains controlled by small-scale dynamo action (Malyshkin & Boldyrev 2008a,b). This competition between large and small-scale dynamos in the kinematic regime will be the theme of our discussion in the next few sections.

5.4. The competition between kinematic large- and small-scale dynamos

The discussion of organised random dynamos in the last section described how the growth rate of magnetic fields is determined by the small-scale dynamo, even if there is a large-scale component. In that simple case, the growth-rate for bound (small-scale) states can be shown to be larger than those for extended (large-scale) solutions.

In a non-random flow, the magnetic field must emerge as an eigenfunction to the kinematic dynamo problem with a well-defined (average) growth rate $\sigma$; hence large and small scales grow at this same rate. It is then appropriate to ask which processes control the growth rate of a dynamo containing large and small scales? Moreover, if this is controlled by the small-scale dynamo action, then what is the role of large-scale (or mean-field) dynamo theory in determining the behaviour of the large scales.

The first of these questions has been addressed over a number of years by a variety of authors. The theoretical considerations of § 4 show that in a flow with coherent structures over spectrum of scales, kinematically the dynamo growth is dominated by dynamo action driven by the dynamo eddy which is just supercritical and has the fastest turnover time. As $Rm$ is increased, this dynamo scale moves further and further down the spectrum generating field at ever faster rates. Moreover, the field that is generated is dominated by energy at scales smaller than the characteristic scale of the eddy. These considerations are all borne out by numerical calculations of kinematic dynamos at high $Rm$ (see e.g. Tobias & Cattaneo 2008a, 2013b). The dynamo growth-rate in the kinematic regime is therefore set by the small-scale dynamo; fundamentally, chaotic stretching on the advective time scale has a tendency to win out over an EMF induced by correlations. In this kinematic regime there is nothing really to stop the generation of exceptionally small scales except for the action of diffusion on the resistive scales. Of course as all scales grow at the same rate this tendency for field to be generated on small scales manifests itself in the form of the eigenfunction, which is dominated by fields close to the resistive scale (all other things being equal).

The answer to the second question, that of the utility of kinematic mean field theory in describing the dynamics of the unfiltered equations, we shall defer to later (§ 5.4.2). We shall first discuss how dynamo theorists have approached the problem of the competition between large and small-scale dynamos in the kinematic regime, and whether there are any strategies for promoting the generation of organised field over the small-scale, random fields in the kinematic regime.

5.4.1. The role of shear: suppressing the small-scale dynamo

Most efforts to promote the large-scale dynamo over the small-scale dynamo involve the addition of a large-scale shear flow. As we have described above, in the absence of shear, the kinematic growth rate is determined by the stretching of the small-scale dynamo. This effect also manifests itself in the distribution of the EMF, which in certain circumstances can be calculated to have a small mean and a large variance caused by the fluctuations. If one is to promote the generation of organised magnetic fields over that of small-scale fields then two strategies immediately suggest themselves; either utilise the shear to boost the mean EMF by introducing new effects (i.e. new terms in the transport coefficients) or utilise the shear to suppress the small-scale dynamo.

Much attention has focussed on the first of these strategies (see e.g. Yousef et al. 2008; Käpylä & Brandenburg 2009; Sridhar & Singh 2010 among others). Indeed it can be shown that the presence of a prescribed large-scale shear flow can lead to extra terms in the transport coefficients and hence in the EMF. Moreover a shear flow can in principle interact with a weak and highly fluctuating EMF to produce large-scale field (Richardson & Proctor 2010). All of these mechanisms can be shown to boost the generation of a mean EMF, in certain circumstances. However, at high $Rm$, if this large-scale field is ever to be seen in the kinematic regime, there must be some mechanism present to suppress the small scales which have a tendency to dominate. It would also help to differentiate the two dynamos if the large-scale dynamo were to exhibit a different dynamics than that of the small scales.

A rather different, and potentially more promising, way to proceed is to recognise that, as indicated by fast dynamo theory, the small-scale dynamo relies on the presence of chaotic stretching to survive at high $Rm$. In order to give the large-scale dynamo a chance (Courvoisier & Kim 2009), it is imperative to suppress the small-scale dynamo by reducing the chaotic stretching in the flow. Such a strategy was proposed by Tobias & Cattaneo (2013b). They argued that the addition of a shear flow to a chaotic (or turbulent) flow will increase the regions of integrability in the flow, thus reducing chaos and hence the effectiveness of the small-scale dynamo at high $Rm$. If such a dynamo also has broken reflectional symmetry and is capable of generating an EMF and hence large-scale fields then it is possible that the suppression of the small-scale dynamo will lead to the emergence of large-scale fields. Tobias & Cattaneo (2013b) and Cattaneo & Tobias (2014) investigated the properties of just such a dynamo. They considered a velocity field in a Cartesian domain that was comprised of a combination of generalised Galloway–Proctor flows imposed at a range of scales together with a large-scale shear flow. In particular they considered a multi-scale generalisation of the cellular Galloway–Proctor flow that is maximally helical and takes the form of an infinite array of clockwise and anti-clockwise rotating helices. The pattern of helices itself rotates in a circle with a scale-dependent frequency and can be decorrelated on a time scale $\tau _d$. The decorrelation time of the flow is also scale-dependent and scales in the same manner as the turnover time. Hence the magnetic Reynolds number is also a function of scale and was chosen to decrease as the spatial scale decreases (see Tobias & Cattaneo 2013b for details).

To this time-dependent helical multi-scale flow, a steady unidirectional large-scale shear of the form

(5.40)\begin{equation} {\boldsymbol u}_s = (V_0 \sin y, 0, 0),\end{equation}

was added and dynamo solutions were sought (in a similar manner to the calculations for the Galloway–Proctor flow of § 4.1) at high $Rm$.

In the absence of shear, the flow is an excellent small-scale dynamo. For $\chi = Rm/Rm_c \gg 1$ all the magnetic energy is, as expected, concentrated at small scales and generated on a typical time scale of a Galloway–Proctor eddy. There is no sign of a large-scale organised component of the magnetic field. The situation changes significantly when shear is added. Although the small-scale dynamo is still operational, its growth-rate is reduced (see Tobias & Cattaneo 2013b). At high $Rm$ in the kinematic regime the role of an imposed shear is always to lower the growth rate of a small-scale dynamo. As explained above, this is due to the introduction of regions of integrability to the otherwise chaotic flow. This statement has caused some controversy in the literature, but is correct at high enough $Rm$ since the operation of such a dynamo relies on the presence of chaos in the flow. Dynamos for which the shear increases the dynamo growth rate are either poor small-scale dynamos to start with, or are simply not being investigated at high enough $\chi = Rm/Rm_c$.

Moreover, once the small-scale dynamo is suppressed, it is reasonable to require that the large-scale dynamo can be detected. Figure 8(a), which shows $\overline {B_x}(\,y,t)$ the average (in $x$) of the toroidal (shear-aligned) field at a given $z$-level, demonstrates that this expectation is met. The average field can be seen to take the form of systematic dynamo waves propagating in the regions of strong shear (of the type envisaged by Parker). In the absence of helicity of the small-scale flow, the shear still amplifies fields in bands but there is no systematic generation of spatio-temporally coherent field. That the shear reduces the fluctuations is clear on examination of the probability distribution (as shown by the p.d.f.) for the EMF in figure 8(b). This shows the distribution of EMF as measured by imposing a mean field. In the absence of shear the distribution is dominated by the fluctuations with a large variance and a small mean (as shown by the black curve). As the strength of the shear is increased the variance of the distribution decreases markedly, although the mean shifts barely at all. This demonstrates (for this flow at least) that the crucial role of shear at high $Rm$ in the kinematic regime is to suppress the fluctuations, rather than modify the mean.

Figure 8. (a) Space–time $(t,y)$ plot of average magnetic field $\overline {B_x}(\,y,t)$ at fixed $z=0$ for dynamo with shear $(V_0=5)$ and small-scale flow with net helicity (after Cattaneo & Tobias 2014); here $0 \le y < 2 {\rm \pi}$ whilst $0 \le t \le 40$. (b) Probability density function of EMF for a fixed small-scale helical flow and a range of strength of shear flow (defined in (5.40)) $V_0 = 0$ (black), $V_0 = 1$ (red), $V_0 = 2$ (yellow), $V_0 = 5$ (green), $V_0 = 10$ (cyan), and $V_0 = 20$ (blue). As the shear is increased the distribution narrows significantly with a small change in the mean EMF (after Cattaneo & Tobias 2014) © AAS. Reproduced with permission.

These results led Cattaneo & Tobias (2014) to propose a general ‘suppression principle’, namely ‘At high $Rm$ large-scale dynamo action can only be observed if there is a mechanism that suppresses the small-scale fluctuations.’ As discussed earlier, the kinematic framework dynamo theory is linear and the solutions are superposable – large-scale dynamo action can therefore only be observed if the small-scale dynamo is suppressed. In the cases discussed above the shear acted so as to suppress the fluctuations. It is possible, however, to envisage other mechanisms that may act in this way. In most circumstances, however, dynamo action proceeds in the nonlinear regime (see § 6). There one can envisage a nonlinear mechanism that leads to either the suppression of the small-scale dynamo or the saturation of the fluctuating field at a reasonable amplitude. In the saturated regime, whether the large-scale field can be observed depends on the ratio of the saturation amplitude for the large and small scales, which is itself a contentious issue in dynamo theory, as we shall see.

5.4.2. What does mean-field theory get right in the kinematic regime?

Mean-field theory describes the evolution of filtered or averaged quantities. In order for this theory to have utility in the kinematic regime, one would expect that application of the filter to the eigenfunctions of the full equations should correspond to the eigenfunctions of the filtered equations. Moreover, the growth rate of the large-scale structure of the solutions should coincide with the growth rate predicted by the filtered equations. As we have argued above, this second expectation is never met at high $Rm$ – the growth rate is determined by the stretching of the small scales (Tobias & Cattaneo 2008a).

However, it is possible to construct more complicated cases, for example with the addition of shear, where the filtered (mean-field) equations have a growth rate and a frequency (Parker 1955; Tobias & Cattaneo 2013b). The shear breaks the isotropy of the filtered equations and the isotropy of the statistics of the velocity in the unfiltered system and therefore one expects a breaking of symmetry of the statistics of the solutions. In that case one might hope that the frequency, as manifested by the speed of propagation of the dynamo waves to be the same both in the solutions of the filtered equations and those of the unfiltered equations in a statistical sense. Here, symmetry breaking acts to perturb the frequency away from zero and is controlled by a change in the symmetry of the large scales. Nigro et al. (2017) considered just such a case based on the Tobias–Cattaneo dynamo of § 5.4.1. They determined a range of solutions with a travelling helical (dynamo) wave, that remains coherent for long periods of time and whose frequency is determined by mean-field effects. The wave could only be identified from the rest of the structure by a persistent phase-coherent signal with the rest of the solution being incoherent in time. They argued from this that it is better to consider a definition of large-scale dynamo action that considers the time scale of evolution of the pattern, rather than one that relies on spatial scales alone, i.e. spatio-temporal averaging is the natural choice for determining organisation. Further, as kinematic large-scale dynamo action consists of a long-lived coherent pattern embedded in a sea of incoherent fluctuations; the filtered equations simply average over the fluctuations. Encouragingly, they predicted that mean-field electrodynamics is capable of predicting the frequency of this coherent signal.

This gives hope for the construction of a sensible framework for deriving a filtered equation in the nonlinear regime. As argued by Nigro et al. (2017), as nonlinearity becomes important different scales may continue to grow at different rates and saturate at relative amplitudes completely different from those of the kinematic phase, controlling the fluctuations relative to the mean. If this were the case then they argue that the filtering technique used to identify spatio-temporally coherent large-scale fields in the kinematic regime may continue to be utilised in the nonlinear regime. We shall therefore finally (and none too soon for the referee!) move into the dynamic regime, and examine solutions to the coupled Navier–Stokes/induction system.

6. Turbulent saturation of the dynamo instability

6.1. Basic considerations

Consider first the ‘simpler’ case where initially the magnetic field is weak, the Lorentz force (being quadratic in the magnetic field) can be neglected, and kinematic theory applies at least to the initial stages of the dynamo evolution. Given the arguments of the previous few sections, at large enough $Rm$ one expects the magnetic field to grow exponentially with a well-defined growth rate, with large and small scales growing at the same rate. The shape of the eigenfunction, and its coherence on time scales longer than typical time scales of the flow, may be controlled by the degree of breaking of reflectional symmetry of the flow. Eventually the kinematic phase must finish, as the Lorentz force becomes significant, and the dynamo instability saturates. Such a dynamo with a well-defined kinematic regime has been termed an essentially kinematic dynamo (Tobias et al. 2011a).

How the dynamo instability saturates depends on the form of the turbulence that underlies the initial instability. Presumably, the dynamo equilibrates on a range of scales with a mean amplitude of the magnetic field that depends on scale (so there will also exist a magnetic energy spectrum alongside the kinetic energy spectrum). However, different considerations are now required, just as for the kinematic regime, to describe the saturation of dynamos for the cases of low and high $Pm$ those cases with or without broken reflectional symmetry.

We shall first of all consider small-scale dynamo saturation, describing the current state of theory and numerical experiments. As might be expected, different considerations come into play for the saturation of organised magnetic fields and we introduce the contentious issue of the generation of large-scale field in the nonlinear regime at high $Rm$ in § 6.5. In order for the instability to saturate in a statistically steady state, on average the inductive terms in the induction equation must balance the diffusive terms. Given that, all other things being equal, the kinematic dynamo puts magnetic energy on the smallest scale possible, it is unlikely that an essentially kinematic dynamo will typically saturate by increasing the diffusion of the magnetic field – since the magnetic field in the kinematic regime is already structured to maximise diffusion; clearly some modification of the inductive nature of the velocity field is needed. Therefore the field usually saturates by primarily modifying the $\boldsymbol {\nabla } \times ({\boldsymbol u} \times {\boldsymbol B})$ term in the induction equation. The simplest way to achieve this is to modify the velocity field via the momentum equation.

The question then arises as to the nature of this modification. Is there an easily identifiable property of the velocity field that in the kinematic phase leads to the exponential growth of the field whilst in the saturated state yields a statistically stationary magnetic field? Discussions of dynamo saturation can be divided into three broad paradigms. The first is the equipartition argument, that the magnetic energy will grow until it is comparable with the kinetic energy. This is a nice criterion as it can be applied with great generality, without knowing much about the details of the dynamo system; however, it suffers from two main problems. It is easy to construct both examples where the saturation magnetic energy is substantially lower than the kinetic energy (Brummell et al. 2001) and where it greatly exceeds it (Stellmach & Hansen 2004). Even in systems where the energies are comparable, they may not be comparable at all scales (Vainshtein & Cattaneo 1992). The second paradigm is a marginal stability argument that the nonlinear effects of the magnetic field are to bring the system back to marginality. In some sense this must be the case as, on average, induction balances diffusion in the induction equation, as argued above. However, this is not an equivalent criterion to reducing $Rm$ to its marginal value by modifying a typical amplitude of the velocity. This is linked with the limitation of the usefulness of the definition of the global $Rm$ based on the system scale discussed in § 2.4.1.

The third paradigm is more sophisticated and invokes some subtle modification of the Lagrangian properties of the flow. Given the arguments of § 4.1 that chaos is required for a dynamo to operate at high $Rm$, it seems plausible that if the magnetic field is able to reduce the level of chaos in the flow then this can lead to saturation. Interestingly, it seems as though the general case is that saturation occurs without suppression of chaos at high $Rm$. Saturated dynamo flows have been shown to continue to stretch exponentially by Cattaneo & Tobias (2009) as evidenced by their ability to amplify exponentially a passive vector field. Indeed the set of vector fields that do not grow exponentially when stretched by the turbulence appears to have measure zero. Unfortunately this implies that there is no general theory for the saturation of the dynamos and that, perhaps unsurprisingly, the physics of the momentum equation will play an important role.

6.2. Saturation of high $Pm$ dynamos

6.2.1. High $Pm$ dynamos: general considerations and random flows

Just as for the kinematic regime, there is a large difference between high and low $Pm$ regimes with regard to the saturation of dynamos. In the high $Pm$ regime the dynamo is operating at scales in the sub-inertial range of the velocity for which the kinetic Reynolds number is small and the inertial term in the momentum equation can be neglected. The velocity can then formally be decomposed into two parts ${\boldsymbol u} = {\boldsymbol u}_{F} + {\boldsymbol u}_M$ where ${\boldsymbol u}_{F}$ is the original velocity driven by the forcing (either mechanical or buoyancy) i.e. it is the velocity of the kinematic dynamo problem and has a characteristic scale $\gg l_{\eta }$, the other ${\boldsymbol u}_M$ is the magnetically driven velocity. These two components then evolve independently (Brummell et al. 2001), as the inertial terms has been neglected.

To fix ideas, consider the case where the velocity is kinematically driven solely on one scale, say the integral scale of the calculation. As noted in § 4.1. this will lead kinematically to a growing dynamo eigenfunction for the magnetic field that is concentrated at the resistive $Rm^{-1/2}$ scale. If $Rm$ is large (i.e. in the kinematic fast dynamo regime) then this is much smaller than the scale of the velocity. These small-scale filamentary fields now drive flows on scales smaller than the velocity via the Lorentz force, eventually generating strong enough magnetic energy at the velocity scale to modify the original flow. Because of the simplicity of system at low $Re$, estimates for the strength of the generated magnetic field can be obtained where $B^2/U^2 \sim Rm^{1/2} Re^{-1}$, where $B$ is measured in units of the Alfvén velocity. This estimate is borne out at moderate $Rm$, but breaks down at higher $Rm$ (Brummell et al. 2001). At such low $Re$ and $Rm$ the saturation mechanism is simply to drive flows via the Lorentz force that can counteract the advective properties of the driven flow.

For multi-scale velocity fields at moderate and large $Re$ the situation is more complicated and saturation can be investigated using semi-analytical models, phenomenological models and numerical experiments. The semi-analytical models are ultimately based on some closure of the MHD equations. For random flows described by the Kazantsev formalism, the magnetically driven velocity produces a change in the velocity correlator, which leads to the nonlinear saturation (Subramanian 1999, 2003). Alternatively a Fokker–Planck equation for the probability distribution function for magnetic fluctuations can be derived. It can be shown that the coefficients of this equation again are determined by the velocity correlation function which can be modified nonlinearly in a similar manner to above (Boldyrev 2001). However, all of these results rely on some assumption of the role of the Lorentz force and that the flows remain $\delta$-correlated in time even in the presence of a saturating magnetic field.

At high $Pm$ an interesting phenomenological description for multi-scale dynamos has been proposed (Schekochihin et al. 2002). The initial kinematic growth of field at the resistive scale (which is much smaller than the viscous scale) is modified when ${\boldsymbol u} \boldsymbol {\cdot } \boldsymbol {\nabla } {\boldsymbol u} \approx {\boldsymbol B} \boldsymbol {\cdot } \boldsymbol {\nabla } {\boldsymbol B}$ at that scale. The left-hand side is easily estimated to be $v^2/l_\nu$. To estimate the right-hand side the authors use the foliated structure of the magnetic field to yield an estimate of $b^2/l_\nu$, where $b$ is a typical field strength at the scale $l_\eta$. The nonlinear saturation therefore begins when the magnetic energy comes into equipartition with the kinetic energy at the viscous scale, and dynamo growth associated with eddies at the viscous scale is halted. This suppression is believed to be the result of a subtle modification of the eddy geometry; the Lorentz force causes the eddies to align with the local magnetic field, thus reducing induction. However, slightly larger eddies can still sustain growth and so growth will continue until the magnetic energy comes into equipartition with the energy of these eddies; this process continues until the magnetic energy reaches equipartition with the kinetic energy at the integral scale. This nonlinear adjustment is characterised by a growth of the magnetic energy on an algebraic (rather than exponential) time with a shift to the generation of magnetic fields on larger scales. It is possible to construct models where the final state only reaches a fraction of global equipartition. This occurs when $Pm$ is large but $Pm \le Re^{1/2}$, and for this regime $B^2/U^2 \approx Pm/Re^{1/2}$ (Schekochihin et al. 2002). In this scenario, however, the characteristic scale of the saturated magnetic field is still smaller than the viscous scale.

6.2.2. High $Pm$ dynamos: numerical experiments

There have been many simulations of dynamo saturation at moderate to high $Pm$ that seem to confirm the phenomenological picture described above (Maron, Cowley & McWilliams 2004). In the kinematic phase the magnetic spectrum appears compatible with the $k^{3/2}$ prediction of the Kazantsev model, whilst the magnetic field does indeed appear foliated. Dynamo saturation occurs when the magnetic and kinetic energies are comparable; as the saturation progresses the magnetic spectrum grows and flattens with the creation of magnetic structures at larger scales. The magnetic energy always exceeds the kinetic energy and the field is more intermittent than the velocity, with the p.d.f. for the velocity field remaining close to that of a Gaussian whilst that for the magnetic field is better described by an exponential. The degree of intermittency is, however, reduced in the saturated state as compared with the kinematic state (Cattaneo 1999). Interestingly, the high $Pm$ formalism seems to persist for $Pm \sim { {O}}(1)$ (see Tobias et al. (2012) for more details) although quite how well the phenomenological picture is replicated is unclear, owing to numerical limitations. It is clear though that at $Pm \sim { {O}}(1)$ the magnetic fields are more intermittent that the velocity and the local magnetic field introduces statistical anisotropy in the flow in the nonlinear regime.

6.3. Saturation of low $Pm$ dynamos

At present, the low $Pm$ regime remains largely unexplored, owing to the computational constraints described earlier for the kinematic problem. Because the dynamo operates in the inertial range of the turbulence inertial terms certainly remain important. Analytical progress can only be made by resorting to standard closure models such as the eddy-damped quasi-normal Markovianised (EDQNM) model (Pouquet, Frisch & Leorat 1976) (see § 6.4 for a discussion of EDQNM and other closures).

However, some progress can be made by adopting dimensional arguments. For non-rotating systems, Fauve & Pétrélis (2007) identify a possible scaling for the level of saturation of dynamos at low $Pm$. For the realistic case where $Rm_c$ becomes independent of $Pm$ as $Pm$ gets small, the scaling is given by

(6.1)\begin{equation} \frac{B^2}{\mu} \sim {\rho U^2}, \end{equation}

so that the level of saturation of the magnetic energy becomes independent of $Pm$. This is in marked contrast to the scaling for large $Pm$ discussed above. Unfortunately, even with current computational power, simulations of dynamos at low $Pm$ are difficult, and it appears as though all simulations (whatever the quoted value of $Pm$) are currently in the $Pm \sim { {O}}(1)$ regime, with all the associated dynamics (see e.g. Mininni 2006; Iskakov et al. 2007). For reasons discussed in Tobias et al. (2012), any calculation that aims to test nonlinear saturation in the low $Pm$ regime (and that resolves the inertial range and the numerically growing eigenfunction) has a requirement in terms of grid points of several thousands (or tens of thousands). Using LES for the velocity should save a factor of ten in resolution (Ponty et al. 2007), but care must be taken in applying sub-grid modelling to MHD flows as emphasised by Miesch (2015). There is, however, a case where dynamo saturation at low $Pm$ can be calculated (Seshasayanan, Gallet & Alexakis 2017b). This utilises an asymptotic reduction of the equations for rapid rotation; as noted earlier rapid rotation has the tendency to make low $Pm$ dynamos act like high $Pm$ dynamos.

6.4. Analytic closure theories for nonlinear dynamos

Analytical progress for nonlinear dynamos is difficult and, just as for hydrodynamic turbulence, almost always relies on closure approximations. It is important to stress that closure theories for MHD are on shakier ground than those for hydrodynamic flows; whereas in hydrodynamic turbulence there is a vast body of experimental evidence and many systematic attempts to test closures and reduction hypotheses, there is no such thing in MHD. Much more needs to be done in evaluating the accuracy and applicability of closure hypotheses for MHD; intuition for such problems is often misleading. One other shortcoming of the procedure is that often expedient assumptions about the nature of the turbulence are made (for example homogeneity and isotropy) in order to make analytical progress in deriving the expressions for the EMF, and hence the transport coefficients. For geophysical and astrophysical situations these assumptions are often poor and so it would be better to proceed without making such assumptions – this is the philosophy behind direct statistical simulation (DSS) as described in § 10.

Here, we do not give details of the various possible closure schemes, as a discussion of these even for the hydrodynamic problem is well beyond the scope here. The interested reader is directed to the excellent review by Yokoi (2019). We only include here a flavour of the assumptions and shortcomings of the various closures.

The simplest closure scheme possible is to consider only quasilinear interactions in both the momentum equation and the induction equations. That is to discard eddy/eddy $\rightarrow$ eddy interactions in both equations. This is expected to work well for the case where there are systematic mean flows and fields, but may not be appropriate when cascade and inverse cascade processes play a dominant role.

For cases where the eddy/eddy $\rightarrow$ eddy interactions are key, more work is needed to capture the nature of the interactions. Following Yokoi (2019) we explain the schemes using a model nonlinear system that schematically represents the hydrodynamic momentum equation in the form

(6.2)\begin{equation} \frac{\partial u}{\partial t} = u u - p + \nu u; \end{equation}

extension of the formalism of the schemes to MHD, although not their applicability, is then straightforward. Here, the $u u$ term is simply a prototypical nonlinear term. Hence

(6.3)\begin{equation} \frac{\partial \langle u \rangle}{\partial t} = \langle u u \rangle - \langle p \rangle + \nu \langle u \rangle, \end{equation}

where $\langle \rangle$ indicates an average quantity. The nonlinearity in this equation of course leads to the problem that the second-order correlation is required to solve for the mean and that the third-order correlation is also needed because

(6.4)\begin{equation} \frac{\partial \langle u u \rangle}{\partial t} = \langle u u u \rangle + \nu \langle u u \rangle. \end{equation}

The fourth-order correlation appears in the equation for the evolution of the third-order correlation because

(6.5)\begin{equation} \frac{\partial \langle u u u \rangle}{\partial t} = \left(\int \langle u u u u\rangle \right) + \langle u u u u\rangle + \nu \langle u u u \rangle, \end{equation}

where the Poisson equation for the pressure has been used.

What is clear is that some closure is needed that can express the fourth-order correlations in terms of lower-order correlations (either via a functional form or by deriving an evolution equation). For inhomogeneous systems this closure can be achieved by setting the values of higher order cumulants (which for the first three correspond to centred moments) to zero (see § 10).

However, if one also makes the assumption of homogeneity, then progress can be made by transforming to spectral (wavenumber) space and deriving the evolution equation for the spectral energy, $E(k)$, in terms of the spectral energy flux, $\varPi (k)$, and the spectral energy transfer, $T(k)$. Because the energy is a measure of the second-order correlations and the flux and energy transfer represent triple correlations, a closure can be achieved by expressing $T(k)$ or $\varPi (k)$ in terms of $E(k)$.

The most basic attempt to achieve this is to consider quasi-normal models (Monin, Yaglom & Lumley 1975) in which the p.d.f. of the velocity field is Gaussian and so the fourth-order correlation may be expressed as a sum of the products of the second-order correlations. However, this severe degree of approximation can lead to lack of realisability in the sense that negative energies may be achieved.

In order to alleviate this problem, which can be traced to the overestimation of triple correlations, Orszag (1970) suggested placing a scale-dependent damping term $\tau (k)$ in the evolution equation for the triple correlations. Together with the assumption that the turbulence does not depend on its history (i.e. lacks memory via Markovianisation) this leads to the EDQNM approximation. The addition of the scale-dependent eddy-damping time can certainly resolve the issue of non-realisability; however, the eddy-damping time depends critically on an assumed eddy-damping rate $\mu ({\boldsymbol k)}$, which is a given parameter. Indeed, simplifications (for example the so-called ‘minimal $\tau$ approximation’, where the eddy-damping time is itself a parameter independent of scale) have also been suggested and explored (Kleeorin et al. 2002).

One of the central uncertainties of utilising closure models for MHD turbulence and dynamos arises because the eddy-damping rate (and hence time) does depend on the magnetic field strength and is not given a priori. Magnetic fields lend memory to the turbulence through the presence of restoring forces and waves, and this should be taken into account in any theory. For all the theories described above, some assumption must be made about this dependence, which does not arise self-consistently from the theory. It is certainly the case that, whereas in hydrodynamic turbulence there is overwhelming evidence that this quantity is of the order of the turnover time, there is no consensus on what this quantity should be in MHD theory. Indeed it may depend sensitively on the field strength and magnetic Reynolds number of the magnetised turbulence.

However, one scheme that does enable some self-consistent derivation of the time scale for the turbulence is the direct-interaction approximation (DIA) introduced by Kraichnan (1959). Here, via Green's function techniques, a closed system of equations comprising an equation for the correlation function and one for the response function is derived. The method of derivation is long and involved (and beyond the scope of this paper – although see Yokoi 2019), taking in such mathematical joys as perturbation expansions, partial summations and renormalisation methods via Feynman diagrams. It goes without saying that, in order to close the system, some assumptions must be made even for this technically challenging procedure. There is much work to be done on understanding the basic response of turbulence in the presence of a magnetic field at different $Rm$ even for the case where the magnetic field has not been self-consistently generated by the turbulence.

Finally for this section, we stress that the closure models described above are all essentially formulated and solved in spectral space, and that, as noted, this approach relies on the turbulence being homogeneous (and for some simplifications isotropic). The astrophysical turbulence that leads to the generation of cosmical magnetic fields is almost never homogeneous (and certainly not isotropic) and so this approach may be limited in applicability. One approach aimed at circumventing this issue is to derive statistical models that do not rely on homogeneity. Perforce, these models are less amenable to analytical descriptions and solution and numerical algorithms and techniques must be developed for their progression (see § 10 for a discussion of DSS). However, some progress can be made in including inhomogeneity via the two-scale direct interaction approximation. This approach combines multiple-scale analysis (of the type used to derive mean-field theory) with the DIA. In short, the method proceeds by introducing two time scales (one order one – called fast for convenience – and one slow) and spatial scales (one small and one order one), so that

(6.6a–d)\begin{equation} {\boldsymbol \xi} = {\boldsymbol x}, \quad {\boldsymbol X} = \epsilon {\boldsymbol x}; \quad \tau = t, \quad T = \epsilon t, \end{equation}

with $\epsilon \ll 1$. All variable are then expanded into mean fields and perturbations so that

(6.7)\begin{equation} f = \overline{f}({\boldsymbol X},T) + f'({\boldsymbol \xi},{\boldsymbol X}; \tau,T). \end{equation}

Progress is now made by assuming that the fluctuation field is homogeneous with respect to the fast space variable, and expanding the fluctuating variables in a power series in the small parameter $\epsilon$, so that $f' = f'_0+\epsilon f'_1+ \cdots$, making assumptions about the statistical properties of $\,f_0$ and calculating the turbulent fluxes (for example Reynolds stresses and EMFs) via DIA. This has the potential to be a formidable approach leading to breakthroughs for certain instability and dynamo problems (see e.g. Yokoi 2018).

6.5. Saturation of systematic fields

In this section we consider probably the most contentious issue in turbulent dynamo theory, the mechanism for the nonlinear saturation of systematic fields. It is fair to say that this issue is technically and computationally challenging and no consensus has emerged. It would be possible to construct a whole review of this topic without really doing it justice, so here I will only recount the salient arguments.

We have demonstrated how prescribed turbulent flows with broken reflectional symmetry can lead to the generation of an EMF that leads to the direct generation of systematic magnetic fields (in a manner analogous to the driving of mean flows via a Reynolds stress in the momentum equation). There are two possible general mechanisms for how systematic fields may continue to be generated in the nonlinear regime. The first is that the kinematic mechanism of direct driving via an EMF continues to proceed in the nonlinear regime. Here we stress that direct driving refers to eddy/eddy $\rightarrow$ mean interactions (in the language of fluid turbulence). If this were the case then it might be possible to develop a theory for the nonlinear behaviour of this driving, including parameterising the role of the large-scale magnetic field in modifying the transport coefficients. A second possibility is that in the nonlinear regime a scale-by-scale inverse cascade proceeds whereby energy is successively transferred locally in scale to larger spatial scales, eventually reaching the system scale of the object. It is fair to say that currently both of these mechanisms (direct driving and inverse cascades) are believed to contribute.

As noted above, even in the nonlinear regime, the induction equation is still formally linear in the magnetic field and so it can be argued that the machinery of kinematic mean-field electrodynamics can, and should, still be of utility. Indeed, the arguments there were based on solution of the induction equation alone, which is an important component of the nonlinear dynamo system.

It should even be possible to continue to couch the direct generation of systematic magnetic fields in the nonlinear regime in terms of the calculation of transport coefficients similar to those outlined in § 5.3.3. What would be needed here is to utilise the saturated velocity field in the equations for the evolution of the fluctuating magnetic field ((5.31)) to calculate the evolution of the EMF. Formally this, of course, would yield the correct prescription, as this is simply a proxy for solving the induction equation. There are two major problems with this approach, however. The first is that there is no general theory for the turbulent response of the fluctuating velocity field to the presence of a magnetic field (either large or small scale) as should have become clear in our discussions earlier. The second is that, even were such a theory possible, we have seen that determination of the transport coefficients can only proceed if certain approximations (for example low $Rm$ or short correlation time, $\tau _c$, for the turbulence) are utilised. So for example, if the turbulence is assumed to maintain a short correlation time then the transport coefficients can be described by equations such as (5.38), where now ${\boldsymbol u}$ and ${\boldsymbol \omega }$ must be interpreted as the saturated velocity and vorticity. Astrophysical dynamos are certainly not at low $Rm$. Moreover it has often been argued that a significant consequence of the presence of magnetic field in turbulence is to add memory to the flow, thereby increasing the correlation time ($\tau _c$) of the turbulence and making the approximations of first-order smoothing less valid in the nonlinear regime (even if they were assumed valid in the kinematic regime). It seems as though the desire to express everything in terms of a mean-field theory based on kinematic considerations is causing considerable complications.

6.5.1. Catastrophic (Vainshtein–Cattaneo) quenching

Perhaps the greatest challenge to mean-field theory in the dynamic regime arises from arguments proposed in the early nineties by Vainshtein and Cattaneo (Cattaneo & Vainshtein 1991; Vainshtein & Cattaneo 1992) and extended by Gruzinov & Diamond (1994). The arguments concern the level of saturation of the organised field that can be generated by an EMF, before the EMF is itself suppressed by the nonlinear effects of Lorentz force in the momentum equation.

One way to proceed would be to assume on energetic grounds that the Lorentz force can only begin to suppress a kinematic effect once the energy in the magnetic field is comparable with that of the turbulence. Traditionally (and naively) it has been assumed that this occurs when the magnetic energy of the mean field is in equipartition with the kinetic energy of the flow, i.e. when

(6.8)\begin{equation} \langle \rho u^2 \rangle \sim B_0^2 / \mu, \end{equation}

where $B_0^2 = \overline{{\boldsymbol B}} \boldsymbol {\cdot } \overline{{\boldsymbol B}}$ is the energy in the mean field. Hence a simple phenomenological model would then be to take the nonlinear dependence of the transport coefficients as

(6.9a,b)\begin{equation} \alpha = \frac{\alpha_0}{1+ B_0^2/ {\mathcal{B}}^2 }, \quad \beta = \frac{\beta_0}{1+ B_0^2/ {\mathcal{B}}^2 },\end{equation}

where $\alpha _0$ and $\beta _0$ are the values of the transport coefficients derived in the kinematic approximation and ${\mathcal {B}}^2$ is the non-dimensional equipartition energy (i.e. normalised with the kinetic energy of the flow). Note that this prescription implicitly assumes that the mean field dominates (or at least is of the same order of magnitude as) the fluctuating field in the saturated regime, and so it is the mean field that leads to saturation of the transport. The expressions in (6.9a,b) are attractive as they lead to a saturation of the mean field for an energy comparable with that of the turbulence. Moreover these expressions rely only on the calculation of the mean or organised fields, which requires little in the way of computational expense.

However, the assumption that the magnetic energy is dominated by the organised component is questionable at high $Rm$ as pointed out clearly by Cattaneo & Vainshtein (1991) and Vainshtein & Cattaneo (1992). The distribution of magnetic energy with spatial scale at high $Rm$ is a key issue for understanding nonlinear dynamo saturation. One can proceed on either physical or mathematical grounds to understand this distribution – computational models, although suggestive and informative, are not at the point of unambiguously settling this issue.

On physical grounds, it is argued that if there is a large-scale component of the field generated kinematically in a turbulent dynamo then turbulence will also amplify magnetic fields on the small scale via the action of eddies with short turnover times. Hence one might expect the small-scale fields (kinematically) to be more energetic than those on large scales. If diffusion plays a role in setting the ratio of these energies then one might expect a relationship of the form

(6.10)\begin{equation} \overline{b^2} \sim Rm^p B_0^2 ,\end{equation}

where $b^2 = {\boldsymbol b'}\boldsymbol {\cdot } {\boldsymbol b'}$ and $0 \le p \le 2$ is a flow- and geometry-dependent coefficient.

If this is the case, and assuming that magnetic fields with this energy do distort the eddies with correlations that lead to the generation of organised field (i.e. those with significant broken reflectional symmetry), then this may indicate that the transport coefficients are altered significantly once this strong small-scale field reaches equipartition with the kinetic energy of the turbulence, i.e. when

(6.11)\begin{equation} \langle \rho u^2 \rangle \sim Rm^p B_0^2 / \mu.\end{equation}

One might then expect that the expressions in (6.9a,b) to be replaced by

(6.12a,b)\begin{equation} \alpha = \frac{\alpha_0}{1+Rm^p B_0^2/ {\mathcal{B}}^2 } , \quad \beta = \frac{\beta_0}{1+Rm^p B_0^2/ {\mathcal{B}}^2 } . \end{equation}

At high $Rm$, one can then argue that the efficiency of the transport coefficients is strongly suppressed, even when the organised field $B_0$ is small; a phenomenon known as ‘catastrophic’ quenching (Cattaneo & Vainshtein 1991; Vainshtein & Cattaneo 1992). It certainly seems plausible that the transport coefficients begin to be suppressed when the mean field is weak. However, the generation of a mean EMF relies on correlations in the turbulence between the fluctuating magnetic field and eddies, which is a subtle interplay, so arguments on physical grounds may be misleading. Hence, other approaches to this issue have been employed.

The first of these is to derive analytical expressions for the EMF (and hence the transport coefficients) utilising closure approximations (such as those discussed in § 6.4) and to then apply conservation laws that are valid in the non-dissipative system. This approach can yield significant insight (in fact even understanding the role played by conservation laws is important) but is only as good as the closure approximations that are utilised. For a more extensive discussion see Diamond, Hughes & Kim (2005) and Brandenburg & Subramanian (2005). We give a brief explanation here.

Let us start by describing a closed domain, with no fluxes in or out of the region. We noted in (2.19) that the total magnetic helicity changes only via dissipation or through surface fluxes from the domain. Hence, for an ideal closed system, magnetic helicity is a quadratic (although not sign-definite) invariant quantity. For a closed dissipative system with no surface fluxes, global magnetic helicity can only be created or destroyed via irreversibility introduced by diffusive processes, i.e.

(6.13)\begin{equation} \frac{\textrm{d}}{\textrm{d} t} \langle {\boldsymbol A}\boldsymbol{\cdot} {\boldsymbol B}\rangle = - 2 \eta \langle\kern1.5pt {\boldsymbol j}\boldsymbol{\cdot} {\boldsymbol B} \rangle. \end{equation}

We recall here that angular brackets refer to volume averages. Note at this point that, for this closed system, this immediately allows some inference about $\langle\kern1.5pt {\boldsymbol j}\boldsymbol {\cdot } {\boldsymbol B}\rangle$ in the steady state, namely that $\langle\kern1.5pt {\boldsymbol j}\boldsymbol {\cdot } {\boldsymbol B}\rangle = 0$.

Furthermore, if one defines an intermediate average (such as coordinate averaging in the mean-field sense as discussed in § 5.1) and denotes it by an overbar then in the steady state

(6.14)\begin{equation} \langle\kern1.5pt {\boldsymbol j}' \boldsymbol{\cdot} {\boldsymbol b}' \rangle + \langle\kern1.5pt \overline{{\boldsymbol j}} \boldsymbol{\cdot} \overline{{\boldsymbol B}} \rangle = 0, \end{equation}

so that the average current helicity from the fluctuations is of the opposite sign to that from the means in the steady state, and if there is no mean current then $\langle\kern1.5pt {\boldsymbol j}' \boldsymbol {\cdot } {\boldsymbol b}' \rangle = 0$.

Indeed, if there is no mean current another exact result emerges from consideration of Ohm's law, giving

(6.15)\begin{equation} {\boldsymbol{\mathcal{E}}} \boldsymbol{\cdot} {\boldsymbol B}_0 = -\frac{1}{\sigma} \langle\kern1.5pt {\boldsymbol j'} \boldsymbol{\cdot} {\boldsymbol b'} \rangle + \langle {\boldsymbol e'} \boldsymbol{\cdot} {\boldsymbol b'} \rangle = \alpha B_0^2,\end{equation}

where ${\boldsymbol B}_0$ is the uniform mean field, with no associated current. This exact result relates the transport coefficient $\alpha$ to the fluctuating current density and the average correlation between the fluctuating magnetic field and EMF.

Finally, for a closed domain one may obtain evolution equations for the magnetic helicity contributions from large and small scales, namely

(6.16)\begin{equation} \frac{\textrm{d}}{\textrm{d} t} \langle \overline{{\boldsymbol A}}\boldsymbol{\cdot} \overline{{\boldsymbol B}}\rangle = 2 \langle \overline{{\boldsymbol {\mathcal{E}}}}\boldsymbol{\cdot} \overline{{\boldsymbol B}} \rangle - 2 \eta \langle\kern1.5pt \overline{{\boldsymbol j}}\boldsymbol{\cdot} \overline{{\boldsymbol B}} \rangle, \end{equation}

and

(6.17)\begin{equation} \frac{\textrm{d}}{\textrm{d} t} \langle {{\boldsymbol a}'}\boldsymbol{\cdot} {{\boldsymbol b}'}\rangle = -2 \langle \overline{{\boldsymbol {\mathcal{E}}}}\boldsymbol{\cdot} \overline{{\boldsymbol B}} \rangle - 2 \eta \langle\kern1.5pt {{\boldsymbol j}'}\boldsymbol{\cdot} {{\boldsymbol b}'} \rangle. \end{equation}

Hence, the $\langle \overline {{\boldsymbol {\mathcal {E}}}}\boldsymbol {\cdot } \overline{{\boldsymbol B}} \rangle$ term acts to transfer magnetic helicity between fluctuating and mean magnetic fields. It is clear then for this closed system that if the large-scale magnetic helicity is to grow it must be at the expense of the small-scale magnetic helicity; recall though that magnetic helicity is not sign definite so it is certainly possible for these to grow together if they are of opposite signs.

We note here that the expressions derived above allow evaluation of an expression for the projection of the average EMF onto the mean field, i.e. ${\boldsymbol {\mathcal {E}} \boldsymbol {\cdot } {\boldsymbol B}_0}$, in terms of correlations between ${\boldsymbol j}'$ and ${\boldsymbol b}'$, instead of correlations between ${\boldsymbol u}'$ and ${\boldsymbol b}'$. However, in the nonlinear regime we know none of ${\boldsymbol u}'$, ${\boldsymbol b}'$ or ${\boldsymbol j}'$, so what has this manipulation gained us?

The answer arises from utilising results from calculations employing the closure approximations described in § 6.4 – specifically calculations employing the EDQNM approximation, and generalisations thereof. These calculations (see e.g. Pouquet et al. 1976; Kleeorin et al. 2002) yield an approximate relationship between the $\alpha$-coefficient for homogeneous isotropic turbulence and the kinetic and current helicity, namely that

(6.18)\begin{equation} \alpha \sim - \frac{\tau_c}{3} \langle {\boldsymbol u}' \boldsymbol{\cdot} {\boldsymbol \omega}' - {\boldsymbol j}' \boldsymbol{\cdot} {\boldsymbol b}' \rangle. \end{equation}

At this point it is worth commenting again on the status and interpretation of this expression. Unlike the earlier global balance relationships, which are exact, it is to be stressed that this expression can only be derived via approximate closure relationships or via linearisation about a pre-existing turbulent MHD state at low $Rm$ and short $\tau _c$ (Proctor 2003). Equation (6.18) has clearly replaced the calculation of correlations between ${\boldsymbol u}'$ and ${\boldsymbol b}'$ with those between ${\boldsymbol j}'$ and ${\boldsymbol b}'$. This is only possible if some dynamical relationship between ${\boldsymbol j}'$ and ${\boldsymbol u}'$ has been approximated.

One key limitation in this expression is the assumption that the correlation time of the turbulence, $\tau _c$ remains unaffected by the presence of large-scale magnetic field, which is not necessarily true at high $Rm$, as discussed earlier. Nonetheless, what is clear is that the expression (6.18) does reduce to that given in (5.38) in the kinematic regime; the correction to the kinematic result often being termed a magnetically driven $\alpha$-effect, $\alpha _M$.

Because, in this approximation, the current helicity from the small scales plays a key role determining the correction to the transport coefficient, it is now possible to combine this with the exact results of (6.17) to derive an evolution equation for the transport coefficient $\alpha$. For details of this see Gruzinov & Diamond (1994), Kleeorin, Rogachevskii & Ruzmaikin (1995) and Blackman & Brandenburg (2002); which include the further assumptions/simplifications (such as relating $\langle {\boldsymbol a}' \boldsymbol {\cdot } {\boldsymbol b}'\rangle$ to $\langle\kern1.5pt {\boldsymbol j}' \boldsymbol {\cdot } {\boldsymbol b}'\rangle$ via a spectral relationship) required to derive the evolution equation.

The evolution equation is given in full in Blackman & Brandenburg (2002) who also show that, in the saturated steady state for large-scale fields with no associated currents (i.e. uniform fields),

(6.19)\begin{equation} \alpha = \frac{\alpha_0}{1+Rm B_0^2/ {\mathcal{B}}^2 }, \end{equation}

a result first derived by Gruzinov & Diamond (1994), which is in agreement with the catastrophic quenching result of Cattaneo & Vainshtein (1991) in (6.12a,b). There is therefore now little debate as to the nature of the quenching of the $\alpha$-effect in closed systems with no large-scale currents; particularly as we shall see that the quenching formula is consistent with the results of numerical experiments (albeit at low and moderate $Rm$.)

Physically this result can be interpreted as the conservation of magnetic helicity placing a topological constraint on the nature of the magnetic field (as discussed in § 2.3.1). In order to grow the large-scale magnetic field, the knotted small-scale field that is naturally generated kinematically must be untangled, which is impossible if the topological constraint is maintained. If the presence of small magnetic diffusion is the only source of irreversibility allowing the untangling of magnetic field, then this places a severe constraint on the level of the mean field that can be generated – or at least the time scale on which it can be generated.

One interpretation of a formula such as that of (6.19) is that large-scale field generation is completely suppressed for mean fields an order of magnitude smaller than equipartition ($B_0 \sim {O}(Rm^{-1/2} \mathcal {B})$). However, consideration of the full evolution equation given in, say, Blackman & Brandenburg (2002) demonstrates that even with this ‘catastrophic’ quenching formula the energy in the large-scale magnetic field continues to grow on a time scale controlled by the irreversible process, in this case magnetic diffusion (an ohmic time). In the long run (for example in the case of the Sun a time scale comparable with the age of the star) reaching a statistically steady state with significant large-scale field may be possible, see also Moffatt & Dormy (2019). However, ‘in the long run we are all dead’ (Keynes 1923).

6.5.2. Large-scale currents and helicity fluxes: can they help?

‘E pur si muove’ Attributed (perhaps incorrectly) to Galileo.

Catastrophic quenching as envisaged by Vainshtein & Cattaneo (1992) and Gruzinov & Diamond (1994) places a severe constraint on the applicability of mean-field theory at high $Rm$. Large-scale organised fields can only be generated from a weak seed field on an ohmic time scale, this being the time taken for the action of ohmic dissipation to untangle the knotted field. However, it may be that other processes can be identified that break reversibility and hence the topological constraints placed on the field, which we shall discuss below. Another option is that the magnetic field may not have been generated all the way from kinematic values by a turbulent EMF, but may simply be sustained by an EMF that arises from an instability of the field itself. These more laminar ‘essentially nonlinear’ dynamos will be described in § 7.

When a mean current $\overline{\boldsymbol{j}}$ is allowed, this naturally leads to the presence of large-scale current helicity $\langle\kern1.5pt \overline {{\boldsymbol j}} \boldsymbol {\cdot } \overline{{\boldsymbol B}}\rangle$. Care does need to be taken here as the presence of a large-scale current leads to both questions of the nature of the separation of scales and the presence of a turbulent diffusion because of the presence of large-scale field gradients. The inclusion of a mean current may seems natural, since it is difficult to see how a dynamo can generate a magnetic field with no mean current. In that case (6.19) (subject to the same assumptions and approximations) is modified to

(6.20)\begin{equation} \alpha = \frac{\alpha_0}{1+Rm B_0^2/ {\mathcal{B}}^2 } + \frac{Rm \, \beta \langle\kern1.5pt \overline{{\boldsymbol j}} \boldsymbol{\cdot} \overline{{\boldsymbol B}}\rangle}{{{\mathcal{B}}^2 }+Rm \langle \overline{{\boldsymbol B}}^2 \rangle}, \end{equation}

as first derived by Gruzinov & Diamond (1994); here, $\beta$ is the magnitude of the turbulent diffusivity. Thus, although the kinematic $\alpha$-effect is catastrophically quenched, if one can (somehow) generate large-scale fields with a non-trivial current helicity then they can be maintained via a magnetically driven EMF, the second term in (6.20) (see also § 7).

However, it is not at all clear how one generates these significant large-scale fields from a weak seed field (e.g. in a galaxy). It may be that shear flows play a major role here in producing a large-scale field strong enough so that nonlinear effects can kick in. However, if one is relying on turbulence (i.e. fluctuation–fluctuation interactions) to do the job then there must be another process breaking reversibility for this to proceed quickly. Blackman & Field (2000) proposed that catastrophic quenching, or the associated generation of large-scale fields on an ohmic time scale, may be circumvented if irreversibility is introduced to the system by allowing magnetic helicity fluxes into and out of the domain. For example if magnetic helicity associated with small scales is preferentially lost on a time scale $\tau _{flux}$ then the rate of generation of large-scale field may occur on this time scale. Mathematically this can be seen by reinstating the helicity flux term into (6.13) and hence also into (6.16) and (6.17). The presence of this term breaks the reliance on ohmic processes to reconfigure (untangle) the magnetic field, and hence may alleviate the catastrophic quenching or allow the more rapid generation of large-scale field. There are two, physically distinct, contributions to the helicity flux; one diffusive and one ideal. Clearly the diffusive flux should become analogously weak to the ohmic dissipation of magnetic helicity in the limit of high $Rm$, so we are relying on the ideal flux to do the job!

There are a number of issues with this picture that require further investigation. The first is that the ideal flux term takes the form of a surface integral; it is not obvious whether a term that relies on losses through a surface can lead to a significant impact on the generation of large-scale fields throughout a large volume. The second is whether the ideal fluxes manage to avoid scaling with the magnetic Reynolds number (not directly, since the diffusivity does not appear explicitly in these terms, but through the form of the magnetic field). This may be checked using numerical calculations as described in the next subsection, at least at moderate $Rm$; this potential lifebelt for the generation of large-scale fields at high $Rm$ should be investigated more closely in the near future. Finally, it may be that it is enough to separate spatially the region of large-scale field generation from that of large-scale field storage. Faced with the proposal of catastrophic quenching, Parker (1993) introduced the concept of interface dynamos, which made use of precisely this idea. More investigation of dynamos with inhomogeneous turbulence and shear flows may reveal whether this model can function at high $Rm$. Shear-enhanced diffusion and other methods of introducing irreversibility on a fast time scale may also play an important role here.

6.5.3. Numerical investigations of catastrophic quenching

The theoretical arguments proposed above can be somewhat tested by devising numerical experiments. These, in general, take two complementary forms: dynamo experiments examine the nature of the generated magnetic field when it is allowed to evolve freely in a turbulent flow, whereas turbulent transport experiments typically impose a large-scale field (perhaps a uniform field or a large-scale field with a large-scale current) and measure the response of the turbulence in creating an EMF. Of course, the first of these approaches is more natural, with the magnetic field self-consistently evolving with the turbulence. The second approach, though, offers more control and the potential to determine the functional dependence of the EMF on parameters such as the imposed field strength, rotation rate, and magnetic Prandtl number.

However, it must be stressed that such numerical simulations are limited in the range of $Rm$ that are accessible and it is difficult to extrapolate the results of these to the astrophysically relevant high $Rm$ regime. For example, a spectral calculation performed with $10^9$ degrees of freedom with a reasonable separation of scales between large and small scales is able to reach magnetic Reynolds numbers of $ {O}(10^3)$ (Hughes & Tobias 2010).

Simulations are usually configured in one of two ways. In the first a flow is driven via a prescribed body force, typically chosen to drive a relatively small-scale flow with significant kinetic helicity (to maximise the chances of driving an EMF and hence a large-scale dynamo) and potentially a shear flow. In the second driving is via buoyancy in thermal convection, where helicity emerges naturally via the interaction with rotation; clearly the second of these is of more relevance astrophysically although there is less control over the form of the flow.

The list of numerical experiments performed over the past twenty or so years is long, and a complete review is well beyond the scope here. The interested reader is directed to Brandenburg & Subramanian (2005) for a thorough review of the early attempts and Brandenburg (2018) for a briefer review of later attempts.

In a closed domain all simulations (albeit at low and moderate $Rm$) are consistent with the picture of catastrophic quenching of the EMF (or diffusive growth of the mean field). For example Cattaneo & Hughes (1996) utilised a forcing that, in the absence of magnetic field, was designed to drive a small-scale Galloway–Proctor flow (here small-scale is defined as $k=1$) and measured the response of the turbulent EMF to an imposed uniform ($k=0$) magnetic field in a triply periodic box. They found results consistent with (6.12a,b). Analogous results were found for dynamo calculations driven by helical turbulence by Blackman & Brandenburg (2002) (and see references therein), where for closed systems the large-scale field only grew ohmically after an exponential kinematic phase. It seems as though there is now reasonable agreement for homogeneous, confined systems; the magnetic field can only emerge after an ohmic time scale.

For open systems, where boundary conditions that allow magnetic helicity to escape or enter the domain, the situation is more controversial, and more detailed calculations at higher $Rm$ are needed. Briefly, it does appear that the dynamos in these open domains do seem to have a different character to those in closed domains – at least at moderate $Rm$. For example, Hubbard & Brandenburg (2010) found that, as $Rm$ is increased both the diffusive volume term and helicity flux contribution to irreversibility decrease, although intriguingly the diffusive volume term decreases more rapidly as a function of $Rm$ (at least for moderate $Rm$). Linear extrapolation would seem to indicate that the boundary term would dominate at sufficiently high $Rm$. Similar results have also been found by Del Sordo, Guerrero & Brandenburg (2013). However, it does seem in all current cases as though the ratio of mean field to fluctuating field continues to decrease with $Rm$ (Brandenburg 2018; Cattaneo, Bodo & Tobias 2020). Clearly this issue is difficult to address computationally, although the status of mean-field theory in the nonlinear regime may become clearer in the next few years.

7. Essentially nonlinear dynamos

Up to this point, we have focussed on dynamos (whether large or small scale) that have a well-defined kinematic regime; the role of the Lorentz force is simply to saturate the dynamo. These types of dynamo have been termed essentially kinematic (whether in the linear or nonlinear regime), and much attention has focussed on whether such dynamos may act to generate organised magnetic fields, as we have seen.

However, another type of dynamo is possible; one where the components of the flow that lead to the generation of the field are themselves driven by the instabilities of, or the intervention of, a finite amplitude magnetic field. Such a dynamo would not exist in the limit of vanishingly small magnetic field (it is necessary to have a finite amplitude magnetic perturbation to sustain such a dynamo). These dynamos, which have been termed essentially nonlinear (Tobias et al. 2011a), often use the presence of a finite amplitude magnetic field to facilitate the extraction of energy from an energy reservoir (typically a shear flow that cannot act as a dynamo in its own right).

It has been argued that such essentially nonlinear dynamos are better candidates for generating organised magnetic fields at high $Rm$ than essentially kinematic dynamos. The line of reasoning here is that, as argued at length above, organised magnetic fields are generated through correlations between the small-scale, usually turbulent, velocity and small-scale magnetic fields. For essentially kinematic dynamos, the small-scale flow exists in the absence of magnetic field and is driven by either an external forcing or a hydrodynamic instability. The magnetic field is then generated by the induction of turbulent motions; at high $Rm$ it is difficult to maintain the correlations between the magnetic field and driving flow that are required for a net EMF. The loss of correlations therefore leaks a weakening of the mean-field generation mechanism.

However, in the contrasting case of an essentially nonlinear dynamo, the finite amplitude magnetic field drives dynamics (either as a result of a hydromagnetic instability or through relaxing a constraint on the flow) yielding magnetic field perturbations and velocity field perturbations that are likely to be correlated (as they have the same source). Of course, this process may be affected adversely by the presence of turbulence that will act so as to destroy these correlations. This is, I believe, an interesting avenue for future research.

Some examples of essentially nonlinear magnetically driven dynamos have been examined in detail. One of particular interest is that dynamo that exists through the interaction of a shear flow with the magnetic buoyancy instability (Cline, Brummell & Cattaneo 2003). In this scenario large-scale toroidal magnetic field is generated from large-scale poloidal field by a shear flow. The regeneration of poloidal field occurs owing to the combined action of magnetic buoyancy and Kelvin–Helmholtz instabilities; this can only occur if the initial magnetic fields exceed a critical threshold. Furthermore, the instabilities lead to fluctuating velocities that are strongly correlated with the unstable magnetic fields – unsurprising as these velocities are themselves driven by the field. This leads to the generation of a large-scale magnetic field and the cycle can begin again. To a fluid dynamicist the description should be reminiscent of a self-sustaining process, particularly that described by Waleffe (1997) for nonlinear transition in wall-bounded shear flows.

This analogy has been carried much further in the description of another essentially nonlinear dynamo in an excellent series of papers by Rincon and collaborators (see e.g. Rincon et al. (2007); Riols et al. (2013) and others). Here the dynamo is driven by the magnetorotational instability, which is a joint instability of differential rotation and magnetic field. Here, a Keplerian shear flow is hydrodynamically stable via Rayleigh's criterion. However, in a finite domain, and in the presence of dissipation, a finite amplitude magnetic field can lead to the destabilisation of the shear and the generation of turbulence. This turbulence may itself lead to the sustainment of the magnetic field sufficient to maintain the instability. Such a ‘bootstrapping’ dynamo bears all the hallmarks of a self-sustaining process, and the strong analogies with nonlinear transition in shear flows have been emphasised, as shown in figure 9, which is adapted from Riols et al. (2013).

Figure 9. (a) Hydrodynamic self-sustaining process (Waleffe 1997) and (b) self-sustaining dynamo processes in shear flows prone to MHD instabilities (Rincon, Ogilvie & Cossu 2007). After Riols et al. (2013).

It remains to be seen whether these mechanisms continue to work at high $Rm$, but this is an extremely promising avenue for research in astrophysical dynamos. Other such dynamos that may benefit from such an analysis include dynamos driven by magnetic buoyancy, current-driven instabilities, joint instabilities of differential rotation and magnetic fields and possibly the geodynamo. Methods from transition, such as the calculation of minimal seeds and exact nonlinear solutions, could be useful in aiding progress in all of these cases.

8. Balances in rapidly rotating magnetised convection and the geodynamo

This perspective has focussed mainly on the key problem of generating magnetic field in the astrophysically relevant limit of high magnetic Reynolds number ($Rm$). However, for the construction of models of planetary magnetic fields, the conceptual difficulty lies not in the solution of the induction equation at high $Rm$, but rather that of the momentum equation in the rapidly rotating limit in the presence of a magnetic field. This difficulty takes the form of the solution of equations with variability on a vast range of temporal scales. Here we shall focus on descriptions of the geodynamo, although of course the theoretical difficulties will propagate through to models of other rapidly rotating solar system planets and exoplanets.

8.1. Time scales in the geodynamo

Earth's magnetic field is believed to be generated in the outer core by the motion of the iron-rich fluid there. The driving mechanism for this motion is still the subject of debate; most models utilise driving by thermal or compositional convection, although the role of instabilities driven by precession or tides are still of interest. The Earth, of course, rotates with a period $(T_\varOmega )$ of one day, which sets the shortest time scale of interest for the geodynamo. Once the magnetic field is generated by the dynamo it allows magnetically mediated waves (so called MAC and torsional waves) with periods of years and decades. It is believed that a typical convective time $T_c$ for the Earth is of the order of $10{-}100$ years, whilst the magnetic field diffuses on a time scale $T_\eta \sim 10^4$ years. Finally the viscous and thermal time scales are extremely long $T_\kappa \sim T_\nu \sim 10^9{-}10^{10}$ years. The ratio of these time scales is encoded in the relevant non-dimensional parameters, i.e.

(8.1a–c)\begin{gather} Ro = \frac{T_\varOmega}{T_c} \sim 10^{-5}{-}10^{-6}, \quad Re = \frac{T_\nu}{T_c} \sim 10^8{-}10^9, \quad Rm = \frac{T_\eta}{T_c} \sim 10^3, \end{gather}

(8.2a–c)\begin{gather}Pm = \frac{T_\eta}{T_{\nu}} \sim 10^{-5}, \quad Pr = \frac{T_\kappa}{T_\nu} \sim 1, \quad E = \frac{T_\varOmega}{T_\nu}\sim 10^{-15}, \end{gather}

where $Ro$ is the Rossby number and $E = Ro/Re$ is the Ekman number. Note that most of these numbers are either very small or very large, reflecting the large separation of time scales. However, the magnetic Reynolds number, the source of all the problems for stellar and galactic dynamos, is not prohibitively large. Therefore we do not anticipate too much difficulty in the solution of the induction equation and our attention should be focussed on the momentum equation.

The vast range of temporal scales provides the most significant obstacle to the construction of a theoretical framework for the geodynamo. However, extra difficulties arise from a lack of knowledge of the strength and the form of the magnetic field in the core of the Earth, which leads to some strong debate about the nature of the important terms in the balance of the momentum equation.

8.2. Geostrophic and magnetostrophic balance

Consider the momentum (2.1), but now in a frame rotating with angular velocity ${\boldsymbol \varOmega }$ and where the body force ${\boldsymbol F}$ arises from buoyancy in a gravity ${\boldsymbol g}$. In this case the equation takes the form

(8.3) \begin{equation} \left.\begin{array}{cccccc@{}} \rho\bigg(\dfrac{\partial {\boldsymbol u}}{\partial t} + {\boldsymbol u} \boldsymbol{\cdot} {\boldsymbol \nabla} {\boldsymbol u}& + 2 {\boldsymbol \varOmega} \times {\boldsymbol u} \bigg) =& - \boldsymbol{\nabla} p& +\ {\boldsymbol j} \times {\boldsymbol B}& +\ \rho \nu \nabla^2 {\boldsymbol u}& +\quad \rho {\boldsymbol g}\\ Inertia& Coriolis& Pressure & Lorentz & Viscous & Buoyancy.\end{array} \right\}. \end{equation}

Here, the pressure terms also includes a contribution from the centrifugal acceleration. In the Earth's interior it is well established that, on time scales of relevance to the dynamo, the inertial and viscous terms do not contribute to the leading-order balance of terms in the momentum equation owing to the smallness of the Rossby and Ekman numbers, respectively. The Coriolis and pressure terms certainly do enter into the leading-order balance, with the buoyancy force acting to drive the dynamics. The debate surrounds the role of the magnetic field in the balance, since this arises from the dynamo generated field and cannot be theoretically determined a priori.

The discussion below follows similar lines of argument to those given in Nataf & Schaeffer (2015) and Aurnou & King (2017) and see the recent papers of Schwaiger, Gastine & Aubert (2019) and Aubert (2019). The importance of the Lorentz force relative to the Coriolis force can be determined by calculating the ratio of the two, given by

(8.4)\begin{equation} \varLambda_{C} = \frac{|\kern1.5pt {\boldsymbol j} \times {\boldsymbol B}|}{| 2 \rho \varOmega \times {\boldsymbol u} |} \sim \frac{B^2}{2 \mu \rho \ell_B \varOmega u}, \end{equation}

where $B$ is a typical magnetic field strength and the magnitude of the current density has been estimated as $|\kern1.5pt {\boldsymbol j}|\sim B/ (\mu \ell _B)$; here, the subscript $C$ is (somewhat provocatively) used for ‘correct’ as this measure gives the correct ratio of the relative importance of the Lorentz and Coriolis forces; $\varLambda _C$ is often termed the modified Elsasser number. Of course the typical amplitude of $B$ and $u$ is a scale-dependent quantity that emerges from the dynamo calculations, but it is of interest to estimate the typical size of this ratio for the measured field (${\sim }40\ \mathrm {\mu }$T) and inferred flow ($10^{-3}\ \textrm {m}\ \textrm {s}^{-1}$) in the Earth's core. For these values $\varLambda _{C} \sim 0.05{-}0.1$. This would seem to indicate that the Lorentz force, at least at the largest scale, is subdominant to the pressure and the Coriolis terms, and that a so-called leading-order geostrophic balance occurs where

(8.5)\begin{equation} \left.\begin{array}{ll@{}} 2 \rho {\boldsymbol \varOmega} \times {\boldsymbol u} \sim & - \boldsymbol{\nabla} p \\ Coriolis& Pressure \end{array}\right\}. \end{equation}

Here, the mathematical symbol $\sim$ is used in the loose sense of being of the same order as (or balancing).

Given the difficulty of estimating the typical amplitudes of the magnetic fields and flows (and indeed the currents) as a function of spatial scale, progress is often made by making one more assumption in the induction equation to link the three. If the fluid is assumed to be at asymptotically low $Rm$ then one may approximate the Lorentz force as

(8.6)\begin{equation} {\boldsymbol F}_L \sim \frac{1}{\mu}\overline{{\boldsymbol B}}\boldsymbol{\cdot} {\boldsymbol \nabla} {\boldsymbol b'}, \end{equation}

where $\overline{{\boldsymbol B}}$ is the large-scale field and ${\boldsymbol b'}$ is the small induced field given by

(8.7)\begin{equation} 0 \sim \overline{{\boldsymbol B}}\boldsymbol{\cdot} {\boldsymbol \nabla} {\boldsymbol u} + \eta \nabla^2 {\boldsymbol b'}, \end{equation}

so that $B \sim \overline{B} u \ell _B/ \eta$. The attractiveness of this approach becomes clear on realising that, on making the approximation, all the unknown terms disappear in the ratio (8.4) and for this limit the ratio of the Lorentz term to the Coriolis term is given by the Elsasser number

(8.8)\begin{equation} \varLambda = \frac{\overline{B}^2}{2 \mu \rho \eta \varOmega }. \end{equation}

Here, $\overline{B}$ is the typical size of the large-scale magnetic field and can be measured, whilst all the other quantities are known (or semi-known) a priori. $\varLambda$ is known as the Elsasser number and is nearly always quoted as the measure of the relative importance of the Lorentz to Coriolis force. In the Earth's core $\varLambda \sim 1$. Taken at face value this would place the Earth's core in magnetostrophic balance at largest scales, i.e.

(8.9)\begin{equation} \left.\begin{array}{ccc@{}} 2 \rho {\boldsymbol \varOmega} \times {\boldsymbol u} \sim & \, - \boldsymbol{\nabla} p & + \, {\boldsymbol j} \times {\boldsymbol B} \\ Coriolis & Pressure & Lorentz \end{array}\right\}. \end{equation}

Even though estimates based on the modified Elsasser number $(\varLambda _C)$ give the leading-order balance as geostrophic at the largest scale, it is clear from (8.4) that the importance of the Lorentz force in the momentum equation is scale dependent. At intermediate or small scales, the precise definition of which depends on the assumed form of the spectrum for the velocity and magnetic field, it is likely that magnetostrophic balance is re-established at a scale for which $Rm \sim 1$ ($\sim 10^4 {-} 10^5$m) in the Earth's core (Nataf & Schaeffer 2015; Aurnou & King 2017). Of course, because $Pm$ is small, this length scale is much larger than the viscous scale (and indeed the Rossby scale). For scales smaller than this one enters a diffusive regime with a strong applied field. In that regime, for scales larger than the Rossby scale the flow is heavily influenced by rotation, whilst for smaller scales the dynamics is prototypical diffusive MHD (Nataf & Schaeffer 2015); the kinetic energy spectrum $E_K(k) \sim k^{-3}$, whilst the magnetic energy spectrum is even steeper $E_M(k) \sim k^{-5}$.

It is worth mentioning here some theoretical arguments in favour of magnetostrophic balance in the geodynamo. Much of these arise from considerations of the linear theory of magnetoconvection, both in a plane layer and spherical shell. For magnetoconvection studies the magnetic field is imposed and the critical Rayleigh number $Ra_c$ is calculated as a function of rotation rate and magnetic field strength. These are measured by the non-dimensional input parameters $E = \nu /2 \varOmega L^2$, $Q = B_0^2/\mu \rho \nu \eta$ or $\varLambda = Q\,E = B_0^2/(2 \mu \rho \eta \varOmega )$. Briefly, as noted in Chandrasekhar's monograph (Chandrasekhar 1961), both rotation and magnetic fields separately act so as to suppress convection, with (in a plane layer) $Ra_c \sim E^{-4/3}$ as $E \rightarrow 0$ so that $Ra_c$ gets large; the associated critical wavenumber $k_c \sim E^{-1/3}$ so that the preferred mode is at small scales. In the non-rotating, magnetised case $Ra_c \sim Q$ as $Q \rightarrow \infty$ so that $Ra_c$ gets large with large field. Here also the associated mode is at small scales so that $k_c \sim Q^{1/6}$. Remarkably when both strong rotation and magnetic fields are present, the preferred mode has $\varLambda \sim 1$ and neither $Ra_c$ nor $k_c$ are asymptotically large. The magnetic field breaks the constraints imposed by rotation to allow efficient convection, there is a sweet spot at $\varLambda \sim 1$ in parameter space for linear convection to onset.

The existence of this sweet spot has formed the basis of the argument for the adoption of magnetostrophy in dynamo calculations. It is argued that, all other things being equal, the dynamo will try to locate itself in a state in which convection is optimised. This is certainly a plausible argument, although it should be stressed that it is not necessarily the case that convection is optimised in the nonlinear regime for the same parameters that optimise it in the linear regime; recall that investigating the linear regime for magnetoconvection is implicitly a low $Rm$ calculation. Moreover the breaking of the rotational constraint does not need to take place at all scales to allow efficient convection, allowing the possibility that geostrophic balance may proceed at large scales whilst magnetostrophy reigns at smaller scales.

At this point we note that it now sometimes traditional to designate a system as being in MAC (Magneto–Archimedes–Coriolis) balance when the primary balance is geostrophic (i.e. Coriolis force balanced by pressure) but the component of the Coriolis force that is not balanced by pressure is balanced by buoyancy and the Lorentz force at the next order, this then takes the form of balancing

(8.10)\begin{equation} \left.\begin{array}{cccc@{}} 2 \rho \,{\boldsymbol \varOmega} \times {\boldsymbol u} \sim & - \, \boldsymbol{\nabla} p & + \, {\boldsymbol j} \times {\boldsymbol B} & + \, \rho {\boldsymbol g}\\ Coriolis & Pressure & Lorentz & Buoyancy \end{array}\right\}. \end{equation}

It is interesting to note that if (8.10) holds then the magnetic field generated should be a function of the modified Rayleigh number ($\tilde {Ra} = (\alpha {\rm \Delta} T L)/\eta \varOmega$) only (Dormy, Oruba & Petitdemange 2018).

Recently Schwaiger et al. (2019) have argued that this balance should be termed QG-MAC (for quasi-geostrophic-MAC) balance to indicate that the primary balance is geostrophic whilst the Lorentz force enters in at the next order, although here we shall use the more usual terminology. However, it is fair to say that there is no overwhelming evidence or argument for assuming a particular force balance in the Earth's core. The investigation of these balances, in the relevant parameter regime of low $Pm$, is a very worthwhile topic of future investigations. Although this is difficult with current computational resources for the dynamo problem, much can be learned by studying balances in the nonlinear rotating magnetoconvection problem (see e.g. Stellmach & Hansen 2004).

8.2.1. Consequences of MAC balance: Taylor's constraint

Enforcing MAC balance imposes a serious constraint on the nature of the possible dynamo solutions. Here, we derive the simplest constraint on the form of the integrated Lorentz force known as Taylor's constraint.

Taking the $\phi$-component (in cylindrical polar coordinates $(s,\phi ,z)$ of (8.10) and integrating over cylinders aligned with the rotation axis one obtains

(8.11)\begin{equation} \int 2 \rho ({\boldsymbol \varOmega} \times {\boldsymbol u})_\phi \,\textrm{d} S \sim - \int (\boldsymbol{\nabla} p)_\phi \,\textrm{d} S + \int (\kern1.5pt {\boldsymbol j} \times {\boldsymbol B})_\phi \, \textrm{d} S + \int (\rho {\boldsymbol g})_\phi \, \textrm{d} S,\end{equation}

where $\textrm {d} S = s \, \textrm {d} \phi \, \textrm {d} z$. As buoyancy is purely radial the final term is zero and the pressure term integrates to zero. Finally the Coriolis term can also be shown to be zero by using the divergence theorem over the cylinder. Hence, the only term to contribute to (8.11) is the Lorentz force and so

(8.12)\begin{equation} \int (\kern1.5pt {\boldsymbol j} \times {\boldsymbol B})_\phi \, \textrm{d} S = \frac{1}{\mu} \int (({\boldsymbol{\nabla} } \times {\boldsymbol B}) \times {\boldsymbol B})_\phi \, \textrm{d} S = 0.\end{equation}

This Taylor constraint ensures that the Lorentz torque on any cylindrical surface parallel to the rotation axis is zero.

It is possible to construct dynamos that solve the induction equation in addition to maintaining Taylor's constraint both in two (Roberts & Wu 2018) and three (Li, Jackson & Livermore 2018) dimensions. This elegant formalism certainly provides a way forward for finding exact solutions that satisfy MAC balance. Numerical computations usually do, however, solve the full momentum equation (albeit in the wrong parameter regime) and utilise the size of the departures from the solution to (8.12) (the degree of ‘Taylorisation’) as a measure of the degree to which they have been successful in obtaining solutions that are in MAC balance.

8.3. Computational models of the geodynamo

For all the reasons described above, computational modelling of the geodynamo remains a formidable venture. The modern era of massively parallel computational geodynamo modelling was started by the seminal work of Glatzmaier & Roberts (1995). They simulated the interaction of rotating, anelastic convection with magnetic fields to produce dynamo models that bore striking resemblance to the observed magnetic field; they even managed to observe reversing magnetic fields. That this tour de force calculation was so successful is surprising given the large discrepancy between the parameters that could be reached computationally and those that pertained to the Earth. With the increase in computational power, a large number of subsequent papers have attempted to refine and improve on the Glatzmaier and Roberts models. The culmination of these efforts so far are two recent high resolution models (Aubert, Gastine & Fournier 2017; Schaeffer et al. 2017).

The first of these models (Aubert et al. 2017) uses large-eddy simulations for the momentum equation (although the induction equation is fully resolved). The parameters are selected to get as close as computationally possible to the values in the Earth whilst maintaining a QG-MAC balance and keeping a constant $Rm$. This is achieved by scaling the input parameters with a parameter $\epsilon$ (in this case the convective power) in a consistent manner – this is similar in spirit to the distinguished limit calculations pioneered by Dormy (2016) discussed below. Here, in the simplest set-up, $Pm \propto \epsilon ^{1/2}$, $E\propto \epsilon$, $Pr \propto 1$, $Rm \propto 1$. Here, $\epsilon \sim 1$ corresponds to the conditions found successfully to produce MAC balance in a moderately rotating model, whilst $\epsilon$ small is appropriate for the Earth. The impressive numerical calculations are performed for a range of parameters on the path, $\epsilon =1$ corresponds to $E \sim 3 \times 10^{-5}$ and $\epsilon = 3 \times 10^{-4}$, gives $E = 10^{-8}$. As expected, the solutions take the form of a MAC balanced dynamo field following diffusivity-free power-based scaling laws. Schaeffer et al. (2017) also consider numerical simulations of dynamo action in a similar parameter regime ($E \sim 10^{-7}$, $Pm \sim 0.1$, $Rm > 500$). They find magnetic fields in MAC balance, with the magnetic energy one order of magnitude larger than the kinetic energy. Interestingly, they find inhomogeneity of the solution with a dynamical contrast between the interior and the exterior of the cylinder parallel to the axis of rotation that circumscribes the inner core, the so-called tangent cylinder. Inside the tangent cylinder, the strong magnetic field is generated by a polar vortex. Outside the tangent cylinder, however, the kinetic energy is mostly non-zonal, with zonal winds being suppressed by the Lorentz force. At different spatial and temporal scales the flow may be either geostrophic (for large-scale and low frequency flows and for small-scale convection in weak field regions) or in MAC balance (for high frequency large-scale modes or convective modes at intermediate scales).

8.3.1. The search for a distinguished limit

As noted above, geodynamo models are usually integrated in completely the wrong parameter range, but sometimes yield solutions that resemble the geomagnetic field. When and why is this the case, and how can one continue to achieve solutions in MAC balance as computers get faster? It is tempting to try to compute numerical models with all the parameters set to be as close as possible to their correct geophysical values (in particular low $E$ and low $Pm$). Is this the most sensible approach, given current computational limitations?

An alternative, and more considered, approach was pioneered by Dormy (2016) who studied the properties of a large set of dynamo models at moderate Ekman number (with $E \approx {O}(10^{-4}$)). He argued that, at such moderate values amenable to rapid computation, it was necessary to increase $Pm$ from its order-one value to achieve the correct balance for the generated magnetic field; indeed it is also of interest to examine rapidly rotating dynamos where inertia is completely suppressed (Hughes & Cattaneo 2016) – in this regime $Pm$ is infinite. Of course, when $E$ is reduced, as computational resources become more lavish, one should move away from these models and $Pm$ should be reduced in tandem with the Ekman number. Dormy argued that a distinguished limit should be found where all parameters, including $Pm$ scale with the Ekman number as $E \rightarrow 0$ to leave one on the balanced dynamo branch. This is also the approach taken by Aubert et al. (2017), as described above, although the scalings selected for the variation of the parameters are slightly different.

8.4. Asymptotic models of the geodynamo

An alternative approach to performing DNS of the geodynamo is to construct asymptotic models valid in the regime of strong rotation (i.e. asymptotically small Rossby number $Ro$). This has much in common with the approach advocated by Dormy (2016) and Aubert et al. (2017). The difference here is that the models are predicated on the fact that the separation of time scales discussed earlier leads naturally to reduced dynamics. There are two reasons why one might seek an asymptotic solution to a complicated problem. The first is that performing an asymptotic expansion often leads to the development of a model where analytic progress is possible and an elegant solution emerges; there are many examples of this use of such models. The second reason is that it is the right thing to do for a given problem; performing the expansion does not lead to a reduction in the complexity of the problem (for example a 3-D partial differential equation may stay a 3-D partial differential equation) but performing the expansion does, by design, lead to a set of reduced equations that no longer contain the small parameters, which is therefore more amenable to solution. Such models have sometimes been viewed with suspicion, ironically sometimes even by researchers happy to utilise the Boussinesq approximation as an asymptotic framework for modelling.

The models for the geodynamo summarised here are of the second kind. The first multi-scale asymptotic models were constructed by Childress & Soward (1972), who performed a weakly nonlinear analysis of Boussinesq convective dynamos in a plane layer; their analysis showed how weak large-scale magnetic fields were readily generated. The Childress–Soward model has been developed, first by Fautrelle & Childress (1982) who extended the analysis to intermediate field strengths, and more recently to include the effects of multiple length scales perpendicular to the rotation axis to enable the description of strongly forced convection. In a series of papers (see Calkins et al. 2015; Calkins, Julien & Tobias 2017; Calkins 2018; Plumley et al. 2018 and the references therein) the asymptotic theory for modelling of quasi-geostropic dynamo models (QGDM) was developed.

Here, we give a brief description of the derivation of the simplest form of the QGDM in a Cartesian plane layer $(x,y,z)$ with gravity pointing in the (negative) $z$-direction. The starting point of the derivation are the non-dimensional equations of Boussinesq magnetised convection in a plane layer of depth $H$, rotating about the vertical $z$-axis. These are given by

(8.13)\begin{gather} \left(\frac{\partial {\boldsymbol u}}{\partial t} + {\boldsymbol u} \boldsymbol{\cdot} \boldsymbol{\nabla} {\boldsymbol u} \right)+\frac{1}{E}\hat {\boldsymbol z} \times {\boldsymbol u} = - \boldsymbol{\nabla} \varPi + \frac{\varLambda}{Pm E} \, {\boldsymbol B} \boldsymbol{\cdot} \boldsymbol{\nabla} {\boldsymbol B} + \frac{Ra}{Pr}\, T + \nabla^2 {\boldsymbol u} , \end{gather}

(8.14)\begin{gather}\frac{\partial T}{\partial t}+{\boldsymbol u} \boldsymbol{\cdot} \boldsymbol{\nabla} T = \frac{1}{Pr}\nabla^2 T, \end{gather}

(8.15)\begin{gather}\frac{\partial {\boldsymbol B}}{\partial t} = \boldsymbol{\nabla} \times ( {\boldsymbol u} \times {\boldsymbol B} )+ \frac{1}{Pm}\nabla^2 {\boldsymbol B}, \end{gather}

(8.16)\begin{gather}\boldsymbol{\nabla} \boldsymbol{\cdot} {\boldsymbol B} = 0, \end{gather}

(8.17)\begin{gather}\boldsymbol{\nabla} \boldsymbol{\cdot} {\boldsymbol u} = 0. \end{gather}

Here, length is non-dimensionalised with respect to $H$, time with the viscous diffusion time across the layer $H^2/\nu$, giving a typical velocity scale of $\nu /H$. In addition magnetic fields are non-dimensionalised with a typical field strength $\overline{B}$ and temperatures with the temperature drop across the layer (${\rm \Delta} T$). In addition to the non-dimensional parameters given in (8.2a–c) and (8.8) the ration of the thermal driving to diffusive processes is given by the familiar Rayleigh number,

(8.18)\begin{equation} Ra = \frac{\alpha {\rm \Delta} T g H^3}{\nu \kappa}, \end{equation}

with $\alpha$ here being the coefficient of thermal expansion.

Just as for the Childress–Soward dynamo model, the anisotropic structure of rotating convection is exploited; $H/L$ is set to be large, where $H$ is the vertical scale of the convection and $L$ is the horizontal scale. We set $H/L = \epsilon ^{-1}$, with $\epsilon$ a small parameter. In general, the degree of anisotropy is expected to be related to the strength of rotation as measured by the Rossby number of the convective turbulence, and so a natural choice is to set $\epsilon = Ro$; this is the scaling utilised by Calkins et al. (2015). However, the Rossby number is an output parameter of a convective system and can only be guaranteed to be small a priori in certain circumstances. It is useful to consider a case where the separation of scales is determined by an input parameter. In the linear regime of rotating convection, anisotropy is governed by the Ekman number and $\epsilon = E^{1/3}$ is a possible choice. The scaling of the anisotropy allows for the use of multiple-scale asymptotics in time and the axial space direction such that

(8.19a,b)\begin{equation} {\partial_z} \rightarrow {\partial_z} + \epsilon {\partial_Z} , \quad {\partial_t} \rightarrow {\partial_t} + \epsilon^{3/2} {\partial_{\tau}} + \epsilon^2 {\partial_{T}}. \end{equation}

Hence $Z = \epsilon z$ is the large-scale vertical coordinate over which convection occurs and $\tau = \epsilon ^{3/2} t$ and $T=\epsilon ^2 t$ are the ‘slow’ mean magnetic and ‘extremely slow’ mean temperature time scales. Note that these time scales are in the correct ordering as described in § 8.1. The slow and fast independent variables are therefore denoted by $(Z,\tau ,T)$ and $( {\boldsymbol x},t )$, respectively. In the more general derivation of Calkins et al. (2015), slow horizontal scales are also utilised.

All of the dependent variables are decomposed into mean and fluctuating variables according to

(8.20)\begin{equation} f = \overline{f}(Z,\tau,T) + f'({\boldsymbol x},Z,t,\tau,T), \end{equation}

with the ‘fast’ averaging operator such that $\overline {f'} \equiv 0$ , Each variable is then expanded in a power series according to

(8.21)\begin{align} f({\boldsymbol x},Z,t,\tau,T) &= \overline{f_0}(Z,\tau,T) + f'_0({\boldsymbol x},Z,t,\tau,T)\nonumber\\ &\quad +\epsilon^{1/2} [\kern1.5pt \overline{f_{1/2}}(Z,\tau,T) + f'_{1/2}({\boldsymbol x},Z,t,\tau,T) ] \nonumber\\ &\quad+ \epsilon [\kern1.5pt \overline{f_{1}}(Z,\tau,T) + f'_{1}({\boldsymbol x},Z,t,\tau,T) ] + O(\epsilon^{3/2}) , \end{align}

and substituted into the governing equations. The equations are then separated into mean and fluctuating components. The leading-order balance in the fluctuating momentum equation is geostrophy and closure of the system is obtained by imposing so-called solvability conditions to eliminate secular growth. This allows the derivation of prognostic equations for the vertical vorticity and vertical velocity, respectively. Similar techniques are utilised for the temperature and induction equations.

In the most realistic of these models both the Rayleigh number and magnetic Prandtl numbers are also scaled so that we set, ${Ra} \sim \epsilon ^{-4}$ and $Pm \sim \epsilon ^{1/2}$. Importantly, the magnetic field is also scaled so that magnetic energy is asymptotically larger than the kinetic energy, which is physically realistic. The method of multiple scales is used to derive 3-D reduced PDEs for large-scale (slow) and small-scale (fast) quantities.

As noted above, the dynamics is such that the QGDM is geostrophic to leading order, with the strong magnetic field entering into the prognostic equation that determines the dynamics. This dynamics in plane-layer models is beginning to be investigated in both the kinematic and nonlinear dynamo regimes. Hydrodynamically, quasi-geostrophic convection in a plane layer passes through four typical configurations, (a) cellular convection, (b) columnar convection, (c) plume convection and (d) geostrophic turbulence as the Rayleigh number increases. Figure 10 shows the vorticity for the final two such asymptotic solutions in the kinematic regime, and the small-scale current generated by the corresponding kinematic dynamo. Interestingly, all the convection regimes act as large-scale dynamos, with the properties of the large-scale magnetic field remaining relatively insensitive to the precise form of the convection. These models have recently been extended into the nonlinear regime by Plumley et al. (2018).

Figure 10. Volumetric renderings of the small-scale vertical vorticity (a,b) and current density (c,d) illustrating two flow regimes observed in the reduced convection simulations, namely the plume regime $Pr=10$, $\tilde {Ra}=200$ (a,c) and the geostrophic turbulence regime (b,d) $Pr=1$, $\tilde {Ra}=100$. Here gravity points downwards whilst rotation is about a vertical axis. After Calkins et al. (2016).

It is to be hoped that these models will be able to give some insight as to the behaviour of rapidly rotating planetary dynamos at low Rossby and magnetic Prandtl number. It will be interesting to utilise other choices for the small parameter that leads to a large anisotropy ratio; these choices may simply be guided by the results of numerical simulations (see e.g. Aubert 2019). Of crucial importance, of course, is to extend these models to more realistic geometries such as spherical shells and thermal annuli, and to investigate the dynamics of reduced models with other dynamic balances.

9. A word on dynamo experiments

Although this perspective is concerned with the theory of turbulent dynamos, it would be remiss to conclude without some discussion of the amazing dynamo experiments that have been performed over the past thirty or so years. Primarily these have been constructed so as to address the question of the nature of dynamo action near to onset, i.e. for $Rm$ close to $Rm_c$.

A typical laboratory dynamo is usually designed to drive a large-scale flow with desirable laminar dynamo properties; i.e. they are usually optimised so that $Rm_c$ is as small as possible; even obtaining a working laboratory dynamo is a formidable challenge. For example, the Riga dynamo experiment (see e.g. Gailitis et al. 2002), is designed to mimic the Ponomarenko dynamo described in § 2.5.1. This dynamo was the first experimental realisation of self-excitation in a liquid metal flow, and took the form of a $3$ metre long helical sodium flow driven by a propeller. The dynamo growth rates and frequencies from the experiments agreed well with numerical predictions (to within $5\,\%\text {--}10\,\%$). Interestingly, computation seems to indicate that de-optimising the azimuthal velocity may lead to the introduction of interesting saturation effects.

The Karlsruhe experiment (see e.g. Stieglitz & Müller 2001) was largely based on the Roberts dynamo of § 2.5.1. Here, a liquid sodium flow was driven in an array of columnar helical pipes, confined in a cylindrical container. The dynamo experiments demonstrated that kinematic dynamo action occurs, with a well-predicted growth rate. The saturated magnetic field oscillated about a well defined mean value for fixed flow rates. The experiments could detect two quasi-dipolar magnetic fields of opposite polarity. For both the Riga and Karlsruhe dynamos, the flows are tightly constrained and remain close to the targeted laminar flows for optimal dynamo action; perhaps it is no surprise that they worked so well.

All of the experimental dynamos utilise liquid metals for which, $Pm \sim { {O}} (10^{-6} {-}10^{-5}) \ll 1$, as discussed earlier. Hence, even though $Rm_{c}$ is moderate. the corresponding Reynolds numbers are vast. For constrained flows, such as those of the Riga and Karlsruhe dynamo experiments, this is not such a problem. However, the presence of turbulence at a high $Re$ does have a significant impact on the experimental search for dynamo action in less constrained flows such as the VKS (von Kármán sodium) experiment (Monchaux et al. 2007). There the flow of liquid sodium was driven in a cylindrical container by rotating propellers. The desired large-scale ‘French washing machine’ (von Kármán) flow is known from computation to be a decent dynamo, and the experiment was optimised so that the (axisymmetric) mean flow provided the lowest possible $Rm_c$. As this mean flow is axisymmetric, the magnetic field generated is expected to have a non-axisymmetric $m=1$ structure. However, because $Re$ is large, this flow is unstable and turbulent fluctuations are found that are comparable in magnitude with the mean flow. As discussed in § 4.2.2, this has a negative effect on the dynamo efficiency of the experiment, increasing $Rm_c$ and rendering dynamo action impossible unless magnetic material is used for the driving disks. Interestingly, when magnetic field growth is sustained, it is an (axisymmetric) $m=0$ magnetic field that dominates; the experiment is a superb example of mean-field electrodynamics in action. Here, a non-axisymmetric turbulent flow driven by vortices near the impellers provides an $\alpha$-effect that interacts with the shear generated near the rotating disks to yield a mean-field dynamo of $\alpha \omega$-type (Ravelet et al. 2012). Although this experiment does rely on the presence of magnetic impellers for dynamo action, it is fair to state that their presence simply acts so as to reduce the $Rm_c$, and that such a set-up would produce a dynamo should higher $Rm$ be reachable. Furthermore, once this material is used the system undergoes a bewildering array of nonlinear dynamics, including the generation and reversals of a large-scale field; much can be learned here of the nonlinear behaviour of large-scale magnetic fields in turbulent systems.

It is important to also mention the ‘whirling dervish’ experiments of the Grenoble group (see e.g. Brito et al. 2011). This liquid sodium experiment investigates magnetised spherical Couette flow, with a strongly magnetised inner core. Although the field is imposed, so that this is not strictly a dynamo experiment, much has been learned about the possible balances and dynamics of waves in rapidly rotating MHD systems from this experiment. Next, and extremely impressively, turbulent transport coefficients have been calculated in liquid metal experiments in a sphere by the Madison group (Rahbarnia et al. 2012). Although of course these experiments are perforce at low $Rm$, they give a crucial insight into the behaviour of driven low $Pm$ fluids.

Dynamo experiments are currently being constructed or modified in both Dresden and Maryland. The Dresden experiment takes the form of a precession driven dynamo, with $8$ m$^3$ of liquid sodium being rotated at angular velocities of up to $10$ Hz, to give $Rm \approx 700$ (Stefani et al. 2015). Incredibly the experiment will generate a gyroscopic torque of $8$ MNm on the foundations of the building. Such experiments are not undertaken lightly.

In general though, it is important to reflect that it is very difficult to have computational and experimental models of dynamos in the same regime. DNS likes to sit close to the $Pm=1$ boundary, whilst liquid metal experiments are constrained to be at low $Pm$. For dynamo theorists interested in data-driven modelling, this may pose some issues, as discussed below.

10. Discussion: what of the future?

‘I may not have gone where I intended to go, but I think I have ended up where I intended to be’.

Douglas Adams

The previous discussion has emphasised that the current status of the theory of turbulent dynamos is uncertain, both for the geodynamo and for its astrophysical counterparts. Despite the heroic efforts to solve the relevant equations numerically, computational efforts are limited to parameter regimes that are far from those that pertain to naturally occurring dynamos. Perhaps this is best emphasised by a recent calculation by P. Käpylä (private communication). He demonstrates that, even if the computational resources were available to undertake numerical simulation of turbulent dynamos, the power required to simulate a star such as the Sun is $10^{22}$ W. This is equivalent to the luminosity of a M9V main sequence red dwarf. Although this is a fairly cool star, the point is well made.

The lack of (computational) power for numerical simulations manifests itself as a different type of deficiency for geodynamo models from astrophysical dynamo models. For astrophysical models, the key point is that the flows that can be simulated are too laminar and too diffusive and there is also not enough separation of scales between the resistive and viscous scales. Hence, for astrophysical dynamo theorists even solving the induction equation is a formidable task, which is beyond the computational power of the near future. For planetary dynamos, such as the geodynamo, the numerical solutions are too viscous and too slowly rotating (although the low $Pm$ problem is also manifest here). It is very difficult to achieve solutions that are in the relevant asymptotic regime of geostrophy or magnetostrophy.

So what is the correct way to proceed? For planetary dynamos, since the viscosity is known to play a minor role, it is useful to examine models that do not rely on viscosity, such as those that are driven by waves (see e.g. Davidson & Ranjan 2018) or those that saturate via the enforcement of Taylor's constraint – or a modified Taylor's constraint that includes the temporal variation of the zonal flow – (see e.g. Wu & Roberts 2015; Li et al. 2018). Another promising avenue to pursue is to determine a distinguished limit (or a path through parameter space) that maintains the correct balances, even though the computational resources are not adequate currently to simulate the vast range of temporal scales; this is the approach taken by Dormy (2016) and Aubert et al. (2017). Furthermore, given the vast separation in spatial and temporal scales for the Earth's interior, a programme based on the derivation and solution of asymptotic models is warranted (Calkins 2018). This programme will require the development of a new class of time stepping methods that are capable of efficiently integrating reduced equations that evolve on very different time scales (i.e. slow–fast dynamics). One issue that remains puzzling has to be the mechanism that leads to the reversals of the Earth's field. Presumably such an abrupt transition relies on nonlinear effects (see e.g. Jones 2008). One cannot rely on inertia, which is small, to be the source of nonlinearity, so that leaves potentially the Lorentz force or the temperature advection and there are constraints on the action of both of these.

For the Sun and other astrophysical bodies, the problem remains the lack of techniques to deal with flows on a vast range of spatial scales. Clearly any numerical code is capable of only resolving an extremely limited range of these scales. There is little prospect in the near future of developing numerical algorithms that are capable of simulating these bodies directly. Even methods designed specifically for simulating a large range of spatial scales, such as those that employ adaptive-mesh refinement, are not suitable for turbulent flows, where high resolution is required globally rather than over a small fraction of the domain. Some procedure must therefore be developed for the parameterisation of the low-order statistics of the scales below the smallest that can be resolved. Such sub-grid models are extremely difficult to construct for MHD turbulence (Miesch 2015). Moreover, we have stressed here that it is key that such parameterisations (even if couched in terms of the construction of turbulent transport coefficients) respect the inherent conservation laws of the original system in the relevant limit. As we have shown, for example, conservation of magnetic helicity plays a vital part in constraining the dynamo solutions. We discuss two such possible methods for making progress in both the planetary and astrophysical dynamo problem here. The first is to utilise DSS via cumulant expansions (see e.g. Marston, Qi & Tobias 2019).

10.1. Direct statistical simulation

The DSS approach involves the direct solution of differential equations derived for the low-order statistics of fluid flows (such as mean flows and two-point correlation functions). In this methodology the equations take the form of a truncated (equal-time) cumulant expansion that respects the inhomogeneity and anisotropy of the underlying dynamo system. This approach is therefore a generalisation of the methods discussed in § 6.4 to anisotropic and inhomogeneous systems. Because of the increased complexity, a computational approach is required to facilitate this approach (see e.g. Tobias, Dagon & Marston 2011b). Details of the technical procedures utilised are given in Marston et al. (2019); here we omit the details, but give a schematic representation.

Briefly, the method proceeds by considering the system as a dynamical system and deriving evolution equations for the cumulants of the probability distribution function of the turbulent state. This is much less computationally expensive than solving the Fokker–Planck equation (or the equivalent Hopf equation) for the whole p.d.f. although that approach may also be worth pursuing (see e.g. Venturi 2018).

Briefly, consider a dynamical system where the state variable $\boldsymbol {q}(\boldsymbol {x},t)$ (in the case of dynamos this may include information about the velocity, pressure, magnetic field, temperature etc of the fluid) evolves via

(10.1)\begin{equation} \boldsymbol{q}_t = \mathcal{L}(\boldsymbol{q}) + \mathcal{N}(\boldsymbol{q,q}) + \boldsymbol{f}(\boldsymbol{x},t). \end{equation}

Here, $\mathcal {L}$ is a linear operator, whilst $\mathcal {N}$ is a nonlinear (in the simplest case quadratic) operator and $\boldsymbol {f}(\boldsymbol {x},t)$ is an (often stochastic) driving term.

The method proceeds by defining the low-order cumulants which, up to third order, are given by

(10.2)\begin{gather} c_1(\boldsymbol{r}) = \langle \boldsymbol{q} \rangle, \end{gather}

(10.3)\begin{gather}c_2(\boldsymbol{r}_1,\boldsymbol{r}_2) = \langle \boldsymbol{q}'(\boldsymbol{r}_1)\, \boldsymbol{q}'(\boldsymbol{r}_2) \rangle, \end{gather}

(10.4)\begin{gather}c_3(\boldsymbol{r}_1,\boldsymbol{r}_2,\boldsymbol{r}_3) = \langle \boldsymbol{q}'(\boldsymbol{r}_1)\, \boldsymbol{q}'(\boldsymbol{r}_2)\, \boldsymbol{q}'(\boldsymbol{r}_3)\rangle, \end{gather}

so that to this order the first cumulant is given by the mean. Here, as before, the averaging process may be taken over a coordinate or be defined as an ensemble mean. The second and third cumulants are the centred moments (or two and three point correlation functions). They may also be thought of as the functional derivatives of the Hopf functional.

Equations for the evolution of the cumulants may be derived either via brute force (potentially utilising symbolic manipulators) or via the Hopf functional technique. Schematically these take the form of a cumulant hierarchy

(10.5)\begin{gather} \frac{\partial {c_1}}{\partial t} = \mathcal{L}(c_1) + \mathcal{N}(c_1^2 +c_2), \end{gather}

(10.6)\begin{gather}\frac{\partial {c_2}}{\partial t} = \mathcal{L}(c_2) + \mathcal{N}(c_1 c_2 +c_3) + \varGamma, \end{gather}

(10.7)\begin{gather}\frac{\partial {c_3}}{\partial t} = \mathcal{L}(c_3) + \mathcal{N}(c_1 c_3 + c_2 c_2 + c_4). \end{gather}

Here, $\varGamma$ encodes the net action of the stochastic driving on the second cumulant (see e.g. Tobias et al. 2011b for more details). As for classical moment hierarchies, the evolution equation for the $n$th cumulant $c_n$ depends on $c_{n+1}$. Hence, as for the models discussed in § 6.4, progress can only be made by truncating the hierarchy. For computational practicality, this is usually done by setting $\mathcal {N}(c_3)=0$, so that the effect of the third cumulant on the evolution equation for the second cumulant is discarded; setting $c_3 = 0$ is certainly sufficient to achieve this. If the cumulant equation is truncated at second order, the approach is known as CE2, which is formally equivalent to the stochastic structural stability theory of Farrell & Ioannou (2003). Truncated at this order the equations are realisable, in the sense that energy densities and probabilities are positive. Moreover, global conservation laws are respected for quadratic invariants (such as magnetic helicity and total energy of the MHD system). This can be shown by deriving the CE2 equations as an exact counterpart of the quasilinear (QL) approximation. CE2 is quasilinear in the sense that it includes interactions of mean quantities with eddies to give eddies and interactions of eddies with eddies to give means (either flows or fields), but neglects interactions of eddies with eddies to give eddies (which we term eddy–eddy scattering). In this sense, CE2 is the self-consistent, nonlinear (i.e. dynamic), inhomogeneous, anisotropic version of first-order smoothing.

Of course it is important to go beyond the traditional QL approximation for DSS. As evidenced by the discussion in § 6.4, a wealth of possible extensions of QL approaches are available.

Perhaps the most obvious route to take is to include one extra evolution equation in the hierarchy for the statistics. This alternative, which has been termed CE3, relies on discarding $c_4$. If this is the route to be taken then some steps must be taken to ensure the realisability of the system. This can be achieved by including an eddy-damping term in the equation for the third cumulant; this form of DSS is then the anisotropic, inhomogeneous version of EDQNM of § 6.4. A more sophisticated approach involves ensuring realisability by projecting out the eigenvectors associated with the unphysical negative eigenvalues of the second cumulant (Marston et al. 2019).

However, a naive implementation of CE3 may be prohibitively computationally expensive owing to the ‘curse of dimensionality’; the statistics typically have a higher dimension than the underlying fields. This limitation has been addressed in a recent paper (Allawala, Tobias & Marston 2020) that performs DSS in a reduced basis; the basis is calculated using a proper orthogonal decomposition of the eigenvectors of the second moment of the fields, retaining only those modes that are important in the evolution of the second cumulant. The procedure can be very effective in reducing the computational cost with limited loss of accuracy. This is a form of unsupervised learning, with training based upon full resolution simulations. An alternative may be to utilise a form of machine learning to represent the high-order statistics based on the state of the first and second cumulants; such an approach would go one step beyond current machine learning approaches that parameterise the relationship between the turbulent transport coefficients and the evolving mean fields.

An alternative to proceeding to higher order is instead to derive statistical equations that are the closures of the generalisation of the QL approximation. The generalised quasilinear (GQL) approximation is achieved by separating all the state variables into large and small scales via a spectral filter rather than by a decomposition into a formal mean and fluctuations (Marston et al. 2016). Nonlinear interactions involving only small scales are then removed, so that the evolution equation for the small scales remains formally linear. This facilitates the closure of the system in a generalised CE2 (or GCE2) method. For a range of fluid and MHD problems it has been shown that GQL performs in a superior fashion to the standard formal QL approximation (see e.g. Tobias & Marston 2017). Further investigation is needed as to the reason for the marked superiority of GQL over QL; however, it has been noted that GQL allows the transfer of energy between fluctuation modes of different scales by scattering off large-scale modes; this mechanism is disallowed in QL.

In the dynamo (and magnetorotational instability) context the simplest form of DSS, namely CE2, has been utilised in the landmark paper of Squire & Bhattacharjee (2015). In that paper the turbulence and dynamo induced by the magnetorotational instability (MRI) were analysed. They considered a case where homogeneous turbulence was susceptible to a large-scale dynamo instability; the instability saturated in an inhomogeneous equilibrium that depends on the value of the magnetic Prandtl number ($Pm$). Examples of the large-scale structure of the magnetic field from the direct statistical simulation, for two different values of $Rm$ is given in figure 11. What is interesting is that the dependence of the angular momentum transport on $Pm$ in the CE2 model is qualitatively similar to the nonlinear MRI turbulence found from DNS. It is also important to note that the use of DSS allowed for the analysis of cases at more extreme values of $Pm$ than were available to DNS. As stated by the authors, this suggests that these ‘$\ldots$ models may be used to gain insight into the astrophysically relevant regimes of very low or high $Pm$’. Finally for statistical models and methods, it is important to derive a firm theoretical foundation on which to build algorithms and computational models. As for pure turbulence and climate research, where statistical models are more widely utilised, there needs to be a program of research that attempts to answer the questions as to whether it possible to construct a non-equilibrium theory that describes the interaction of mean flows and magnetic fields with turbulence in inhomogeneous, anisotropic dissipative systems (such as dynamo systems). For example, is there a well-defined procedure for connecting fluctuations to dissipation in turbulent magnetised flows that are not homogeneous and isotropic (and may be stratified and rotating)? These will probably take the form of generalised fluctuation–dissipation theorems; progress in this area could lead to the development of meaningful statistical theories for the turbulent generation of magnetic field.

Figure 11. Evolution of the large-scale mode as a function of $z$ and $t$; at $Pm=4$ for (a) $Rm = 4500$, and (b) $Rm = 12\ 000$. In (a) the solution reaches a steady state, whilst in (b) DSS leads to a time-dependent state. After Squire & Bhattacharjee (2015).

10.2. A data-driven approach to dynamo theory?

In the past few years it has been suggested that progress in modelling many physical systems can be accelerated by utilising data-driven methods (see e.g. Spears et al. 2018; Brunton, Noack & Koumoutsakos 2020). These methods take advantage of the rapid advances in computational power and algorithmic development and can be used either independently or in conjunction with a physics-based approach. As noted by Kutz (2017), the two cultures of statistics and data science, namely machine learning and statistical learning, may both be of interest for the modelling of turbulent flows. In general, the former is concerned with prediction whilst the latter is concerned with ‘inference of interpretable models from data’ (sometimes termed model reduction).

How might these methods be brought to bear on the dynamo problems at hand? It is perhaps easier to visualise a route to success for the second approach. In fact, dynamo theorists through the application and evaluation of mean-field theory have been trying to construct reduced models for dynamos for over fifty years. As we have seen, in mean-field theory, the idea is to parameterise the effects of all the unresolved scales using transport coefficients in reduced models for the mean fields. Here there is a strict separation of scales and separate evolution equations can be determined. There has recently been an interesting application of machine learning to the problem of constructing turbulent transport coefficients by Nauman & Nättilä (2019). Here DNS of a helically forced dynamo system was performed. The data from this simulation were compared with a mean-field model, with the EMF parameterised as either a linear or nonlinear functional of the mean field. Other analytical models for the EMF were also compared with the data using Bayesian inference with Markov chain Monte Carlo sampling. From training on data, it was possible to construct the transport coefficients, although it was shown to be extremely difficult to constrain the formulae for the quenching of the transport coefficients by the mean field. Interestingly, the data-driven approach could also show that non-locality in time could play an important role in relating the transport coefficients to the mean field.

A central problem for the application of a data-driven approach to the dynamo problem arises quite simply because of the lack of available data to constrain the models. As we have seen, it is simply not possible to perform laboratory or numerical experiments at (or even close to) the relevant conditions for an astrophysical object (or indeed the Earth's dynamo). So what is the way forward with current computational resources?

Computational models are certainly capable of faithfully and efficiently simulating a range of spatial and temporal scales, although calculations will necessarily miss contributions to the evolution from the scales above and below those being simulated. Typically, as we have seen, the effect of the smaller spatial scales can be modelled utilising their statistical properties since these evolve on a faster time scale than the larger scales. A hierarchical multi-scale data-driven approach (HMDD) may therefore be effective. For this approach, the starting point is a simulation of the smallest scales at parameter values that are comfortable for the available computational resources; the availability of these resources will certainly limit the range of scales that can be simulated. These simulations will be conducted for a range of imposed larger-scale conditions (e.g. magnetic field and shear strength and orientation) and the data from these will be used for training. In this case the output of the training will be the statistics of the small scales for all larger scale conditions. Once this training has been achieved, the process can be repeated with a range of scales on a slightly larger length scale and longer time scale, with the statistics of the small scales parameterised by the output of the training. Of course, the parameters pertaining to the evolution of these larger scales will be different, as will the prevailing balances, so the statistics of these scales can only be found via the data from this simulation being utilised for training again. This process is repeated until a global-scale calculation is possible (with the accurate, learned statistics from the smaller scales repeatedly being utilised). Although, this approach may require many levels of training for a realistic geophysical or astrophysical problem, at each level of training the computational costs are designed not to be restrictive. This HMDD programme of research is technically challenging, but perhaps is able to marry the better aspects of physics driven and data-driven approaches together without ad hoc assumptions. It is clear that data-driven methods will play a role in advancing our understanding of dynamos, but there are so little data (here laboratory experiments and astrophysical observations will play a key role in validation and verification of the procedures) that one must proceed with caution.

Finally, for data-driven methods in dynamos it is intriguing to note that machine learning is now at the point where learning techniques can be utilised to solve equations efficiently after training (Li et al. 2020) and that there are now techniques for constructing models by recovering reduced equation sets from data (Brunton, Proctor & Kutz 2016). Both of these advances, in conjunction with the necessary understanding of the allied physics and mathematical structure of dynamos, should lead to a significant leaps in our understanding and ability to model dynamos.

10.3. Final thoughts

Whatever the successful approach turns out to be, I believe it will not emerge from a blind application of computing power. However, that is not to say that there is no role for direct numerical solution of dynamo systems. Such an approach should be regarded as the construction of silicon-based thought experiments, rather than the direct modelling of the astrophysical objects. Much can be learned from the investigation of such models; though, given the restrictions and the sensitivity of the dynamo system, extreme care must be taken in extrapolating conclusions to the relevant regimes for geophysics or astrophysics.