2 Mathematical formulation

We consider thermally driven convection in a three-dimensional, Cartesian layer of incompressible (Boussinesq) fluid rotating about the vertical. The fluid layer has depth $d$ , angular velocity $\unicode[STIX]{x1D6FA}$ , density $\unicode[STIX]{x1D70C}$ , kinematic viscosity $\unicode[STIX]{x1D708}$ , thermal diffusivity $\unicode[STIX]{x1D705}$ and magnetic diffusivity $\unicode[STIX]{x1D702}$ . Following standard practice, we adopt the layer depth $d$ , the thermal relaxation time $d^{2}/\unicode[STIX]{x1D705}$ and the temperature drop across the layer $\unicode[STIX]{x0394}T$ as the units of length, time and temperature respectively. We scale magnetic field intensities with $(2\unicode[STIX]{x1D6FA}\unicode[STIX]{x1D705}\unicode[STIX]{x1D707}_{0}\unicode[STIX]{x1D70C})^{1/2}$ , where $\unicode[STIX]{x1D707}_{0}$ is the magnetic permeability of the medium. With these units, and in standard notation, the magnetohydrodynamic (MHD) equations read

(2.1) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D70E}^{-1}(\unicode[STIX]{x2202}_{t}\boldsymbol{u}+\boldsymbol{u}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\boldsymbol{u})+E^{-1}\boldsymbol{e}_{z}\times \boldsymbol{u}=-\unicode[STIX]{x1D735}\tilde{p}+E^{-1}\boldsymbol{J}\times \boldsymbol{B}+Ra\,\unicode[STIX]{x1D703}\boldsymbol{e}_{z}+\unicode[STIX]{x1D6FB}^{2}\boldsymbol{u}, & \displaystyle\end{eqnarray}$$

(2.2) $$\begin{eqnarray}\displaystyle & \displaystyle (\unicode[STIX]{x2202}_{t}-q^{-1}\unicode[STIX]{x1D6FB}^{2})\boldsymbol{B}+\boldsymbol{u}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\boldsymbol{B}=\boldsymbol{B}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\boldsymbol{u}, & \displaystyle\end{eqnarray}$$

(2.3) $$\begin{eqnarray}\displaystyle & \displaystyle (\unicode[STIX]{x2202}_{t}-\unicode[STIX]{x1D6FB}^{2})\unicode[STIX]{x1D703}+\boldsymbol{u}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\unicode[STIX]{x1D703}=w, & \displaystyle\end{eqnarray}$$

(2.4) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{B}=\unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{u}=0, & \displaystyle\end{eqnarray}$$

where $\boldsymbol{J}=\unicode[STIX]{x1D735}\times \boldsymbol{B}$ is the current density, $\unicode[STIX]{x1D703}$ denotes the temperature fluctuations relative to a linear background profile and the velocity $\boldsymbol{u}=(u,v,w)$ . Four dimensionless numbers appear explicitly: the Rayleigh number $Ra$ , the Ekman number $E$ , the Prandtl number $\unicode[STIX]{x1D70E}$ and the Roberts number $q$ ; these are defined by

(2.5a-d ) $$\begin{eqnarray}Ra=\frac{g\unicode[STIX]{x1D6FC}\unicode[STIX]{x0394}Td^{3}}{\unicode[STIX]{x1D705}\unicode[STIX]{x1D708}},\quad E=\frac{\unicode[STIX]{x1D708}}{2\unicode[STIX]{x1D6FA}d^{2}},\quad \unicode[STIX]{x1D70E}=\frac{\unicode[STIX]{x1D708}}{\unicode[STIX]{x1D705}},\quad q=\frac{\unicode[STIX]{x1D705}}{\unicode[STIX]{x1D702}},\end{eqnarray}$$

where $g$ is the gravitational acceleration and $\unicode[STIX]{x1D6FC}$ is the coefficient of thermal expansion.

The infinite Prandtl number limit of the momentum equation (2.1) results in the equation

(2.6) $$\begin{eqnarray}\boldsymbol{e}_{z}\times \boldsymbol{u}=-\unicode[STIX]{x1D735}p+\boldsymbol{J}\times \boldsymbol{B}+R\unicode[STIX]{x1D703}\boldsymbol{e}_{z}+E\unicode[STIX]{x1D6FB}^{2}\boldsymbol{u},\end{eqnarray}$$

where $R$ is the rotational Rayleigh number defined by

(2.7) $$\begin{eqnarray}R=Ra\,E=\frac{g\unicode[STIX]{x1D6FC}\unicode[STIX]{x0394}Td}{2\unicode[STIX]{x1D6FA}\unicode[STIX]{x1D705}}.\end{eqnarray}$$

The new pressure $p$ is given by $p=E\tilde{p}$ ; its scaling though is immaterial since for Boussinesq convection the role of the pressure is simply to keep the flow solenoidal.

It is worth pointing out that the limit of infinite Prandtl number is a delicate one, and should not be thought of as being obtained either by letting $\unicode[STIX]{x1D708}\rightarrow \infty$ or $\unicode[STIX]{x1D705}\rightarrow 0$ with all the other parameters fixed. Indeed, $\unicode[STIX]{x1D708}$ and $\unicode[STIX]{x1D705}$ appear in $R$ , $E$ and $q$ , all of which must remain finite. The most glaring manifestation of taking this limit is the removal of the inertia terms from the momentum equation. This is not a mere simplification of the system: it has both mathematical and physical consequences. Mathematically, the momentum equation has become a linear diagnostic equation that defines the velocity in terms of the thermal and magnetic variables. Physically, the removal of the inertia terms has profound implications for the interpretation of the behaviour of the system. In the absence of buoyancy forces and dissipation, the original system (with inertia) supports two families of linear waves. For moderate rotation rates and strong background magnetic field, the two families correspond to rotationally modified Alfvén waves. On the other hand, if rotation dominates, then one family consists of modified inertial waves, whereas the other consists of hydromagnetic inertial waves; for these latter waves, both inertial and magnetic effects come into play, with the wave frequency given by $\unicode[STIX]{x1D714}^{2}\approx \unicode[STIX]{x1D714}_{A}^{4}/\unicode[STIX]{x1D714}_{I}^{2}$ , where $\unicode[STIX]{x1D714}_{A}$ and $\unicode[STIX]{x1D714}_{I}$ are the Alfvén and inertial frequencies respectively. In the new (inertialess) system, with momentum equation (2.6), the inertial waves acquire an infinite propagation speed. Only the hydromagnetic inertial waves remain; for these, both $\unicode[STIX]{x1D714}_{A}$ and $\unicode[STIX]{x1D714}_{I}$ tend to infinity, but such that $\unicode[STIX]{x1D714}_{A}^{4}/\unicode[STIX]{x1D714}_{I}^{2}$ is finite.

Thus the equations with infinite Prandtl number describe a rapidly rotating, strongly magnetised system. It is for this reason that we chose to scale the magnetic field with $(2\unicode[STIX]{x1D6FA}\unicode[STIX]{x1D705}\unicode[STIX]{x1D707}_{0}\unicode[STIX]{x1D70C})^{1/2}$ ; had we chosen the conventional scaling of magnetic field with $(\unicode[STIX]{x1D707}_{0}\unicode[STIX]{x1D70C})^{1/2}\unicode[STIX]{x1D705}/d$ then the field can be brought into play as $\unicode[STIX]{x1D70E}\rightarrow \infty$ only by a further scaling with $\unicode[STIX]{x1D70E}^{1/2}$ . Although there is no ambiguity in determining $\boldsymbol{u}$ and $\boldsymbol{B}$ from the inertialess equations, a difficulty arises in comparing $|\boldsymbol{u}|^{2}$ and $|\boldsymbol{B}|^{2}$ . If we denote the dimensional quantities by $\hat{\boldsymbol{u}}$ and $\hat{\boldsymbol{B}}$ then

(2.8) $$\begin{eqnarray}\frac{|\hat{\boldsymbol{B}}|^{2}/\unicode[STIX]{x1D707}_{0}\unicode[STIX]{x1D70C}}{|\hat{\boldsymbol{u}}|^{2}}=\frac{2\unicode[STIX]{x1D6FA}d^{2}}{\unicode[STIX]{x1D705}}\frac{|\boldsymbol{B}|^{2}}{|\boldsymbol{u}|^{2}}=\frac{\unicode[STIX]{x1D70E}}{E}\frac{|\boldsymbol{B}|^{2}}{|\boldsymbol{u}|^{2}}.\end{eqnarray}$$

Thus, in this system, with $\unicode[STIX]{x1D70E}\rightarrow \infty$ and $E$ finite, comparing kinetic and magnetic energies is meaningless.

This raises an interesting issue about how to interpret magnetic field intensities. Traditionally, one uses the equivalent Alfvén speed to compare the magnetic energy with the kinetic energy of the moving fluid; here, however, the Alfvén speed is formally infinite and so this is of no use. Instead, one way to proceed is to note that whereas equation (2.1) expresses the conservation of momentum, equation (2.6) expresses a force balance requirement. With this new system, momentum and energy are therefore not that well defined; indeed, it is now velocities, not accelerations, that are related to forces. Interpretation of the meaning of $\boldsymbol{B}$ must therefore be through comparison with the Coriolis force. That said, and as we shall see, it is still helpful to use the average values of $|\boldsymbol{u}|^{2}$ and $|\boldsymbol{B}|^{2}$ to measure the magnitudes of $\boldsymbol{u}$ and $\boldsymbol{B}$ relative to themselves, but not to each other.

In the horizontal directions we assume that all fields are periodic with periodicity $\unicode[STIX]{x1D706}$ – the aspect ratio. In the vertical we consider standard illustrative boundary conditions, namely that the boundaries are perfectly conducting, both thermally and electrically, impermeable and stress free. Formally these correspond to

(2.9) $$\begin{eqnarray}\unicode[STIX]{x1D703}=w=\unicode[STIX]{x2202}_{z}u=\unicode[STIX]{x2202}_{z}v=B_{z}=\unicode[STIX]{x2202}_{z}B_{x}=\unicode[STIX]{x2202}_{z}B_{y}=0\quad \text{at }z=0,1.\end{eqnarray}$$

We solve (2.2)–(2.4) numerically by standard pseudo-spectral methods. Details concerning the numerical methods can be found in Cattaneo, Emonet & Weiss (Reference Cattaneo, Emonet and Weiss2003). The numerical scheme requires at every time step a knowledge of the velocity $\boldsymbol{u}$ , which can be computed from (2.6). The details of how this can be performed are contained in appendix A. The implementation of the momentum equation was checked against the linear stability results of Roberts & Jones (Reference Roberts and Jones2000). The numerical resolution and parameter values for all the simulations presented in this paper are summarised in table 1.

Table 1. Summary of the parameter values and numerical resolution for the simulations. An asterisk in the $q$ column indicates that the dynamo failed, starting from a weak seed field.

3 Hydrodynamic convection

Before tackling the dynamo problem, it is useful to present some of the basic features of rotating non-magnetic convection at infinite Prandtl number that may play a role in the process of magnetic field generation.

Stability to infinitesimal perturbations that vary as

(3.1) $$\begin{eqnarray}\exp (st+\text{i}(\boldsymbol{k}_{H}\boldsymbol{\cdot }\boldsymbol{x}_{H}))\frac{\sin }{\cos }\unicode[STIX]{x03C0}z\end{eqnarray}$$

is governed by the dispersion relation

(3.2) $$\begin{eqnarray}s=\frac{REk_{H}^{2}k^{2}}{E^{2}k^{6}+\unicode[STIX]{x03C0}^{2}}-k^{2},\end{eqnarray}$$

where $k_{H}$ is the horizontal wavenumber and $k^{2}=k_{H}^{2}+\unicode[STIX]{x03C0}^{2}$ (cf. Roberts & Jones Reference Roberts and Jones2000). Although the standard dispersion relation for rotating Boussinesq convection is cubic in the growth rate $s$ (see, for example, Chandrasekhar Reference Chandrasekhar1961), here the inertial modes are filtered out, leading to the linear relation (3.2). It follows from (3.2) that convection sets in as a steady bifurcation when

(3.3) $$\begin{eqnarray}R=R_{0}=\frac{E^{2}k^{6}+\unicode[STIX]{x03C0}^{2}}{Ek_{H}^{2}}.\end{eqnarray}$$

Since the onset of convection is via a time-independent mode, it follows that the marginal stability criteria for the inertialess system must be identical to those of the system with inertia (which includes the $\unicode[STIX]{x2202}_{t}\boldsymbol{u}$ term and the nonlinear $\boldsymbol{u}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\boldsymbol{u}$ term; see Chandrasekhar (Reference Chandrasekhar1961)). In particular, for any Ekman number $E$ , there is a unique value of the horizontal wavenumber, $k_{H}=k_{H}^{crit}$ say, for which $R_{0}$ takes its minimum value, $R_{0}^{crit}$ . As $E\rightarrow 0$

(3.4a,b ) $$\begin{eqnarray}k_{H}^{crit}\sim \left(\frac{\unicode[STIX]{x03C0}^{2}}{2}\right)^{1/6}E^{-1/3},\quad R_{0}^{crit}\sim 3\left(\frac{\unicode[STIX]{x03C0}^{2}}{2}\right)^{2/3}E^{-1/3}.\end{eqnarray}$$

Figure 1 plots $k_{H}^{crit}$ and $R_{0}^{crit}$ as functions of $E$ , together with the asymptotic results (3.4). Clearly, for $E\lesssim 10^{-3}$ , the system is in the rapidly rotating asymptotic regime.

Figure 1. $k_{H}^{crit}$ and $R_{0}^{crit}$ as functions of $E^{-1}$ (solid lines), together with the asymptotic results (3.4) (dotted lines).

As $R$ is increased above $R_{0}^{crit}$ , the convection becomes nonlinear and time dependent, with associated changes in the planform. Figure 2 shows density plots of the temperature in vertical planes and horizontal planes near the surface, for three increasing values of $R$ at $E=10^{-4}$ , for which $R_{0}^{crit}\approx 187$ . Near onset, the convection assumes a columnar structure with horizontal scale controlled by $k_{H}^{crit}$ (figure 2 a,d, with $R=500$ ). Within each column, the fluid spins one way at the top and the other way at the bottom. As $R$ is increased, the convective pattern becomes asymmetric. Near the upper boundaries, the cold downflowing fluid forms a network surrounding the hot isolated upflows. Typically, the upflow centres develop a shallow cold return flow in which the fluid spins in the opposite direction to the hot fluid surrounding it (figure 2 b,e, with $R=10^{3}$ ). Because of the up–down symmetry of the Boussinesq equations, the reverse is true in the vicinity of the lower boundary. At yet higher values of $R$ , the convective cells are sheared out, leading to a network of thin downflows, winding itself round the upflows (figure 2 c,f, with $R=10^{4}$ ). Although the plots in figure 2 are for the specific value of $E=10^{-4}$ , we have found the behaviour described to be generic, provided that the Ekman number is sufficiently small.

For convection in rotating systems, we expect a correlation between flows and circulation. This manifests itself in a lack of reflectional symmetry, which can be measured by the kinetic helicity ${\mathcal{H}}=\langle \boldsymbol{u}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\times \boldsymbol{u}\rangle$ , where $\langle \cdot \rangle$ denotes an average over the fluid volume. In the Boussinesq approximation, the helicity density is antisymmetric with respect to the midplane of the layer, and so it is the distribution of helicity density in each half-layer that is of importance. Figure 3 shows the relative helicity, defined by

(3.5) $$\begin{eqnarray}h(z)=\frac{\langle \boldsymbol{u}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\times \boldsymbol{u}\rangle _{H}}{\langle \boldsymbol{u}^{2}\rangle _{H}^{1/2}\langle (\unicode[STIX]{x1D735}\times \boldsymbol{u})^{2}\rangle _{H}^{1/2}},\end{eqnarray}$$

where $\langle \cdot \rangle _{H}$ denotes an average over horizontal planes, for the three cases of figure 2. It can be seen that even for the case with the most vigorous convection, the influence of rotation is still strong.

Figure 2. Density plots of temperature fluctuations, in which light (dark) tones correspond to hot (cold) fluid: (a–c) show the $xy$ -plane near the upper surface ( $z=0.95$ ); (d–f) show the $xz$ -plane at the midpoint value of $y$ (i.e.  $y=\unicode[STIX]{x1D706}/2$ ). In all cases, $E=10^{-4}$ ; (a,d) $R=500$ , (b,e) $R=10^{3}$ , (c,f) $R=10^{4}$ . For (a,d,b,e), the aspect ratio $\unicode[STIX]{x1D706}=5$ ; for (c,f), $\unicode[STIX]{x1D706}=2$ .

Figure 3. Snapshots of relative helicity $h(z)$ versus $z$ for the three cases of figure 2: $E=10^{-4}$ ; $R=500$ (solid line); $R=10^{3}$ (dashed line); $R=10^{4}$ (dot-dashed line). Exact antisymmetry about the mid-plane ( $z=0.5$ ) would be recovered by time averaging.

In general, as $R$ increases at fixed $E$ , the rotational constraint is weakened and hence the characteristic horizontal scale of the convection increases. Nevertheless, since the inertial term is absent, the characteristic scale of the convection still decreases as $E^{-1/3}$ even in the nonlinear regime, as can be seen from simple considerations of force balance in the vorticity equation. This is illustrated by figure 4, which shows the horizontal power spectra of the velocity for three convective flows at different values of $E$ with $R$ chosen so that they all have the same degree of supercriticality. For fixed $E$ , the efficiency of the convection increases with $R$ . This is illustrated in figure 5, which shows the horizontally averaged temperature as a function of depth, for the same parameter values as figure 2. At the highest value of $R$ , the convection is so efficient that the layer is isothermal except in very thin boundary layers.

Figure 4. Horizontal power spectra of the kinetic energy for three cases of nonlinear convection with the same degree of supercriticality. Solid line: $E=10^{-3}$ , $R=500$ ; dashed line: $E=10^{-4}$ , $R=1077$ ; dot-dashed line: $E=10^{-5}$ , $R=2321$ . In all cases, the wavenumber refers to the case of $\unicode[STIX]{x1D706}=5$ . The spectra were computed over the interior regions of the simulations ( $0.1<z<0.9$ ) and time averaged. The location of the spectral peaks is consistent with a scaling in which the characteristic convective length scale decreases as $E^{1/3}$ .

Figure 5. Horizontally averaged temperature profiles for the three cases of figure 2: $E=10^{-4}$ ; $R=500$ (solid line); $R=10^{3}$ (dashed line); $R=10^{4}$ (dot-dashed line).

4 Dynamo action

In order to study the influence of the rotational constraint on convectively driven dynamos, we have investigated three different values of $E$ ( $10^{-3}$ , $10^{-4}$ , $10^{-5}$ ); for each of these we have considered a range of values of $R$ , from a few times supercritical up to vigorous convection. In each case, the convection is evolved from small random perturbations of the basic state, to an eventual stationary state; these purely hydrodynamic convective states are what we shall refer to as the kinematic flows. For each such flow, we introduce a small magnetic perturbation, at $t=0$ , which is then evolved for a range of values of $q$ , effectively increasing $Rm$ by increasing $q$ . Each simulation is evolved either until the field decays or until it grows and saturates, with the convection attaining a new, magnetised, stationary state. In table 1 we indicate which of the hydrodynamic flows are able to maintain dynamo action.

It is natural to ask under what conditions dynamo action can take place. It is not surprising to find that dynamo action is possible once the convection is sufficiently vigorous and the magnetic Reynolds number is large enough. A more challenging question is to ask what determines the saturation amplitude of the magnetic field and what are the mechanisms that govern the saturation process. In general, this is a very difficult question, not just to answer, but even to formulate precisely. Here, because of the absence of inertia in the momentum equation, we can definitely provide a precise formulation of the problem of dynamo saturation, and are also able, to some degree, to provide an answer. The crucial property of the formulation of the convection problem here is that the momentum equation is linear in the velocity. This allows us to decompose the velocity into the sum of two components, which we choose to have meaningful physical interpretations. Thus we write the total velocity $\boldsymbol{u}$ as

(4.1) $$\begin{eqnarray}\boldsymbol{u}=\boldsymbol{u}_{T}+\boldsymbol{u}_{M},\end{eqnarray}$$

where $\boldsymbol{u}_{T}$ and $\boldsymbol{u}_{M}$ satisfy, respectively, the equations

(4.2) $$\begin{eqnarray}\displaystyle & \displaystyle \boldsymbol{e}_{z}\times \boldsymbol{u}_{T}=-\unicode[STIX]{x1D735}p_{T}+R\unicode[STIX]{x1D703}\boldsymbol{e}_{z}+E\unicode[STIX]{x1D6FB}^{2}\boldsymbol{u}_{T}, & \displaystyle\end{eqnarray}$$

(4.3) $$\begin{eqnarray}\displaystyle & \displaystyle \boldsymbol{e}_{z}\times \boldsymbol{u}_{M}=-\unicode[STIX]{x1D735}p_{M}+\boldsymbol{J}\times \boldsymbol{B}+E\unicode[STIX]{x1D6FB}^{2}\boldsymbol{u}_{M}, & \displaystyle\end{eqnarray}$$

with $\unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{u}_{T}=\unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{u}_{M}=0$ and $p=p_{T}+p_{M}$ .

The component $\boldsymbol{u}_{T}$ is driven by Coriolis and buoyancy forces, whereas the component $\boldsymbol{u}_{M}$ is driven by Coriolis and magnetic forces. The velocity $\boldsymbol{u}_{T}$ exists even in the absence of magnetic fields, and defines the kinematic problem; $\boldsymbol{u}_{M}$ , on the other hand, exists only by virtue of the presence of the magnetic field. One can think of $\boldsymbol{u}_{M}$ as the means by which the system reacts to the presence of the magnetic field, leading eventually to the saturation of magnetic field growth. If the buoyancy force were to be prescribed, then $\boldsymbol{u}_{T}$ would always be independent of the magnetic field. As things are, in the saturated regime, the temperature distribution depends on both $\boldsymbol{u}_{T}$ and $\boldsymbol{u}_{M}$ , via (2.3), and hence $\boldsymbol{u}_{T}$ does not remain as the kinematic velocity.

It is not unreasonable to conceive of a situation in which the temperature distribution is not strongly contaminated by $\boldsymbol{u}_{M}$ . In this scenario, buoyancy drives the convection, the convection generates a magnetic field by dynamo action, and the magnetic field drives a velocity $\boldsymbol{u}_{M}$ , which slightly modifies $\boldsymbol{u}$ so as to saturate the dynamo growth. Here, $\boldsymbol{u}$ and $\boldsymbol{u}_{T}$ are similar, with the magnetic field having a weak effect on the velocity – we categorise such solutions as weak field solutions. That said, one can envisage a different scenario in which, when saturation occurs, a substantial fraction of the total velocity is driven magnetically. In this case, magnetic forces do not contribute just to the saturation of the dynamo, but to its actual driving. Here, the velocity in the saturated state is completely unrelated to that of the kinematic state; indeed, the kinematic state may not even act as a dynamo. In this case, the magnetic field has a strong effect on the velocity – we thus classify such solutions as strong field solutions. The justification for the consideration of strong field solutions is that at small $E$ the convection is strongly rotationally constrained. There are circumstances in which the presence of the magnetic field can alleviate this constraint, and possibly lead to more vigorous convection, which, in turn, can support stronger fields. This is well known in the related, but simpler, problem of rotating convection, in which the presence of an imposed magnetic field reduces considerably the critical Rayleigh number for the onset of convection (Chandrasekhar Reference Chandrasekhar1961; Eltayeb & Roberts Reference Eltayeb and Roberts1970; Eltayeb Reference Eltayeb1972). The dynamo problem is more complicated, in that the magnetic field is self-generated; nonetheless, it is not unreasonable to expect that a similar mechanism may be operating here. That being the case, one would expect strong field solutions to be favoured for small values of $E$ and moderate values of $R$ . As we shall see presently, these considerations are borne out by the results.

4.1 Weak field solutions

As an illustrative example of a weak field solution, we consider the case of $E=10^{-4}$ , $R=10^{4}$ , $q=1$ , whose convective pattern in the kinematic regime is shown in figure 2(c). Figure 6 shows the evolution of $\langle \boldsymbol{u}^{2}\rangle$ and $\langle \boldsymbol{B}^{2}\rangle$ , following the introduction of a small seed field. There follows an exponential growth in $\langle \boldsymbol{B}^{2}\rangle$ , accompanied by a decrease in $\langle \boldsymbol{u}^{2}\rangle$ , followed by a nonlinear saturation in which the system settles to a stationary state.

Figure 6. Time traces of $\langle \boldsymbol{u}^{2}\rangle$ and $\langle \boldsymbol{B}^{2}\rangle$ for the weak field dynamo with $E=10^{-4}$ , $R=10^{4}$ , $q=1$ . The amplitude of $\langle \boldsymbol{B}^{2}\rangle$ has been multiplied by $10^{5}$ in order that the two can be represented clearly on one plot.

Figure 7 shows snapshots of the vertical velocity near the top of the domain, for the kinematic and dynamic regimes. It can be seen that whereas the overall structure of the convective patterns remains roughly the same, the characteristic scale of the convection increases slightly in the dynamic regime. This increase in scale is a manifestation of the magnetic field relaxing the rotational constraint; as $R$ is increased, rendering the convection more supercritical, the increase in the convective scale in the saturated regime diminishes. This is confirmed by figure 8, which shows the velocity spectra for the two regimes; except for a slight shift towards smaller wavenumbers, the two spectra are quite similar.

Figure 7. Density plots of vertical velocity close to the upper boundary ( $z=0.95$ ), for ( $a$ ) the kinematic and ( $b$ ) the dynamic regime, for $E=10^{-4}$ , $R=10^{4}$ , $q=1$ ; light (dark) tones correspond to upward (downward) moving fluid. The plots are scaled independently.

Figure 8. Horizontal power spectra of the kinetic energy for the kinematic regime (solid line) and dynamic regime (dashed line) for the weak field dynamo: $E=10^{-4}$ , $R=10^{4}$ , $q=1$ , $\unicode[STIX]{x1D706}=2$ . The spectra were computed over the interior regions ( $0.1<z<0.9$ ) of the simulations and time averaged.

The magnetic field structure can be gauged from figure 9, which shows the vertical component of the electric current near the top of the layer, and the vertical component of the magnetic field in the middle of the layer. The current has a filamentary structure, with the filaments outlining the downflows of the convective patterns. The vertical magnetic field is concentrated in sheets, by and large, with a characteristic scale comparable with that of the convective flow. As can be seen from figure 6, the flow and field evolve on essentially the same time scale.

Figure 9. Density plots of (a) the vertical current near the upper boundary ( $z=0.95$ ) and (b) the vertical magnetic field at the midplane.

One of the key ideas of mean field electrodynamics is that the generation of large-scale (mean) magnetic fields, on scales much larger than that of the underlying flow, is associated with a lack of reflectional symmetry in the flow. Although, from inspection of figure 9, no large-scale structure in the field is immediately apparent, it is nonetheless natural to look more closely into this issue, given that the flow is undeniably helical. Given our geometry and boundary conditions, we can look for the presence of a mean field that is horizontal, but whose strength and orientation may vary with height $z$ (cf. Childress & Soward Reference Childress and Soward1972). Figure 10(a) shows the time average of $\langle \boldsymbol{B}_{H}\rangle _{H}$ , where $\boldsymbol{B}_{H}$ is the horizontal magnetic field and $\langle \cdot \rangle _{H}$ denotes an average over the horizontal plane. Figure 10(b) shows the quantity

(4.4) $$\begin{eqnarray}\unicode[STIX]{x1D6E4}(z)=\frac{|\langle \boldsymbol{B}_{H}\rangle _{H}|^{2}}{\langle |\boldsymbol{B}_{H}|^{2}\rangle _{H}},\end{eqnarray}$$

at several instances, approximately one turnover time apart, and as a time average. It can be seen that $\unicode[STIX]{x1D6E4}$ is a measure of the strength of the large-scale field relative to the total strength (Cattaneo & Hughes Reference Cattaneo and Hughes2006). Three things should be noted: $\unicode[STIX]{x1D6E4}$ is very small, it has no apparent organisation in $z$ and its temporal coherence is no longer than a turnover time. These points taken together suggest that $\boldsymbol{B}_{H}$ is nothing more than the horizontal projection of a random field; we find no evidence of any large-scale organisation of the magnetic field.

Figure 10. (a) Snapshot of the mean horizontal field as a function of depth for the weak field dynamo. (b) Snapshots of $\unicode[STIX]{x1D6E4}(z)$ versus $z$ at three times (the dashed, dotted and dot-dashed lines), together with the time average of $\unicode[STIX]{x1D6E4}(z)$ calculated over ten snapshots (solid line). The plot in (a) corresponds to the case shown as a dot-dashed line.

We can obtain some insight into the mechanism of field saturation by exploiting the velocity decomposition (4.1). Figure 11 shows the vertical components of $\boldsymbol{u}$ , $\boldsymbol{u}_{T}$ and $\boldsymbol{u}_{M}$ near the upper boundary, at the same time as for the plots in figure 9. The most striking feature is the remarkable similarity between $w$ and $w_{T}$ . Taken at face value, one might expect $\boldsymbol{u}_{M}$ to be miniscule. However, that is not entirely correct; for instance, at this particular time, $\langle \boldsymbol{u}_{M}^{2}\rangle /\langle \boldsymbol{u}_{T}^{2}\rangle \approx 0.3$ . It is though difficult to detect the features of $\boldsymbol{u}_{M}$ that are ultimately responsible for the dynamo saturation. For instance, the two marked boxes in figure 11 outline features in $\boldsymbol{u}_{M}$ . In the black box, the component of $w_{T}$ is stronger than that in the full velocity $w$ ; the effect of the Lorentz force is thus to slow down the vertical velocity. By contrast, in the green box, the opposite is true, with the total velocity being stronger then the thermally driven flow; the effect of the Lorentz force is thus to speed up the velocity. The conclusion is that the changes required in the velocity to saturate the dynamo are tremendously subtle.

Figure 11. Density plots, scaled consistently, of the vertical components of (a) the full velocity, (b) the thermal velocity and (c) the magnetic velocity. The time and location ( $z=0.95$ ) are identical to those of the plot of the vertical component of the current in figure 9.

A better approach comes by noting that from the induction equation, magnetic field is generated by field line stretching (the $\boldsymbol{B}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\boldsymbol{u}$ term); hence it is velocity gradients, not the velocity per se, that are of interest. This motivates looking at the vorticity field. Figure 12 shows the vertical vorticity corresponding to $\boldsymbol{u}$ , $\boldsymbol{u}_{T}$ and $\boldsymbol{u}_{M}$ at the same time and location as in figure 9. Unlike their velocity counterparts, the difference between $\unicode[STIX]{x1D74E}$ and $\unicode[STIX]{x1D74E}_{T}$ is now more striking; indeed, $\unicode[STIX]{x1D74E}_{T}$ looks like a blurred version of $\unicode[STIX]{x1D74E}$ . The difference $\unicode[STIX]{x1D74E}_{m}$ therefore consists mainly of these missing sharp vorticity features. The mechanism of saturation therefore clearly requires highly localised changes in the velocity gradients. Unlike for the magnetically driven velocity, changes in the magnetically driven vorticity are reflected in the overall vorticity. Nevertheless, the action of the magnetic torques is not simply to remove or decrease regions of high vorticity, but, somewhat counterintuitively, instead acts to sharpen the vorticity distribution of the thermally driven component. It is easy to verify that the sharpening of the vorticity field is a magnetic phenomenon by noting that the vorticity distribution for the kinematic velocity (figure 12 d) exhibits the same degree of smoothness as the thermally driven velocity in the saturated regime.

Figure 12. Density plots of the vertical components of (a) the full vorticity, (b) the thermal vorticity and (c) the magnetic vorticity, corresponding to the velocities shown in figure 11. For comparison, a representative plot of the vertical component of the vorticity in the kinematic regime is shown in (d). Plots (a), (b) and (c) are scaled consistently; plot (d) is scaled independently.

4.2 Strong field solutions

As an illustrative example of a strong field solution, we first consider the case of $E=10^{-4}$ , $R=500$ , $q=20$ , whose convective pattern in the kinematic regime is shown in figure 2(a). Figure 13 shows the evolution of $\langle \boldsymbol{u}^{2}\rangle$ and $\langle \boldsymbol{B}^{2}\rangle$ , following the introduction of a weak seed magnetic field. As in the weak field case, the magnetic energy grows exponentially and then saturates. However, the effect of this field growth on the kinetic energy is radically different. Following the saturation of the field, the time scale of variation of the kinetic energy changes dramatically; it is clear that a much shorter time scale has been introduced. We note also that the growth and saturation of the magnetic field actually causes an abrupt increase in the kinetic energy. This second feature is found in most, but not all, of our simulations of strong field solutions; the marked change in temporal variations is, however, common to all such solutions.

Figure 13. Time traces of $\langle \boldsymbol{u}^{2}\rangle$ and $\langle \boldsymbol{B}^{2}\rangle$ for the strong field dynamo with $E=10^{-4}$ , $R=500$ , $q=20$ . The amplitude of $\langle \boldsymbol{B}^{2}\rangle$ has been multiplied by $5\times 10^{3}$ in order that the two can be represented clearly on one plot.

Figure 14 shows snapshots of the vertical velocity near the top of the domain, for the kinematic and dynamic regimes. Whereas for the weak field solutions, except for a slight change in scale, the nature of the convection is the same, here not only the scale, but the very organisation of the convection has been changed. These changes are made apparent by inspection of the velocity spectra in the two regimes, as shown in figure 15. In the kinematic regime, there is a sharp peak at the convective scale, with steep fall-offs to both larger and smaller scales. By contrast, in the dynamical regime, there is a flattening in the spectrum, with more energy at both smaller and larger scales than in the original convective flow.

Figure 14. Density plots of vertical velocity close to the upper boundary ( $z=0.95$ ), for (a) the kinematic and (b) the dynamic regime, for $E=10^{-4}$ , $R=500$ , $q=20$ ; light (dark) tones correspond to upward (downward) moving fluid. The plots are scaled independently.

Figure 15. Horizontal power spectra of the kinetic energy for the kinematic regime (solid line) and dynamic regime (dashed line) for the strong field dynamo: $E=10^{-4}$ , $R=500$ , $q=20$ , $\unicode[STIX]{x1D706}=5$ . The spectra were computed over the interior regions ( $0.1<z<0.9$ ) of the simulations and time averaged.

The structure of the magnetic field in the middle of the layer and that of the vertical current near the top of the domain are shown in figure 16. The thing to note is its remarkable similarity to the corresponding image for the weak field, shown in figure 9. Even though one is a strong field solution and the other a weak field solution, from a magnetic point of view they look the same. We shall return to this point in the following subsection. We note here that in order to make the comparison easier, we are showing only a part of the horizontal plane for the strong field solution with the same aspect ratio ( $\unicode[STIX]{x1D706}=2$ ) as the weak field solution of figure 9.

Figure 16. Density plots of (a) the vertical current near the upper boundary ( $z=0.95$ ) and (b) the vertical magnetic field at the midplane for the strong field case of $E=10^{-4}$ , $R=500$ , $q=20$ . Note that in order to make a comparison with figure 9, only a $2\times 2$ section of the domain is shown; the simulation has aspect ratio $\unicode[STIX]{x1D706}=5$ .

Figure 17 shows the velocity decomposition (4.1), analogous to figure 11 for the weak field solution. As seen in the velocity spectra, this figure gives a visual confirmation that in the strong field solutions the kinematic and saturated velocities are dramatically different. Any prominent feature in $\boldsymbol{u}$ can be seen to originate either from $\boldsymbol{u}_{T}$ or $\boldsymbol{u}_{M}$ . Therefore, whereas in the weak field case, $\boldsymbol{u}$ and $\boldsymbol{u}_{T}$ are nearly equal, here $\boldsymbol{u}_{T}$ and $\boldsymbol{u}_{M}$ both contribute significantly to the full velocity. In fact, for smaller values of $E$ , it is possible to find strong field solutions in which the contribution to the full velocity from $\boldsymbol{u}_{M}$ becomes the dominant one. This is illustrated in figure 18, which shows the velocity decomposition for a strong field solution with $E=10^{-5}$ , $R=10^{3}$ , $q=10$ . It is striking that here $\boldsymbol{u}$ and $\boldsymbol{u}_{M}$ are virtually identical, the very opposite of the weak field case. The reason for this is that strong field solutions require both a strong rotational constraint and also sufficiently vigorous convection in order to generate a magnetic field that can break that constraint. As $E$ becomes smaller, we anticipate therefore that the magnetically-driven component of the velocity will become more and more dominant.

Figure 17. Density plots, scaled consistently, of the vertical components of (a) the full velocity, (b) the thermal velocity and (c) the magnetic velocity for the strong field case of $E=10^{-4}$ , $R=500$ , $q=20$ . The time and location ( $z=0.95$ ) are identical to those of the plot of the vertical component of the current in figure 16.

Figure 18. Density plots at $z=0.95$ , scaled consistently, of the vertical components of ( $a$ ) the full velocity, ( $b$ ) the thermal velocity and ( $c$ ) the magnetic velocity for the ‘super strong’ field case of $E=10^{-5}$ , $R=10^{3}$ , $q=10$ .

4.3 Hysteresis and solution branches

In the previous subsections we have identified two types of solutions: weak and strong. Two issues naturally arise. One is the relation between these two types of solution as parameters are varied. The other, given the fact that the strong field solutions are clearly of finite amplitude, is whether these solutions exhibit hysteretical behaviour. We shall address these points by considering how the nature of the solutions changes for a fixed value of $E$ ( $E=10^{-4}$ ) as $R$ is varied; physically this corresponds to fixing the rotation rate and varying the intensity of the convective driving.

As a starting point, we know that at $R=500$ there is a strong field solution. By changing $R$ in small steps, and keeping $q$ fixed at $q=10$ , we are able to examine how the strong field solution changes as $R$ is both increased and decreased. The results are summarised in figure 19, which plots the root mean square value of the magnetic field strength versus $R$ ; each asterisk corresponds to an actual simulation that used a neighbouring simulation as initial condition. Solutions of finite amplitude can be traced back until $R\approx 330$ ; for smaller values of $R$ the field decays to zero and no dynamo solutions can be found. Furthermore, we have investigated the onset of dynamo action from infinitesimal magnetic field perturbations; for this value of $q$ , we find instability for $R\gtrsim 420$ . Taken together, these results suggest a bifurcation diagram of the form sketched in figure 20(a). For $R<R_{0}$ there is no dynamo solution; for $R_{0}<R<R_{3}$ , dynamo solutions can be found only if the amplitude of the initial perturbation is sufficiently large; and for $R>R_{3}$ all initial conditions are unstable to dynamo action. We anticipate that the trivial solution is connected to the stable finite amplitude branch via a subcritical branch of unstable solutions. This picture should be contrasted with another possible bifurcation diagram often discussed in the dynamo literature (see Roberts Reference Roberts, Roberts and Soward1978; Fautrelle & Childress Reference Fautrelle and Childress1982; Roberts & Soward Reference Roberts and Soward1992; St-Pierre Reference St-Pierre, Proctor and Gilbert1994), in which there is a supercritical branch of stable solutions that intersects the unstable subcritical branch, as sketched in figure 20(b). In this case, there is a range of $R$ ( $R_{1}<R<R_{2}$ ) in which there exist two non-trivial stable solutions. This leads naturally to the classification of solutions in terms of weak and strong field branches. In our system, interpreted as in figure 20(a), there is only one non-trivial stable branch. The distinction we make between weak and strong field solutions is in terms of how the saturated velocity relates to the kinematic velocity. Recently, the scenario described in figure 20(b) has been verified by the numerical investigation by Dormy (Reference Dormy2016) of dynamo action in a rapidly rotating spherical shell with finite inertia. This leads us to conjecture that the absence of the stable supercritical branch in our system is related to the absence of inertia. That being the case, it is reasonable to expect that the extent of the supercritical weak branch, or, equivalently, the distance between $R_{1}$ and $R_{2}$ , decreases with decreasing $E$ .

Figure 19. Time-averaged values of $\langle B^{2}\rangle ^{1/2}$ versus $R$ , with $E=10^{-4}$ , $q=10$ . The dynamo fails for $R\lesssim 330$ ; the zero field solution is unstable for $R\gtrsim 420$ .

Figure 20. Sketches of the amplitude of the magnetic field versus $R$ ; stable (unstable) solutions are denoted by solid (dashed) lines. In (a), infinitesimal magnetic fields are stable for $R<R_{3}$ . All stable dynamo solutions lie on the strong field branch, at finite amplitude, and these can persist down to $R=R_{0}$ . In (b), a strong field branch again exist for $R>R_{0}$ . Now though, infinitesimal magnetic fields are stable only for $R<R_{1}$ ; for $R_{1}<R<R_{2}$ there is also a stable weak field branch.

Given that in our system there appears to be only one branch of non-trivial solutions up to moderate values of $R$ , it is natural to ask what is the relation between our strong and weak field solutions. In particular, is there a gradual transformation from one to the other as $R$ is increased? The evidence suggests that this is indeed the case. As $R$ increases, the velocity changes gradually between the strong field and weak field solutions. However, as noted already, the magnetic field remains essentially unchanged. Thus, in some sense, the system produces only one type of dynamo, namely that in which the influence of rotation is minimised. At high values of $R$ , the influence of rotation is small, and so this is achieved by a slight increase in spatial scale in the velocity and where $\boldsymbol{u}\approx \boldsymbol{u}_{T}$ . By contrast, at smaller values of $R$ , the influence of rotation is strong; dynamo action succeeds by a dramatic increase in spatial scale in the velocity, driven almost totally by the magnetic component $\boldsymbol{u}_{M}$ .

In order to gain more physical insight into the nature of this system it is instructive to consider the simpler problem of rotating convection in the presence of an imposed magnetic field. It is well known that the effect of rotation is always to stabilise the system, in the sense that the critical value of the Rayleigh number for the onset of convection increases with increasing rotation. The same is true for a magnetic field, in the absence of rotation. However, as noted originally by Chandrasekhar (Reference Chandrasekhar1961), in the presence of rotation a magnetic field can actually be destabilising, in the sense that the critical value of $Ra$ decreases with increasing field strength, up to some point. Similar results apply to rotating convection in the nonlinear regime; for a fixed value of $Ra$ the convective efficiency can increase with imposed field strength, again, up to some point (Stellmach & Hansen Reference Stellmach and Hansen2004). These considerations motivate us to conceive of a gedankenexperiment in which, for fixed rotation, and for each Rayleigh number, the strength of the imposed magnetic field is tuned such that the convection is as efficient as possible, as measured by the Nusselt number, say. If one were then to plot these maximised Nusselt numbers as a function of $Ra$ , we would envisage a plot reminiscent of figure 19, with a monotonic increase of Nusselt number with $Ra$ . There will be a value of $Ra$ below which convection is impossible, no matter what the strength of the magnetic field; above that, there will be a range of $Ra$ for which convection will be possible only if a magnetic field is present; and above that, there will be a range of $Ra$ for which convection would have occurred anyway, but for which the efficiency is enhanced by the magnetic field. Eventually, of course, at very high $Ra$ , the magnetic field will cease to have a beneficial effect. The important thing to note is that at moderate values of $Ra$ , the effect of the optimal magnetic field is always to make the system more supercritical. In our dynamo model, the effect is the same except that the magnetic field is not externally imposed, but is self-determined. Thus the dynamo seems to operate by removing as much of the rotational constraint as possible; at high $R$ there is little constraint to remove, whereas at low values of $R$ the constraint to be removed is significant.

The computations leading to figure 19 are all performed at the fixed value of $q=10$ . However, the weak and strong field solutions exhibited in §§ 4.1 and 4.2 have different values of $q$ ( $q=1$ and $q=20$ respectively). These choices were motivated chiefly by numerical considerations, since dynamos driven by very vigorous convection and with high values of $q$ becomes prohibitively expensive. Hence it is natural to decrease the value of $q$ as $R$ increases. However, such variations will not influence the fundamental nature of the solutions; for example, if we had run the case of § 4.1, but with $q=20$ , the solution would still be a weak field solution, but with magnetic boundary layers that would be extremely expensive to compute.

5 Force balance considerations

Having identified two different mechanisms of saturation, leading to weak and strong solutions, it is important to determine the nature of the force balance in each case. There are four dynamical ingredients: buoyancy, the Coriolis force, viscous stresses and magnetic forces. Following common practice, we shall refer to these by A (for Archimedean), C, V and M. For purely hydrodynamic convection, the balance is between buoyancy, Coriolis and viscous forces (VAC); simple considerations of vorticity production lead to the estimate that the critical scale of convection $\ell$ scales as $\ell \sim E^{-1/3}$ . When the dynamo operates, there are more possibilities: the key question is whether the resulting force balance remains as VAC, or whether the magnetic field is important, in which case the balance can be either MAC or VMAC. In general, it is not straightforward to classify the solution in terms of these categories. However, here, because we are dealing with an inertialess system, we can make progress by again exploiting the decomposition (4.1).

By definition, our weak field solutions are ones for which the kinematic and dynamic velocities are similar. As such, since the hydrodynamic balance is VAC, then so must be the balance in the MHD regime. We have already argued in § 4 that the magnetic field enters only in a very subtle way in order to saturate the dynamo growth. The more complicated issue is to determine what happens in the case of the strong field solutions. Clearly magnetic forces enter the picture in a significant manner. We have already noted that for small $E$ the magnetic velocity $\boldsymbol{u}_{M}$ is almost identical to the full velocity, and the latter has a characteristic horizontal scale that is much larger than that in the corresponding hydrodynamic regime. Inspection of (4.3) suggests that to a first approximation, magnetic forces almost balance the Coriolis force, so that the rotational constraints are greatly reduced. However, it is rather difficult to have convection without some buoyancy and so (4.2) must come into play. Indeed it does, but with larger effective values of $R$ and $E$ for the given rotation rate, leading to solutions similar to the weak field cases obtained at much lower rotation rates. One could therefore argue that the strong field solutions could be classified as MC–VAC; in other words, the strong field solutions organise themselves to resemble weak field solutions by doing whatever it takes to reduce the influence of the Coriolis force.

We have made progress through the introduction of the velocities $\boldsymbol{u}_{T}$ and $\boldsymbol{u}_{M}$ , necessitating the solution of an auxiliary equation, either (4.2) or (4.3). More typically, the approach is to solve the governing MHD equations (2.1)–(2.4), and, for those cases in which there is reason to believe that inertia is negligible, to address the issue of force balance by the application of what is known as Taylor’s constraint (Taylor Reference Taylor1963). If inertia can be neglected, then the momentum equation reduces to (2.6). Taking the $z$ -component of the curl of (2.6) then gives

(5.1) $$\begin{eqnarray}\frac{\unicode[STIX]{x2202}w}{\unicode[STIX]{x2202}z}=(\unicode[STIX]{x1D735}\times (\boldsymbol{J}\times \boldsymbol{B}))\boldsymbol{\cdot }\boldsymbol{e}_{z}+E(\unicode[STIX]{x1D6FB}^{2}\unicode[STIX]{x1D74E})\boldsymbol{\cdot }\boldsymbol{e}_{z},\end{eqnarray}$$

where the vorticity $\unicode[STIX]{x1D74E}=\unicode[STIX]{x1D735}\times \boldsymbol{u}$ . On integrating from $z=0$ to $z=1$ , the Coriolis term (i.e. the left-hand side) vanishes as a result of the impermeable boundary; this leads to the exact result

(5.2) $$\begin{eqnarray}\int _{0}^{1}(\unicode[STIX]{x1D735}\times (\boldsymbol{J}\times \boldsymbol{B}))\boldsymbol{\cdot }\boldsymbol{e}_{z}\,\text{d}z=-E\int _{0}^{1}(\unicode[STIX]{x1D6FB}^{2}\unicode[STIX]{x1D74E})\boldsymbol{\cdot }\boldsymbol{e}_{z}\,\text{d}z.\end{eqnarray}$$

In the limit of vanishing $E$ , assuming that the vorticity remains sufficiently regular, the integral on the right-hand side vanishes. This gives Taylor’s constraint in the form

(5.3) $$\begin{eqnarray}\int _{0}^{1}(\unicode[STIX]{x1D735}\times (\boldsymbol{J}\times \boldsymbol{B}))\boldsymbol{\cdot }\boldsymbol{e}_{z}\,\text{d}z=0,\end{eqnarray}$$

which has to be satisfied pointwise for each value of $x$ and $y$ (Jones & Roberts Reference Jones and Roberts2000). A dynamo system that generates a solution with strong magnetic field, no inertia and negligible viscous stresses (i.e. MAC balance) must satisfy the constraint (5.3). Thus it can be used to test whether a solution with these properties has been found. It is important to realise that the vanishing of the integral is due to cancellation and not simply to vanishing of the integrand. This is typically addressed by normalising the integral through consideration of the ‘Taylorisation’ function

(5.4) $$\begin{eqnarray}{\mathcal{T}}(x,y)=\frac{\displaystyle \int _{0}^{1}(\unicode[STIX]{x1D735}\times (\boldsymbol{J}\times \boldsymbol{B}))\boldsymbol{\cdot }\boldsymbol{e}_{z}\,\text{d}z}{\displaystyle \int _{0}^{1}|(\unicode[STIX]{x1D735}\times (\boldsymbol{J}\times \boldsymbol{B}))\boldsymbol{\cdot }\boldsymbol{e}_{z}|\,\text{d}z}.\end{eqnarray}$$

Note, given the exact result (5.2), one might consider the alternative normalisation defined by

(5.5) $$\begin{eqnarray}\widehat{{\mathcal{T}}}=\frac{\displaystyle \int _{0}^{1}\unicode[STIX]{x1D6FB}^{2}\unicode[STIX]{x1D714}_{z}\,\text{d}z}{\displaystyle \int _{0}^{1}|\unicode[STIX]{x1D6FB}^{2}\unicode[STIX]{x1D714}_{z}|\,\text{d}z}.\end{eqnarray}$$

If $|{\mathcal{T}}(x,y)|\ll 1$ for all values of $x$ and $y$ then the solution should have the properties listed above. In our case, we know that the solutions are not affected by inertia and that the magnetic field plays a significant role in our strong field solutions. We could thus use (5.4) as a measure of the importance of viscous stresses in determining the solution. Figure 21 shows the results of calculating ${\mathcal{T}}(x,y)$ , at two different times, for the strong field solution with $E=10^{-5}$ , $R=10^{3}$ , $q=10$ , for which $\boldsymbol{u}$ , $\boldsymbol{u}_{T}$ and $\boldsymbol{u}_{M}$ are shown in figure 18. Clearly ${\mathcal{T}}(x,y)$ fluctuates very rapidly on small scales and attains $O(1)$ values over much of the domain. What therefore should one make of this?

Figure 21. Density plots of the Taylorisation function ${\mathcal{T}}(x,y)$ for the case of $E=10^{-5}$ , $R=10^{3}$ , $q=10$ , $\unicode[STIX]{x1D706}=2$ . The panels correspond to two different times.

Taken at face value, this result suggests that our strong field solution is not in MAC balance, and hence that viscous stresses must play a role; a similar conclusion was reached by Jones & Roberts (Reference Jones and Roberts2000). Indeed, this is not inconsistent with our idea that the strong field solution is MC–VAC. Having said that, it is perhaps useful to introduce a cautionary note. The calculation of ${\mathcal{T}}$ or $\widehat{{\mathcal{T}}}$ is tricky and prone to the introduction of high frequency errors. The evaluation of ${\mathcal{T}}$ involves a product and two derivatives; the former introduces aliasing errors and the latter amplifies numerical errors. Alternatively, the evaluation of $\widehat{{\mathcal{T}}}$ does not involve a product, but does introduce a further derivative of the primitive variable $\boldsymbol{u}$ . To give a measure of the problems involved, we note that the three derivatives acting on the highest frequencies give rise to an amplification factor of $O(10^{7})$ for the run shown in figure 18. Bearing this in mind, it natural to ask if the fine structure and rapid oscillations present in figure 21 might not simply be a numerical artefact. Partly to address this issue we have increased the resolution of the solution to $1024\times 1024\times 257$ prior to the calculation of ${\mathcal{T}}(x,y)$ . This was achieved by interpolating the solution to a finer grid and integrating the evolution equations for $10^{6}$ steps at the higher resolution. The results have remained largely the same, which leads us to conclude that the lack of Taylorisation as evidenced in figure 21 is physical and not a numerical artefact. However, for simulations that are only marginally resolved, the possibility of numerical contamination should be taken into consideration before drawing any conclusion.

6 Concluding remarks

We have considered dynamo action in thermally driven convection under the influence of rotation but with no inertia. Our primary motivation for considering an inertialess system is that it allows a useful decomposition of the velocity into thermal and magnetic components, which can be used to gain valuable insights into the processes that lead to dynamo saturation. We have identified two types of dynamo solution, weak and strong, based on how close the saturated velocity remains to the kinematic one. For the weak field solutions, the saturated velocity remains close to the kinematic velocity, with just a slight increase in scale and a slight decrease in convective efficiency. The velocity $\boldsymbol{u}$ in the saturated regime bears a remarkable similarity to the thermal velocity $\boldsymbol{u}_{T}$ ; the role of $\boldsymbol{u}_{M}$ is to introduce very subtle changes, sufficient to halt the dynamo growth. By contrast, for the strong field solutions, the saturated velocity bears no resemblance to the kinematic flow. There is a marked increase in scale and the convective efficiency increases substantially (the heat flux increases by $50\,\%$ for the case shown in figure 14). As $E$ decreases, the velocity $\boldsymbol{u}$ becomes closer to the magnetic component $\boldsymbol{u}_{M}$ .

Remarkably, the magnetic field for the weak and strong field solutions look very similar. In other words, although the method by which dynamo action is achieved is very different, the final product is essentially the same. This has led us to speculate that weak and strong field solutions lie on the same dynamo branch. This is a finite amplitude branch of solutions, presumably connected to the non-dynamo solution by an unstable, subcritical branch, as sketched in figure 20(a). We have found no evidence for a finite amplitude, supercritical branch of solutions. This should be contrasted with the traditional picture in which both branches coexist, as shown in figure 20(b). We conjecture therefore that the nature of the bifurcation diagram is related to the importance of inertia. In the absence of rotation, for which inertia is important, we believe that only the lower branch exists; by contrast, in our system, when rotation is dominant and inertia is negligible, there is no lower branch of solutions. It is only in cases in which rotation is important but inertia still plays a role that the upper and lower branches can coexist. This would go some way towards explaining why it has been difficult to find the upper branch of solutions by simulating the full system and reducing the effect of inertia simply by reducing $E$ (see, for example, Soderlund, King & Aurnou Reference Soderlund, King and Aurnou2012; King & Buffett Reference King and Buffett2013). The problem is that, in the full system, $E$ probably has to be extremely small before the upper branch manifests itself. This can, to some extent, be circumvented by considering the case of high Prandtl number; in the current study, $Pr$ is formally infinite, inertia drops out of the equation, and the upper branch comes shining through. It is however of interest to note that in the plane layer dynamo simulations of Stellmach & Hansen (Reference Stellmach and Hansen2004), a transition to large-scale convection occurs at low Ekman number. Indeed, the horizontal convection scale increases until it is comparable with the size of the domain, which for computational reasons was constrained to be fairly small. Possibly because of this restriction, the large-scale convective flow is not capable of maintaining the strong field indefinitely; the field decays and the flow reverts to being small scale. Field amplification can then restart, leading eventually to another transition to large-scale convection. Thus the amplitude of the generated magnetic field fluctuates strongly in time.

It is important to make the distinction between weak and strong field branches and what we refer to as weak and strong field solutions. It is slightly ironic that our weak and strong field solutions lie on the one branch that is traditionally called the strong field branch. Our nomenclature addresses the issue of the force balance that leads to dynamo saturation rather than the branch on which the solutions lie. As such, we found that our weak solutions may be characterised as VAC, and our strong field solution as MC–VAC, in which the principal role of the magnetic field is to remove rotational constraints.

The fact that both our weak and strong solutions have viscous stresses playing a role raises the interesting question of what happens to these inertialess solutions as $E$ is decreased to yet smaller values. This is not only of intrinsic mathematical interest, but is pertinent to the application of some of these ideas to planetary dynamos, which are believed to operate at very small values of $E$ (see, for example, Dormy Reference Dormy2016; Hughes & Cattaneo Reference Hughes and Cattaneo2016). This might be usefully addressed through consideration of how the bifurcation diagram in figure 20(a) might change as $E$ is decreased. We recall that the solutions in $R_{0}<R<R_{3}$ are definitely strong, whereas those with $R\gg R_{3}$ are weak. As $E$ is decreased, the value of $R_{3}$ , denoting the onset of dynamo action from infinitesimal perturbations, will increase. This implies that by decreasing $E$ , a weak field solution would be transformed into a strong field solution. A more fundamental issue is to address what happens to $R_{0}$ as $E$ is decreased. Here there are two physically interesting possibilities. One is that the location of $R_{0}$ remains unchanged. This would imply that the range of values of $R$ for which strong field solutions exist continues to increase, and that any finite amplitude strong field solution cannot be stabilised by increasing the rotation. Presumably, this would be the case if the solutions were in MAC balance; the value of $E$ is irrelevant and the solution continues to operate for arbitrarily small values of $E$ . The other possibility is that $R_{0}$ increases as $E$ decreases, but at a rate no faster than that of $R_{3}$ . Then the domain of existence of strong field solutions increases, but any strong field solution will ultimately be stabilised at sufficiently high rotation. It would be of tremendous interest to distinguish between these possibilities. However, we admit that this would be computationally extremely challenging. To give an idea of the scale of the problem, decreasing $E$ from $10^{-4}$ to $10^{-5}$ in our simulations led to over an order of magnitude increase in computational requirements.