1. Introduction
The shallow-water equations are widely used as an idealised model of stratified fluid dynamics in a thin layer, as generically occurs in planetary atmospheres and oceans (e.g. Zeitlin Reference Zeitlin2018). In their simplest incarnation with no bottom topography, the equations describe the motion of an inviscid fluid of constant density occupying
$0\lt z\lt h(\boldsymbol {x},t)$
, beneath an overlying quiescent fluid of negligible density; here
$\boldsymbol {x}$
is the horizontal position, and
$z$
is an upwards vertical coordinate. When the fluid depth
$h(\boldsymbol {x},t)$
is much smaller than the horizontal length scale of the flow, the hydrostatic approximation can be made, and solutions exist with the horizontal flow
$\boldsymbol {u}$
independent of
$z$
(e.g. Gill Reference Gill1982). This leads to the coupled equations
\begin{align} \partial _t \boldsymbol {u} + \boldsymbol {u} \cdot {\boldsymbol {\nabla }} \boldsymbol {u} = - g {\boldsymbol {\nabla }} h + \boldsymbol {F}, \end{align}
\begin{align} \partial _t h + {\boldsymbol {\nabla }} \cdot (h {\boldsymbol {u}}) = 0, \end{align}
where
$g$
is the acceleration due to gravity, and
$\boldsymbol {F}$
is any
$z$
-independent forcing or dissipation. There is an obvious extension including background rotation, as in the equations originally derived by Laplace (Reference Laplace1776). Although the shallow-water equations have direct applications to barotropic flow in the ocean, they are often used with a reduced gravity
$g'$
to model upper oceanic flows above a deep quiescent layer of larger density, or as a quasi-two-dimensional
$({\boldsymbol {x}},t)$
idealisation of three-dimensional
$({\boldsymbol {x}},z,t)$
baroclinic dynamics in a continuously stratified flow, perhaps using the idea of equivalent depth (e.g. Gill Reference Gill1982; Zeitlin Reference Zeitlin2018).
For numerical solutions of the shallow-water equations in a strongly nonlinear regime, a scale-selective dissipation term is usually included in
$\boldsymbol {F}$
. An obvious choice is to set
$\boldsymbol {F}=\nu \nabla ^2 \boldsymbol {u}$
in (1.1), where
$\nabla ^2$
is the horizontal Laplacian operator. But this choice is undesirable: it does not lead to negative definite energy dissipation, and it violates angular momentum conservation (Gent Reference Gent1993; Schär & Smith Reference Schär and Smith1993; Shchepetkin & O’Brien Reference Shchepetkin and O’Brien1996; Ochoa et al. Reference Ochoa, Sheinbaum and Jiménez2011). Two approaches have been used to generate alternative forms of the dissipation that are consistent with the fundamental physical principle that it be the divergence of a symmetric tensor (Batchelor Reference Batchelor1967). In the first approach, Shchepetkin & O’Brien (Reference Shchepetkin and O’Brien1996) and Gilbert et al. (Reference Gilbert, Riedinger and Thuburn2014) set
\begin{equation} F_i = \frac {1}{h} \frac {\partial {}}{\partial {x_j}} (h \sigma _{ij}), \; \; \; \sigma _{ij} = \nu \left ( \frac {\partial {u_i}}{\partial {x_j}} + \frac {\partial {u_j}}{\partial {x_i}} - \varsigma \delta _{ij} \frac {\partial {u_k}}{\partial {x_k}} \right ), \end{equation}
for some parameter
$\varsigma$
, building on the study of Schär & Smith (Reference Schär and Smith1993) with
$\varsigma =1$
. The factors of
$h$
and the symmetric form of
$\sigma _{ij}$
ensure conservation of angular momentum, and Gilbert et al. (Reference Gilbert, Riedinger and Thuburn2014) proved negative semidefinite energy dissipation provided
$\varsigma \leqslant 1$
. However, this approach does not uniquely determine a value of
$\varsigma$
. The second approach is to develop an asymptotic reduction of the full three-dimensional Navier–Stokes equations as
$\varepsilon \to 0$
(Marche Reference Marche2007), where
$\varepsilon$
is the aspect ratio of the flow (i.e. a characteristic depth divided by a horizontal length scale). The leading-order momentum balance is then
$\partial ^2 {\boldsymbol {u}} / \partial {z}^2 =0$
; however, applying zero tangential stress at the top and bottom of the fluid layer, the leading-order flow
$\boldsymbol {u}$
is undetermined and independent of
$z$
, consistent with the standard shallow-water hypothesis. At the next order as
$\varepsilon \to 0$
, the shallow-water equations (1.1)–(1.2) emerge with a viscous term involving horizontal derivatives of the leading-order flow
$\boldsymbol {u}$
. Indeed, the viscous term that emerges is simply (1.3) with
$\varsigma =-2$
.
These modelling strategies can be extended to thin stratified layers with magnetic fields, as often occur in planetary and stellar atmospheres and interiors. Motivated by considerations of the solar tachocline, the equations of shallow-water magnetohydrodynamics (SWMHD) were introduced by Gilman (Reference Gilman2000). For an inviscid and perfectly conducting fluid, he showed that the extension of the system (1.1)–(1.2) is
\begin{align}\partial _t \boldsymbol {u} + \boldsymbol {u} \cdot {\boldsymbol {\nabla }} \boldsymbol {u} = {\boldsymbol {b}} \cdot {\boldsymbol {\nabla }} {\boldsymbol {b}} - g {\boldsymbol {\nabla }} h + \boldsymbol {F},\end{align}
\begin{align} \partial _t {\boldsymbol {b}} + \boldsymbol {u} \cdot {\boldsymbol {\nabla }} {\boldsymbol {b}} = {\boldsymbol {b}} \cdot {\boldsymbol {\nabla }} \boldsymbol {u},\end{align}
\begin{align} \partial _t h + {\boldsymbol {\nabla }} \cdot (h {\boldsymbol {u}}) = 0, \end{align}
where
${\boldsymbol {b}}(\boldsymbol {x},t)$
is the horizontal magnetic field (measured in units of the Alfvén speed), which, like
$\boldsymbol {u}(\boldsymbol {x},t)$
, can be taken to be independent of
$z$
. Then integrating the three-dimensional solenoidal condition across the fluid layer gives
\begin{equation} {\boldsymbol {\nabla }} \cdot \left ( h {\boldsymbol {b}} \right ) = 0, \end{equation}
upon assuming that the free surface is composed of magnetic field lines, and that there is no normal (i.e. vertical) field at the flat bottom. As shown by Dellar (Reference Dellar2002), (1.4)–(1.6) may be cast in a conservative form for the variables
$h {\boldsymbol {u}}$
,
$h {\boldsymbol {b}}$
and
$h$
. In this form (see below), it is immediately clear that (1.7) is consistent with (1.5) and (1.6); that is, if
${\boldsymbol {\nabla }} \cdot (h {\boldsymbol {b}})=0$
holds initially, then it will remain so. The equations of SWMHD have been used to model waves and instabilities in various geophysical and astrophysical settings (e.g. Schecter et al. Reference Schecter, Boyd and Gilman2001; Gilman & Dikpati Reference Gilman and Dikpati2002; Zaqarashvili et al. Reference Zaqarashvili, Oliver, Ballester and Shergelashvili2008; Hunter Reference Hunter2015; Mak et al. Reference Mak, Griffiths and Hughes2016; Márquez Artavia et al. Reference Márquez Artavia, Jones and Tobias2017), although, since none of these settings involve a free surface, either
$g$
in (1.4) should be interpreted as a reduced gravity
$g'$
(Gilman Reference Gilman2000), or the layer depth should be interpreted as an equivalent depth, as in Mak et al. (Reference Mak, Griffiths and Hughes2016).
Just as the hydrodynamic shallow-water equations have been extended to include a diffusive term to account for viscosity, it is natural to ask how the equations of SWMHD can be extended to include a diffusive term to account for finite conductivity. Indeed, the means by which magnetic (Ohmic) diffusion is implemented is arguably more important than how viscous diffusion is implemented, because (1.4) could involve balances between any combination of advection, pressure gradients, the Lorentz force and possibly Coriolis terms, with viscous diffusion playing a minor role. However, the extended shallow-water induction equation would involve only advection and diffusion, and so the consequences of implementing either of these terms erroneously could be serious. In particular, one might be concerned how the form of a magnetic diffusion term influences dynamo action in SWMHD.
To be precise, we introduce a dissipative term
$\boldsymbol {d}(\boldsymbol {x},t)$
in (1.5), as
\begin{equation} \partial _t {\boldsymbol {b}} + \boldsymbol {u} \cdot {\boldsymbol {\nabla }} {\boldsymbol {b}} = {\boldsymbol {b}} \cdot {\boldsymbol {\nabla }} \boldsymbol {u} + \boldsymbol {d} \end{equation}
or, equivalently, as
\begin{equation} \partial _t (h {\boldsymbol {b}}) = {\boldsymbol {\nabla }} \times (\boldsymbol {u} \times h{\boldsymbol {b}}) + h \boldsymbol {d}, \end{equation}
using (1.6). When
$\boldsymbol {d}=0$
, (1.9) is in the form given by Hunter (Reference Hunter2015), and can be reduced to equation (18c) of Dellar (Reference Dellar2002). When
$\boldsymbol {d} \neq 0$
, it provides an immediate constraint on the form of
$\boldsymbol {d}$
, since taking the divergence of (1.9) and using (1.7) gives
\begin{equation} {\boldsymbol {\nabla }} \cdot (h \boldsymbol {d}) = 0. \end{equation}
A second constraint can be derived by considering the domain-integrated energy equation
\begin{equation} \frac {\textrm {d} {E}}{\textrm {d} {t}} = \int h \boldsymbol {u} \cdot {\boldsymbol F} \, \textrm {d} S + \int h {\boldsymbol {b}} \cdot \boldsymbol {d} \, \textrm {d} S, \; \; \mbox {with} \; \; E=\tfrac {1}{2} h \bigl ( \boldsymbol {u}^2 + {\boldsymbol {b}}^2 \bigr ) + \tfrac {1}{2} g h^2, \end{equation}
where
$\textrm {d}S$
is the two-dimensional area element, and we have taken the boundary energy fluxes to vanish, which is guaranteed for appropriate lateral boundary conditions, or for an unbounded flow with
$|\boldsymbol {u}| \to 0$
and
$|{\boldsymbol {b}}| \to 0$
as
$|\boldsymbol {x}| \to \infty$
. We require the Ohmic dissipation to be negative semidefinite, i.e.
\begin{equation} \int h {\boldsymbol {b}} \cdot \boldsymbol {d} \, \textrm {d} S \leqslant 0. \end{equation}
As noted by Mak (Reference Mak2013), for the straightforward choice
$\boldsymbol {d}=\eta \nabla ^2 \boldsymbol {b}$
one cannot prove that either (1.10) or (1.12) is satisfied. The former failure is particularly significant: setting
$\boldsymbol {d}=\eta \nabla ^2 \boldsymbol {b}$
introduces a fundamental inconsistency in the SWMHD formulation, since the constraint (1.7) is not satisfied. This simple diffusion was used in the numerical simulations of SWMHD by Lillo et al. (Reference Lillo, Mininni and Gómez2005), whose results should be treated with caution: in particular, the SWMHD dynamo action they reported could be unphysical.
What forms of
$\boldsymbol {d}$
are consistent with (1.10) and (1.12)? Some first steps in this direction were taken by Mak (Reference Mak2013), who noted that
$\boldsymbol {d}=\eta h^{-1} \nabla ^2 (h {\boldsymbol {b}})$
satisfies (1.10), but also that (1.12) will not be satisfied, in general. One can do better by considering the form
\begin{equation} \boldsymbol {d} = -\frac {1}{h} \, {\boldsymbol {\nabla }} \times \bigl [ \eta h^p {\boldsymbol {\nabla }} \times ( h^q {\boldsymbol {b}} ) \bigr ], \end{equation}
for some
$p$
and
$q$
, and where
$\eta$
may vary horizontally. By construction, this automatically satisfies (1.10). Then, again assuming that the lateral boundary fluxes vanish (e.g. by
$|{\boldsymbol {b}}| \to 0$
as
$|\boldsymbol {x}| \to \infty$
), the Ohmic dissipation
\begin{equation} \int h {\boldsymbol {b}} \cdot \boldsymbol {d} \, \textrm {d} S = -\int \eta h^p ({\boldsymbol {\nabla }} \times {\boldsymbol {b}}) \cdot ({\boldsymbol {\nabla }} \times h^q {\boldsymbol {b}}) \, \textrm {d} S. \end{equation}
Thus (1.12) is certainly satisfied when
$q=0$
. Just as in the case of the viscous diffusion ansatz (1.3), which satisfies the necessary physical constraints when
$\varsigma \leqslant 1$
, we now have a magnetic diffusion ansatz (1.13) that satisfies the necessary physical constraints when
$q=0$
with
$p$
arbitrary. If
$\eta$
has dimensions of
$\textrm {L}^2 \, \textrm {T}^{-1}$
(i.e. it is a diffusivity), then we would need to take
$p=1$
on dimensional grounds, giving
\begin{equation} \boldsymbol {d} = -\frac {1}{h} \, {\boldsymbol {\nabla }} \times \left ( \eta h {\boldsymbol {\nabla }} \times {\boldsymbol {b}} \right ). \end{equation}
So, starting from the ansatz (1.13), we have argued for a plausible form (1.15) for
$\boldsymbol {d}$
. Our main aims here are to show that (1.15) can also be derived systematically by an asymptotic analysis of the three-dimensional induction equation, and to explore some implications of this form for the equations of SWMHD, particularly with dynamo action in mind.
We start, in § 2, by returning to the straightforward choice
$\boldsymbol {d}=\eta \nabla ^2 {\boldsymbol {b}}$
, and investigating the possibility of SWMHD dynamo action. This straightforward choice was adopted by Lillo et al. (Reference Lillo, Mininni and Gómez2005), who considered the SWMHD evolution of forced helical turbulent flows. Here, in order to isolate and understand more clearly any dynamo action in the SWMHD system, we consider the simpler case of the shallow-water analogue of the circularly polarised (CP) flow of Galloway & Proctor (Reference Galloway and Proctor1992) – a flow that has received considerable attention in dynamo studies. Using numerical simulations, we show that SWMHD dynamo action is indeed possible for a range of
$\eta$
. Furthermore, we are able to make comparison with the corresponding magnetohydrodynamical dynamo resulting from the Galloway & Proctor (Reference Galloway and Proctor1992) flow. Whether or not the SWMHD dynamo action is physically realistic is another matter. In § 3, we return to the full three-dimensional induction equation with a three-dimensional Laplacian diffusion, and perform an asymptotic analysis for a thin fluid layer with appropriate conditions on the magnetic field at the free surface and bottom. The ideas here are analogous to those used by Marche (Reference Marche2007) to derive a physically consistent viscous diffusion term for the hydrodynamic shallow-water equations. The outcome of our calculation is a set of equations for SWMHD with an expression for
$\boldsymbol {d}$
that is consistent with both the shallow-water solenoidal constraint (1.10) and the requirement of negative semidefinite Ohmic dissipation (1.12). In § 3.2, we set out some properties of the magnetic diffusion term in more detail, and establish a simple type of antidynamo theorem, thus confirming that the SWMHD dynamo action reported in § 2 is spurious, and arises solely owing to the choice
$\boldsymbol {d}=\eta \nabla ^2 {\boldsymbol {b}}$
. In § 3.3, we revisit the Galloway–Proctor flow numerically, but now with the correct form of the magnetic diffusion; in stark contrast to the exponential growth of magnetic energy with
$\boldsymbol {d}=\eta \nabla ^2 {\boldsymbol {b}}$
, the magnetic energy now decays exponentially. In § 4, we give detailed expressions for the components of the physically consistent magnetic diffusion term in spherical geometry, given the importance of this for astrophysical applications. We conclude in § 5.
2. Shallow-water ‘dynamo action’
As discussed in the introduction, one might be tempted to include magnetic diffusion in the SWMHD induction equation simply through the addition of an
$\eta \nabla ^2 {\boldsymbol {b}}$
term, thus mimicking the diffusion term in the full induction equation. This is the form adopted by Lillo et al. (Reference Lillo, Mininni and Gómez2005), who considered, as a basic state flow, a highly time-dependent hydrodynamical shallow-water flow driven by a large-scale helical forcing. They then showed that the introduction of a weak seed field leads to the growth and subsequent saturation of magnetic energy. It is though hard to draw any detailed conclusions about this particular SWMHD dynamo, since the values of the key parameters, the fluid and magnetic Reynolds numbers, are not provided. In this section, therefore, we look in more detail at the evolution of the magnetic field under the assumption that the magnetic diffusion takes the form
$\eta \nabla ^2 {\boldsymbol {b}}$
. Incompressible, two-dimensional planar flows cannot support dynamo action (Zeldovich Reference Zeldovich1957). Thus, to exhibit dynamo action in the SWMHD equations requires flows with a possibly appreciable variation in height; attaining numerical stability is then not straightforward, but is more readily achieved for unsteady flows. To make contact with classical investigations of dynamo action in incompressible fluids, we shall therefore consider an unsteady, forced shallow-water flow related to a particular incompressible flow widely used in dynamo studies. In § 2.1 we describe briefly the kinematic dynamo properties resulting from solution of the full (three-dimensional) induction equation; in § 2.2 we describe the kinematic properties of what might be regarded as the analogous SWMHD dynamo.
2.1. Classical dynamo action driven by a two-dimensional flow
The kinematic dynamo problem – in which the flow is prescribed and the field evolves solely under the induction equation – is simplified by considering two-dimensional flows – i.e. flows that are invariant in one Cartesian direction. For such flows, as we shall see presently, it is possible to draw an analogy with shallow-water ‘dynamo action’. If the velocity is incompressible, it may be expressed as
\begin{equation} \tilde {\boldsymbol {u}} = \widetilde {\boldsymbol {\nabla }} \times ( \psi {\hat {\boldsymbol {z}}} )+ w {\hat {\boldsymbol {z}}}, \end{equation}
where
$\psi$
and
$w$
are functions of
$x$
,
$y$
and
$t$
. Here we use a tilde to denote three-dimensional vector fields; unless otherwise stated, unadorned quantities represent vector fields with components only in the
$(x,y)$
-plane, as in § 1. Likewise we have
${\boldsymbol {\nabla }} = \hat {\boldsymbol {x}} \partial _x + \hat {\boldsymbol {y}} \partial _y$
as the planar operator and
$\widetilde {\boldsymbol {\nabla }} = \hat {\boldsymbol {x}} \partial _x + \hat {\boldsymbol {y}} \partial _y + \hat {\boldsymbol {z}} \partial _z$
in three dimensions.
A widely studied example of the form (2.1) is the unsteady flow introduced by Galloway & Proctor (Reference Galloway and Proctor1992), in their study of fast dynamo action, with
\begin{equation} \psi = w = A \bigl ( \cos \left ( x + \cos t \right ) + \sin \left ( y + \sin t \right ) \bigr ). \end{equation}
We note that the vorticity is parallel to the velocity: the flow is said to be Beltrami, or maximally helical. For incompressible flows, the induction equation, in dimensionless form, may be written as
\begin{equation} \frac {\partial \tilde {\boldsymbol {b}}}{\partial t} + {\tilde {\boldsymbol {u}}} \cdot \widetilde {\boldsymbol {\nabla }} {\tilde {\boldsymbol {b}}} = {\tilde {\boldsymbol {b}}} \cdot \widetilde {\boldsymbol {\nabla }} {\tilde {\boldsymbol {u}}} + \hat \eta {\widetilde \nabla }^2 {\tilde {\boldsymbol {b}}}, \end{equation}
where
$\hat \eta$
is the (constant) dimensionless magnetic diffusivity, which is inversely proportional to the magnetic Reynolds number
$Rm$
. In the kinematic regime, for flows that are independent of
$z$
, the magnetic field may be expressed in the form
\begin{equation} {\tilde {\boldsymbol {b}}}(x,y,z,t) = {\hat {\boldsymbol {b}}}(x,y,t) \exp (ikz). \end{equation}
For a given wavenumber
$k$
, therefore, the problem involves only two spatial dimensions,
$x$
and
$y$
. The induction equation (2.3) is solved numerically as an initial value problem, using a pseudospectral spatial representation in conjunction with second-order exponential time differencing with Runge–Kutta time stepping (scheme ETD2RK from Cox & Matthews (Reference Cox and Matthews2002)). After any initial transient, the magnetic field grows or decays, with an accompanying oscillation, with growth rate
$s$
. For the particular case of
$A=1.5$
and
$\hat \eta ^{-1} = 100$
, the mode of maximum growth rate has wavenumber
$k=0.61$
and dynamo growth rate
$s=0.38$
. Contours of the
$z$
-components of the magnetic field and the electric current (
${\tilde {\boldsymbol {j}}} =\widetilde {\boldsymbol {\nabla }} \times {\tilde {\boldsymbol {b}}}$
) are shown in figure 1, highlighting their fine-scale structure.
Contour plots on a plane
$z = \textrm {constant}$
of the long-term kinematic solutions for (a)
${\tilde {\boldsymbol {b}}} \cdot {\hat {\boldsymbol {z}}}$
and (b)
${\tilde {\boldsymbol {j}}} \cdot {\hat {\boldsymbol {z}}}$
, for the flow (2.2), with
$A=1.5$
,
$\hat \eta ^{-1} = 100$
and wavenumber
$k=0.61$
. The values are normalised such that
$\mathrm {max} |{\tilde {\boldsymbol {b}}} \cdot {\hat {\boldsymbol {z}}} | =1$
. The calculation was performed with
$256$
Fourier modes in each direction.
2.2. Shallow-water Galloway–Proctor dynamo
For comparison, we now address the kinematic evolution of the magnetic field in a forced, dissipative shallow-water system. We solve, numerically, (1.4) with the addition of forcing and viscous terms to the right-hand side but excluding the Lorentz force, (1.5) with the addition of a magnetic diffusion term to the right-hand side, and (1.6). As discussed above, we are here exploring the implications of expressing the magnetic diffusion term as a Laplacian. For simplicity, and also because it is widely adopted in shallow-water studies, we choose chiefly to employ a two-dimensional Laplacian operator also for the viscous diffusion. Since our focus in this paper is on the evolution of the magnetic field, we do not anticipate that the particular choice of diffusion for the velocity will be a critical factor. We shall, however, briefly address the case when the viscous dissipation takes the form (1.3), with
$\varsigma =-2$
.
We thus first consider the equations
\begin{align} \partial _t \boldsymbol {u} + \boldsymbol {u} \cdot {\boldsymbol {\nabla }} \boldsymbol {u} = - g {\boldsymbol {\nabla }} h + {\boldsymbol {P}} + \nu \nabla ^2 {\boldsymbol {u}}, \end{align}
\begin{align} \partial _t {\boldsymbol {b}} + \boldsymbol {u} \cdot {\boldsymbol {\nabla }} {\boldsymbol {b}} = {\boldsymbol {b}} \cdot {\boldsymbol {\nabla }} \boldsymbol {u} + \eta \nabla ^2 {\boldsymbol {b}}, \end{align}
\begin{align} \partial _t h + {\boldsymbol {\nabla }} \cdot (h {\boldsymbol {u}}) = 0 , \end{align}
where
$\boldsymbol {P}$
denotes the forcing term and
$\nu$
and
$\eta$
denote the (constant) kinematic viscosity and magnetic diffusivity. In dimensionless form, on scaling velocities and horizontal lengths with representative values
$U$
and
$L$
, and fluid depth with the undisturbed depth
$H$
, these may be written as
\begin{align} \partial _t \boldsymbol {u} + \boldsymbol {u} \cdot {\boldsymbol {\nabla }} \boldsymbol {u} = - F^{-2} {\boldsymbol {\nabla }} h + {\boldsymbol {P}} + \hat \nu \nabla ^2 {\boldsymbol {u}}, \end{align}
\begin{align} \partial _t {\boldsymbol {b}} + \boldsymbol {u} \cdot {\boldsymbol {\nabla }} {\boldsymbol {b}} = {\boldsymbol {b}} \cdot {\boldsymbol {\nabla }} \boldsymbol {u} + \hat \eta \nabla ^2 {\boldsymbol {b}}, \end{align}
\begin{align} \partial _t h + {\boldsymbol {\nabla }} \cdot (h {\boldsymbol {u}}) = 0, \end{align}
where
$F=U/\sqrt {gH}$
is the Froude number,
$\hat \nu = \nu /UL$
and
$\hat \eta =\eta /UL$
are scaled diffusivities (inversely proportional to the Reynolds number
$Re$
and magnetic Reynolds number
$Rm$
, respectively), and
$\boldsymbol {P}$
is now the dimensionless forcing.
To draw an analogy with the dynamo described in § 2.1, we suppose that the system is forced by the horizontal projection of the body force that in an incompressible fluid would (at least for sufficiently small fluid Reynolds number) lead to the Galloway–Proctor flow (2.2). Since the flow is incompressible and maximally helical (thus with
$\tilde {\boldsymbol {u}} \cdot \widetilde {\boldsymbol {\nabla }} \tilde {\boldsymbol {u}} = \frac {1}{2} \widetilde {\boldsymbol {\nabla }} \tilde {\boldsymbol {u}}^2$
), it is driven by the forcing
$\widetilde {\boldsymbol {P}} = (\partial _t - \hat \nu \widetilde {\nabla }^2) \tilde {\boldsymbol {u}}$
(see e.g. Cattaneo & Hughes Reference Cattaneo and Hughes1996). Thus, for the shallow-water system, we adopt the forcing
${\boldsymbol {P}} = (P_x, P_y) =(\widetilde P_x, \widetilde P_y)$
using the horizontal components of
$\widetilde {{\boldsymbol {P}}}$
given by
\begin{align} \widetilde P_x &= A \bigl ( (-\cos t \sin (\sin t) + \hat \nu \cos (\sin t)) \cos y - (\cos t \cos (\sin t) + \hat \nu \sin (\sin t)) \sin y \bigr ), \end{align}
\begin{align} \widetilde P_y &= A\bigl ((-\sin t \cos (\cos t) + \hat \nu \sin (\cos t)) \cos x + (\sin t \sin (\cos t) + \hat \nu \cos (\cos t)) \sin x \bigr ). \end{align}
Starting from an initial condition of uniform depth
$h$
(
$\equiv 1$
), zero velocity and zero magnetic field, (2.8) and (2.10) are first evolved in time, on a
$2 \pi \times 2 \pi$
domain, until a stationary, purely hydrodynamic state is attained. As an illustrative example, we again consider the specific case of
$A=1.5$
, for comparison with the Galloway–Proctor dynamo discussed in § 2.1, and take
$F=\sqrt {2/3}$
,
$\hat \nu =0.1$
. We again employ a pseudospectral Fourier representation with ETD2RK time stepping, now with
$512$
Fourier modes in each direction. The flow evolves to a periodic state, with
$\langle h^2 \rangle ^{1/2} = 1.19$
,
$\langle {\boldsymbol {u}}^2 \rangle = 2.09$
,
$\langle h {\boldsymbol {u}}^2 \rangle = 1.89$
, where angle brackets denote an average over
$x$
,
$y$
and
$t$
. Snapshots of the
$z$
-component of the vorticity,
$q$
, say, and the height
$h$
in the hydrodynamic stationary state are shown in figure 2.
Snapshots of contours of (
$a$
)
$q$
and (
$b$
)
$h$
, in the stationary shallow-water hydrodynamic state resulting from the forcing (2.11) with
$A=1.5$
,
$F=\sqrt {2/3}$
,
$\hat \nu =0.1$
.
Long-term kinematic evolution of
$\langle h {\boldsymbol {b}}^2 \rangle$
for the hydrodynamic flow resulting from the forcing (2.11) with
$A=1.5$
,
$F=\sqrt {2/3}$
,
$\hat \nu =0.1$
, with Laplacian viscosity and with Laplacian diffusion for the magnetic field. The different curves are for
$(a)$
$\hat \eta ^{-1}=5$
,
$(b)$
$\hat \eta ^{-1}=10$
,
$(c)$
$\hat \eta ^{-1}=20$
,
$(d)$
$\hat \eta ^{-1}=100$
.
Snapshots of contours of the (exponentially growing) (a)
${\tilde {\kern-2pt \boldsymbol j}} \cdot {\boldsymbol {\hat z}}$
and (b)
$q$
, for the kinematic field evolution driven by the stationary hydrodynamic flow resulting from the forcing (2.11) with
$A=1.5$
,
$F=\sqrt {2/3}$
,
$\hat \nu =0.1$
,
$\hat \eta =0.1$
, and with Laplacian diffusion for the magnetic field. In (
$a$
), the values have been normalised; the values themselves are immaterial in a kinematic field evolution.
Long-term kinematic evolution of
$\langle h {\boldsymbol {b}}^2 \rangle$
for the hydrodynamic flow resulting from the forcing (2.11) with
$A=1.5$
,
$F=\sqrt {2/3}$
,
$\hat \nu =0.05$
, and with viscous diffusion given by (1.3), with
$\varsigma =-2$
. In (a)
$\hat \eta ^{-1}=30$
; in (b)
$\hat \eta ^{-1}=50$
.
To explore the kinematic evolution of the magnetic field, we introduce a seed field of zero mean into the hydrodynamic flow and solve (2.8)–(2.10). The long-time behaviour is characterised by exponential (and oscillatory) growth or decay. Figure 3 shows the exponential growth of magnetic energy versus time for a range of values of
$\hat \eta ^{-1}$
; note that the dependence of the growth rate on
$\hat \eta$
is non-monotonic. As a comparison with the Galloway–Proctor dynamo described in § 2.1, the dynamo growth rate (half the growth rate of the magnetic energy) for
$\hat \eta ^{-1}=10$
is given by
$s=0.11$
, and for
$\hat \eta ^{-1}=100$
,
$s=0.022$
. Snapshots of the
$z$
-components of the electric current and the vorticity for the case of
$\hat \eta ^{-1}=10$
are shown in figure 4. As noted above, with Laplacian diffusion for the magnetic field, the constraint
${\boldsymbol {\nabla }} \cdot (h {\boldsymbol {b}})=0$
is not satisfied; thus, for the shallow-water dynamos shown in figure 3,
${\boldsymbol {\nabla }} \cdot (h {\boldsymbol {b}})$
grows exponentially in time.
To confirm our belief that shallow-water dynamo action is not dependent on the precise form of viscous dissipation adopted – particularly since the motions are driven by an arbitrary forcing – but is a consequence of the combination of the height and induction equations, we have also explored the magnetic field evolution when the flow is again driven by the forcing (2.11), but now with the dissipative term given by (1.3), with
$\varsigma =-2$
. Figure 5 shows the magnetic energy, plotted logarithmically, versus time for
$\hat \nu =0.05$
and for the two cases of
$\hat \eta ^{-1}=30$
and
$\hat \eta ^{-1}=50$
. The magnetic energy, which is oscillatory, exhibits clear exponential growth, again demonstrating shallow-water dynamo action.
Figures 3 and 5 are indeed reminiscent of plots of kinematic dynamo action, showing the exponential amplification of an infinitesimally weak magnetic field. This shallow-water dynamo is, however, a very different beast to its classical counterpart, as can be seen by comparison of the induction equations (2.3) and (2.9). In (2.3),
$\tilde {\boldsymbol {b}}$
is solenoidal and magnetic field growth depends crucially on the field being three-dimensional; if
$k=0$
, then, by a Cartesian analogue of Cowling’s theorem forbidding dynamo-generated axisymmetric fields (Cowling Reference Cowling1933), the magnetic energy can only decay. By contrast, in (2.9),
${\boldsymbol {b}}=(b_x,b_y)$
is not solenoidal and has no
$z$
-dependence; the means of field amplification is clearly therefore very different in the two cases. Whereas the term
${\widetilde \nabla }^2 \tilde {\boldsymbol {b}}$
in (2.3) is always dissipative, there is no such guarantee for the corresponding term in (2.9). Can field growth thus be attributed exclusively to the form of the ‘dissipative’ term adopted in (2.9)? It is clearly important therefore to establish precisely what form this term should take, and then to understand its implications. This is our next aim.
3. Asymptotic reduction of the three-dimensional induction equation
In this section, we derive a physically consistent magnetic diffusion term for SWMHD, by performing an asymptotic analysis of the full three-dimensional diffusive induction equation as the aspect ratio
$\varepsilon \to 0$
. Even though we need not consider the hydrodynamic aspects of the flow in detail, it is useful to sketch how the corresponding hydrodynamic analysis as
$\varepsilon \to 0$
leads to a physically consistent viscous diffusion term in the shallow-water equations (Marche Reference Marche2007); also see the analysis of Levermore & Sammartino (Reference Levermore and Sammartino2001) for a closely related system under the rigid-lid approximation. The hydrodynamic analysis has three key requirements, namely that (i) there is zero tangential stress at the free surface; (ii) there is zero tangential stress at the bottom; (iii) the Reynolds number
$Re$
(based on the horizontal length scale) is of order unity as
$\varepsilon \to 0$
. Requirements (ii) and (iii) are generally inappropriate for oceanic flows, where there will be no slip at the bottom, and
$Re \gg 1$
. However, requirements (i) and (ii) are essential for the leading-order horizontal momentum balance
$\partial ^2 {{\boldsymbol {u}}} / \partial {z}^2 =0$
to have a non-trivial solution that is independent of
$z$
(required for a shallow-water-like outcome), whilst requirement (iii) ensures that a viscous diffusion term appears at the next order (in the physically desirable form (1.3), with
$\varsigma =-2$
), alongside the standard terms of the shallow-water momentum equation. Even though the analysis only formally holds for
$Re$
of order unity as
$\varepsilon \to 0$
, this is really just a convenient way of generating a physically consistent diffusion term, and in practice one might still deploy it in numerical simulations at high
$Re$
.
Here we adopt a similar philosophy for the problem of magnetic diffusion in SWMHD. We will thus need boundary conditions on the magnetic field that allow the leading-order equations to have a non-trivial solution that is independent of
$z$
, and assume that the magnetic Reynolds number
$Rm$
is of order unity, even though we might eventually deploy the resulting magnetic diffusion term in numerical simulations at high
$Rm$
.
3.1. Derivation of the magnetic diffusion term
Without approximation, the induction equation for an incompressible flow, the diffusion term and solenoidal condition may be written as
\begin{align} \partial _t \tilde {\boldsymbol {b}} + \tilde {\boldsymbol {u}} \cdot \widetilde {\boldsymbol {\nabla }} \tilde {\boldsymbol {b}} - \tilde {\boldsymbol {b}} \cdot \widetilde {\boldsymbol {\nabla }} \tilde {\boldsymbol {u}} = \tilde {\boldsymbol {d}}, \end{align}
\begin{align} \tilde {\boldsymbol {d}} = - \widetilde {\boldsymbol {\nabla }} \times ( \eta \widetilde {\boldsymbol {\nabla }} \times \tilde {\boldsymbol {b}}) , \end{align}
\begin{align} \widetilde {\boldsymbol {\nabla }} \cdot \tilde {\boldsymbol {b}}= 0, \end{align}
where, as in § 2, we use a tilde to denote three-dimensional vector fields and operators. We allow a spatially dependent magnetic diffusivity, but take this to be independent of the vertical coordinate, i.e.
$\eta = \eta (x,y)$
. Equations (3.1)–(3.3) are to be solved in a plane layer of fluid,
$0\leqslant z \leqslant h(x,y,t)$
.
The boundary conditions on
$\tilde {\boldsymbol {b}}$
at
$z=0$
and
$z=h(x,y,t)$
depend upon the assumed form of
$\tilde {\boldsymbol {b}}$
and the electric field
$\tilde {\boldsymbol {E}}$
outside the fluid layer. We assume a perfectly conducting exterior with zero magnetic field, in which case
$\tilde {\boldsymbol {b}}=0$
and
$\tilde {\boldsymbol {E}}=0$
for both
$z\lt 0$
and
$z\gt h(x,y,t)$
. The boundary conditions then follow upon integrating
$\widetilde {\boldsymbol {\nabla }} \cdot \tilde {\boldsymbol {b}}=0$
over a pillbox sitting along the boundary, and applying Faraday’s Law to a thin rectangular contour straddling the boundary. At
$z=0$
, the result is standard:
$\hat {\boldsymbol {z}} \cdot \tilde {\boldsymbol {b}}$
and
$\hat {\boldsymbol {z}} \times \tilde {\boldsymbol {E}}$
both vanish, where
$\hat {\boldsymbol {z}}$
is a unit vector in the vertical. However, the calculation is more subtle at
$z=h(x,y,t)$
, since the integrals must be performed in a frame moving with the interface. Denoting values in this moving frame with primes, and using square brackets to denote a change across the interface, we obtain
\begin{equation} [ \tilde {\boldsymbol {n}} \cdot \tilde {\boldsymbol {b}}' ]=0, \; \; \; [ \tilde {\boldsymbol {n}} \times \tilde {\boldsymbol {E}}' ]=0, \end{equation}
where
$\tilde {{\boldsymbol {n}}}$
is any vector normal to the interface (e.g. Roberts Reference Roberts1967). From Ohm’s law, we can write
$\tilde {\boldsymbol {E}}'=\eta \widetilde {\boldsymbol {\nabla }} \times \tilde {\boldsymbol {b}}' - \tilde {\boldsymbol {u}}' \times \tilde {\boldsymbol {b}}'$
, and since
$\tilde {\boldsymbol {u}}' \cdot \tilde {\boldsymbol {n}}=0$
(the frame moves with the interface), (3.4) implies
\begin{equation} [ \tilde {\boldsymbol {n}} \cdot \tilde {\boldsymbol {b}}' ]=0, \; \; \; \eta \tilde {\boldsymbol {n}} \times \bigl [ \widetilde {\boldsymbol {\nabla }} \times \tilde {\boldsymbol {b}}' \bigr ] = (\tilde {\boldsymbol {n}} \cdot \tilde {\boldsymbol {b}}') \, \bigl [ \tilde {\boldsymbol {u}}' \bigr ]. \end{equation}
But
$\tilde {\boldsymbol {b}}'=\tilde {\boldsymbol {b}}$
(it is frame independent), and, for a perfectly conducting exterior with zero magnetic field, (3.5) reduces to
$\tilde {\boldsymbol {n}} \cdot \tilde {\boldsymbol {b}}=0$
and
$\eta \tilde {\boldsymbol {n}} \times \bigl ( \widetilde {\boldsymbol {\nabla }} \times \tilde {\boldsymbol {b}} \bigr )=0$
at the interface. These are just standard conditions of zero normal field and zero tangential current (the latter can also be demonstrated by integrating (3.1) across the interface and using the Reynolds transport theorem). When
$\eta \neq 0$
, we thus solve (3.1)–(3.3) subject to
\begin{align} \hat {\boldsymbol {z}} \cdot \tilde {\boldsymbol {b}} = 0 , \quad \hat {\boldsymbol {z}} \times ( \widetilde {\boldsymbol {\nabla }} \times \tilde {\boldsymbol {b}} ) = 0 \; \; \mbox {on} \; \; z=0, \end{align}
\begin{align} {\tilde {\boldsymbol {n}}} \cdot \tilde {\boldsymbol {b}} = 0 , \quad {\tilde {\boldsymbol {n}}} \times ( \widetilde {\boldsymbol {\nabla }} \times \tilde {\boldsymbol {b}} ) = 0 \; \; \mbox {on} \; \; z=h(x,y,t). \end{align}
We now consider the shallow-water limit: after an appropriate rescaling based on a fluid depth scale
$H$
and horizontal length scale
$L$
with
$H/L = \varepsilon \ll 1$
, the fluid is confined in the layer with
$ 0 \leqslant z \leqslant h(x,y,t)$
, where
$h$
is the original layer depth scaled by
$H$
. The three-dimensional flow
$\tilde {\boldsymbol {u}}$
and magnetic field
$\tilde {\boldsymbol {b}}$
(both scaled by a representative speed
$U$
) and gradient operator
$\widetilde {\boldsymbol {\nabla }}$
take the form
\begin{equation} \tilde {\boldsymbol {u}} = {\boldsymbol {u}} + \varepsilon w \, \hat {\boldsymbol {z}}, \; \; \; \tilde {\boldsymbol {b}} = {\boldsymbol {b}} + \varepsilon c \, \hat {\boldsymbol {z}}, \; \; \; \widetilde {\boldsymbol {\nabla }} = {\boldsymbol {\nabla }} + \varepsilon ^{-1} \hat {\boldsymbol {z}} \, \partial _z . \end{equation}
Here, as before,
$\boldsymbol {u}$
,
$\boldsymbol {b}$
and
$\boldsymbol {\nabla }$
are the horizontal components of the flow, field and gradient operator, whilst
$\varepsilon w$
,
$\varepsilon c$
and
$\varepsilon ^{-1} \partial _z$
are the vertical components. We take the (surface) normal vector field as
\begin{equation} \tilde {\boldsymbol {n}} = - \varepsilon {\boldsymbol {\nabla }} h + \hat {\boldsymbol {z}}. \end{equation}
Note that
$\boldsymbol {u}$
,
$\boldsymbol {b}$
,
$w$
and
$c$
depend on all of
$(x,y,z,t)$
at the outset. When we expand in powers of
$\varepsilon$
, it will be the leading-order horizontal terms
$\boldsymbol {u}_0$
and
$\boldsymbol {b}_0$
that are
$z$
-independent and which will constitute the fields governed by the SWMHD system.
The three-dimensional induction equation (3.1) and solenoidal condition (3.3) become
\begin{align} \left ( \partial _{t} + {\boldsymbol {u}} \cdot {\boldsymbol {\nabla }} + w \, \partial _{z} \right ) \tilde {\boldsymbol {b}} & = \left ( {\boldsymbol {b}} \cdot {\boldsymbol {\nabla }} + c \, \partial _{z} \right ) \tilde {\boldsymbol {u}} + \tilde {\boldsymbol {d}}, \end{align}
\begin{align} {\boldsymbol {\nabla }} \cdot \boldsymbol {b} + \partial _z c & = 0, \end{align}
where, in (3.10), time has been scaled by the advective time scale
$L/U$
, and
$\tilde {\boldsymbol {d}}$
is the scaled version of the magnetic diffusion term (3.2). This can be expressed in terms of
$\boldsymbol {b}$
and
$c$
by using (3.8) to write
\begin{align} \widetilde {\boldsymbol {\nabla }} \times \tilde {\boldsymbol {b}} & = ( {\boldsymbol {\nabla }} + \varepsilon ^{-1} \hat {\boldsymbol {z}} \partial _{z} ) \times ( {\boldsymbol {b}} + \varepsilon c \hat {\boldsymbol {z}} ) = \varepsilon ^{-1} \hat {\boldsymbol {z}} \times \partial _z{{\boldsymbol {b}}} + {\boldsymbol {\nabla }} \times {\boldsymbol {b}} + \varepsilon {\boldsymbol {\nabla }} c \times \hat {\boldsymbol {z}} , \end{align}
\begin{align} \implies \tilde {\boldsymbol {d}} & = \varepsilon ^{-2} \hat \eta \, \partial _z^2 {{\boldsymbol {b}}} - \varepsilon ^{-1} \hat {\boldsymbol {z}} {\boldsymbol {\nabla }} \cdot (\hat \eta \partial _{z} {\boldsymbol {b}} ) - {\boldsymbol {\nabla }} \times \left ( \hat \eta {\boldsymbol {\nabla }} \times {\boldsymbol {b}} \right ) - \hat \eta \partial _{z} {\boldsymbol {\nabla }} c + \varepsilon \hat {\boldsymbol {z}} {\boldsymbol {\nabla }} \cdot ( \hat \eta {\boldsymbol {\nabla }} c) , \end{align}
where
$\hat \eta (x,y)=\eta /UL$
is the scaled magnetic diffusivity, as in (2.9). Since the first, third and fourth terms on the right-hand side of (3.13) are horizontal whilst the second and fifth terms are vertical, we can split (3.10) into its horizontal and vertical components:
\begin{align} \left ( \partial _{t} + {\boldsymbol {u}} \cdot {\boldsymbol {\nabla }} + w \partial _{z} \right ) {\boldsymbol {b}} & = \left ( {\boldsymbol {b}} \cdot {\boldsymbol {\nabla }} + c \partial _{z} \right ) {\boldsymbol {u}} + \varepsilon ^{-2} \hat \eta \, \partial _z^2 {{\boldsymbol {b}}} - {\boldsymbol {\nabla }} \times \left ( \hat \eta {\boldsymbol {\nabla }} \times {\boldsymbol {b}} \right ) - \hat \eta \partial _{z} {\boldsymbol {\nabla }} c, \end{align}
\begin{align} \left ( \partial _{t} + {\boldsymbol {u}} \cdot {\boldsymbol {\nabla }} + w \partial _{z} \right ) c & = \left ( {\boldsymbol {b}} \cdot {\boldsymbol {\nabla }} + c \partial _{z} \right ) w - \varepsilon ^{-2} {\boldsymbol {\nabla }} \cdot (\hat \eta \partial _z \boldsymbol {b}) + {\boldsymbol {\nabla }} \cdot (\hat \eta {\boldsymbol {\nabla }} c). \end{align}
We turn now to the boundary conditions (3.6) and (3.7), the scaled versions of which are
\begin{align} c &= 0 \; \; \mbox {on} \; \; z =0,\end{align}
\begin{align} - \partial _z \boldsymbol {b} + \varepsilon ^2 {\boldsymbol {\nabla }} c & = 0 \; \; \mbox {on} \; \; z =0,\end{align}
\begin{align} c - \boldsymbol {b} \cdot {\boldsymbol {\nabla }} h & = 0 \; \; \mbox {on} \; \;z=h(x,y,t),\end{align}
\begin{align} - \partial _z \boldsymbol {b} - \varepsilon \hat {\boldsymbol {z}} {\boldsymbol {\nabla }} h \cdot \partial _z \boldsymbol {b} + \varepsilon ^2 {\boldsymbol {\nabla }} c - \varepsilon ^2 {\boldsymbol {\nabla }} h & \times ({\boldsymbol {\nabla }} \times \boldsymbol {b}) + \varepsilon ^3 \hat {\boldsymbol {z}} {\boldsymbol {\nabla }} h \cdot {\boldsymbol {\nabla }} c = 0 \; \; \mbox {on} \; \; z=h(x,y,t),\end{align}
using (3.9). Equation (3.19) can also be split into horizontal and vertical components:
\begin{align} \partial _z \boldsymbol {b} = \varepsilon ^2 \bigl ( {\boldsymbol {\nabla }} c - {\boldsymbol {\nabla }} h \times \left ( {\boldsymbol {\nabla }} \times \boldsymbol {b} \right ) \bigr ) \; \; \mbox {on} \; z=h(x,y,t),\end{align}
\begin{align} {\boldsymbol {\nabla }} h \cdot ( \partial _z \boldsymbol {b} - \varepsilon ^2 {\boldsymbol {\nabla }} c ) = 0 \; \; \mbox {on} \; z=h(x,y,t).\end{align}
All the above is exact, albeit rescaled. We now consider the shallow-water limit, i.e.
$\varepsilon \to 0$
. Although
$\hat \eta$
could, in principle, be chosen to depend upon
$\varepsilon$
as this limit is taken, the natural way for second-order horizontal derivatives in the diffusion term to enter into a shallow-water-like balance of (3.14) is with
$\hat \eta$
independent of
$\varepsilon$
. We thus consider the limit
$\varepsilon \to 0$
, with
$\hat \eta$
of order unity (or equivalently
$Rm$
of order unity). The governing equations are (3.14)–(3.15), with boundary conditions (3.16)–(3.18) and (3.20)–(3.21). Noting that the small parameter in this system is
$\varepsilon ^2$
rather than
$\varepsilon ^1$
, we introduce expansions
\begin{equation} {\boldsymbol {b}} = {\boldsymbol {b}}_0 + \varepsilon ^2 {\boldsymbol {b}}_1 + \cdots , \qquad c=c_0+\varepsilon ^2 c_1 + \cdots . \end{equation}
The hydrodynamic expansions are well known to occur in the same way, i.e.
${\boldsymbol {u}} = {\boldsymbol {u}}_0 + \varepsilon ^2 {\boldsymbol {u}}_1 + \cdots$
and
$h = h_0 + \varepsilon ^2 h_1 + \cdots$
. As is standard in shallow-water systems, the hydrodynamic equations (which we do not give here) may be satisfied by taking
\begin{equation} \partial _z {{\boldsymbol {u}}_0} = 0, \end{equation}
so that incompressibility implies
\begin{equation} w_0 = -z {\boldsymbol {\nabla }} \cdot {\boldsymbol {u}}_0, \end{equation}
having applied
$\tilde {\boldsymbol {u}} \cdot \hat {\boldsymbol {z}}=0$
at
$z=0$
. Then the kinematic condition at
$z=h$
implies
\begin{equation} \partial _t {h_0} + {\boldsymbol {\nabla }} \cdot (h_0 {\boldsymbol {u}}_0) = 0. \end{equation}
Introducing expansions of the form (3.22) into the horizontal induction equation (3.14), the leading-order terms yield
$0=\hat \eta \partial ^2_z{{\boldsymbol {b}}_0}$
. Since
$\partial _z {{\boldsymbol {b}}_0}=0$
at
$z=0$
by (3.17) and at
$z=h$
by (3.20), it follows that
\begin{equation} \partial _z {{\boldsymbol {b}}_0} = 0 \; \; \mbox {for all}\,z. \end{equation}
That is, the leading-order horizontal field
$\boldsymbol {b}_0 = \boldsymbol {b}_0(x,y,t)$
is independent of
$z$
, as is the case for
$\boldsymbol {u}_0$
from (3.23). Then, from (3.11), which implies
$ \partial _z {c_0} = -{\boldsymbol {\nabla }} \cdot {\boldsymbol {b}}_0$
, and (3.16), which implies
$c_0 =0$
on
$z=0$
, we obtain
\begin{equation} c_0 = -z {\boldsymbol {\nabla }} \cdot {\boldsymbol {b}}_0. \end{equation}
Since (3.18) implies
$c_0 = \boldsymbol {b}_0 \cdot {\boldsymbol {\nabla }} h_0$
on
$z=h_0$
, combining with (3.27) yields the appropriate divergence free condition for magnetic field,
\begin{equation} {\boldsymbol {\nabla }} \cdot (h_0 {\boldsymbol {b}}_0)=0 . \end{equation}
At order
$\varepsilon ^0$
, (3.14) yields
\begin{equation} \left ( \partial _{t} + {\boldsymbol {u}}_0 \cdot {\boldsymbol {\nabla }} \right ) {\boldsymbol {b}}_0 = {\boldsymbol {b}}_0 \cdot {\boldsymbol {\nabla }} {\boldsymbol {u}}_0 + \hat \eta \partial _z^2 {{\boldsymbol {b}}_1} - {\boldsymbol {\nabla }} \times \left (\hat \eta {\boldsymbol {\nabla }} \times {\boldsymbol {b}}_0 \right ) - \hat \eta \partial _{z} {\boldsymbol {\nabla }} c_0 , \end{equation}
where we have also used (3.23). There are two distinct ways to proceed at this point. The first approach is to integrate (3.29) over the layer depth to obtain
\begin{equation} h_0 \left ( \partial _{t} + {\boldsymbol {u}}_0 \cdot {\boldsymbol {\nabla }} \right ) {\boldsymbol {b}}_0 = h_0 {\boldsymbol {b}}_0 \cdot {\boldsymbol {\nabla }} {\boldsymbol {u}}_0 - h_0 {\boldsymbol {\nabla }} \times \left ( \hat \eta {\boldsymbol {\nabla }} \times {\boldsymbol {b}}_0 \right ) + \hat \eta \Bigl [ \partial _z {{\boldsymbol {b}}_1} - {\boldsymbol {\nabla }} c_0 \Bigr ]_{z=0}^{h_0}\, . \end{equation}
The terms in the square bracket can be evaluated using the
$O(\varepsilon ^2)$
terms of (3.17) and (3.20), which are
\begin{equation} \partial _z {\boldsymbol {b}}_1 = {\boldsymbol {\nabla }} c_0\; \mbox {at} \; z=0, \end{equation}
\begin{equation} \partial _z {\boldsymbol {b}}_1 = {\boldsymbol {\nabla }} c_0 - {\boldsymbol {\nabla }} h_0 \times \left ( {\boldsymbol {\nabla }} \times {\boldsymbol {b}}_0 \right ) \;\mbox {at} \; z=h_0. \end{equation}
Substituting in (3.30) and combining terms gives
\begin{equation} \left ( \partial _{t} + {\boldsymbol {u}}_0 \cdot {\boldsymbol {\nabla }} \right ) {\boldsymbol {b}}_0 = {\boldsymbol {b}}_0 \cdot {\boldsymbol {\nabla }} {\boldsymbol {u}}_0 - h_0^{-1} {\boldsymbol {\nabla }} \times (\hat \eta h_0 {\boldsymbol {\nabla }} \times {\boldsymbol {b}}_0 ). \end{equation}
This is the key result and goal of this paper, namely the induction equation governing the leading-order horizontal fields
$\boldsymbol {b}_0(x,y,t)$
,
$\boldsymbol {u}_0(x,y,t)$
and
$h_0(x,y,t)$
as
$\varepsilon \to 0$
, with
$\hat \eta$
of order unity. Dropping the zero subscript and returning to unscaled variables, this provides the shallow-water form of the induction equation, namely
\begin{equation} \partial _{t} {\boldsymbol {b}} + {\boldsymbol {u}} \cdot {\boldsymbol {\nabla }} {\boldsymbol {b}} = {\boldsymbol {b}} \cdot {\boldsymbol {\nabla }} {\boldsymbol {u}} + \boldsymbol {d}, \end{equation}
with the physically consistent diffusion term
\begin{equation} \boldsymbol {d} = - h^{-1} {\boldsymbol {\nabla }} \times ( \eta h {\boldsymbol {\nabla }} \times \boldsymbol {b}) , \end{equation}
as in (1.15).
The second approach to deriving (3.33) from (3.29) is to recognise that there is a hidden consistency requirement in the above analysis. This can be made explicit by noting that, with the exception of
$\hat{\eta} \partial _z^2 {\boldsymbol {b}}_1$
, all terms of (3.29) have already been found to be independent of
$z$
. It follows that
$\partial _z^2 {{\boldsymbol {b}}_1}$
must also be independent of
$z$
, so that
$\partial _z {{\boldsymbol {b}}_1}$
is linear in
$z$
. Using (3.31) and (3.32) it follows that
\begin{equation} \partial _z {{\boldsymbol {b}}_1} = {\boldsymbol {\nabla }} c_0 \big |_{z=0} \, \left ( 1 - {z}/{h_0} \right ) + \bigl [ {\boldsymbol {\nabla }} c_0 \big |_{z=h_0} - {\boldsymbol {\nabla }} h_0 \times \left ( {\boldsymbol {\nabla }} \times {\boldsymbol {b}}_0 \right ) \bigr ] \, ({z}/{h_0}) , \end{equation}
and so
\begin{equation} \partial _z^2 {\boldsymbol {b}}_1 = h_0^{-1} \bigl [ {\boldsymbol {\nabla }} c_0 \bigr ]_{z=0}^{h_0} - h_0^{-1} {\boldsymbol {\nabla }} h_0 \times \left ( {\boldsymbol {\nabla }} \times {\boldsymbol {b}}_0 \right ) = \partial _{z} {\boldsymbol {\nabla }} c_0 - h_0^{-1} {\boldsymbol {\nabla }} h_0 \times \left ( {\boldsymbol {\nabla }} \times {\boldsymbol {b}}_0 \right ), \end{equation}
since
$c_0$
is also linear in
$z$
from (3.27). It is then easily checked that substituting (3.37) into (3.29) once more gives (3.33).
Finally, we also need to verify that the vertical component of the induction equation, i.e. (3.15), is satisfied to the same degree of approximation. On substituting expansions of the form (3.22), the leading-order,
$O ( \varepsilon ^{-2} )$
, term of (3.15) is zero as
$\boldsymbol {b}_0$
is independent of
$z$
. At the next order in
$\varepsilon$
, we find
\begin{equation} \left ( \partial _{t} + {\boldsymbol {u}}_0 \cdot{\boldsymbol {\nabla }} + w_0 \partial _{z} \right ) c_0 = \left ( {\boldsymbol {b}}_0 \cdot {\boldsymbol {\nabla }} + c_0 \partial _{z} \right ) w_0 - {\boldsymbol {\nabla }} \cdot (\hat \eta \partial _z \boldsymbol {b}_1) + {\boldsymbol {\nabla }} \cdot (\hat \eta {\boldsymbol {\nabla }} c_0). \end{equation}
We will omit the details, but it can be checked that this equation is satisfied identically. This can be done by taking the divergence of (3.29), using (3.24) and (3.27), and noting that the combination
$\partial _z \boldsymbol {b}_1 - {\boldsymbol {\nabla }} c_0$
is linear in
$z$
with (3.31) holding.
3.2. Properties of the magnetic diffusion term
Having established, from the thin layer approximation to the full three-dimensional system, that a physically consistent diffusion term is (3.35) for the shallow-water induction equation written in the form (3.34), we now check that evolving quantities such as the magnetic energy and magnetic flux have the properties we would expect. Since we have confirmed the magnetic diffusion in the form of (1.15), or (1.13) with
$p=1$
,
$q=0$
, the solenoidal condition
$\boldsymbol {\nabla } \cdot (h\boldsymbol {b})=0$
is preserved in time, while for magnetic energy we have
\begin{equation} \frac {\mathrm d E_M}{{\mathrm d}t} \equiv \frac {\mathrm d}{{\mathrm d}t} \int \tfrac {1}{2} h \boldsymbol {b}^2 \, {{\mathrm d}S} = \int h \boldsymbol {b} \cdot ({\boldsymbol {\nabla }} \boldsymbol {u}) \cdot \boldsymbol {b} \, {{\mathrm d}S} - \int \eta h ({\boldsymbol {\nabla }} \times \boldsymbol {b})^2 \, {{\mathrm d}S} . \end{equation}
Here we adopt the boundary conditions that there is no normal component of
$\boldsymbol {u}$
or
$\boldsymbol {b}$
, and no tangential component of the current
$\eta {\boldsymbol {\nabla }} \times \boldsymbol {b}$
to any curve bounding the region containing fluid in the
$(x,y)$
-plane (exterior perfect conductor). So, in agreement with (1.14), the Ohmic dissipation term is negative semidefinite, as desired.
The diffusion term may be expanded to see its structure; it is convenient to add a term that is zero (from (1.7)) and take
$\eta$
constant to write
\begin{align} \eta ^{-1} \boldsymbol {d} & = {\boldsymbol {\nabla }} \bigl [h^{-1} {\boldsymbol {\nabla }} \cdot (h\boldsymbol {b})\bigr ] - h^{-1} {\boldsymbol {\nabla }} \times ( h {\boldsymbol {\nabla }} \times \boldsymbol {b}) \end{align}
\begin{align} & = \nabla ^2 \boldsymbol {b} + {\boldsymbol {\nabla }} ( \boldsymbol {b} \cdot h^{-1} {\boldsymbol {\nabla }} h) + ({\boldsymbol {\nabla }} \times \boldsymbol {b} ) \times h^{-1} {\boldsymbol {\nabla }} h, \end{align}
which, in components with
$\boldsymbol {b} = b_x {\hat {\boldsymbol {x}}} + b_y {\hat {\boldsymbol {y}}}$
, amounts to
\begin{align} \eta ^{-1} d_x & = \nabla ^2 b_{x} + h^{-1} ( \partial _x h \, \partial _x + \partial _y h \, \partial _y ) b_{x} + \partial _x (h^{-1}\partial _x h) b_x + \partial _x (h^{-1}\partial _y h) b_y , \end{align}
\begin{align} \eta ^{-1} d_y & = \nabla ^2 b_{y} +h^{-1}( \partial _x h\, \partial _x + \partial _y h\, \partial _y ) b_{y} + \partial _y(h^{-1} \partial _x h) b_x + \partial _y(h^{-1} \partial _y h) b_y . \end{align}
We have the usual Laplacian terms plus coupling of the components through the height field.
A more compact formulation is to use the divergence-free condition (1.7) to introduce a flux function
$A$
for the magnetic field, defined by
\begin{equation} h\boldsymbol {b} = {\boldsymbol {\nabla }} \times (A\hat {\boldsymbol {z}}) = (\partial _y A, -\partial _x A, 0) , \end{equation}
and having the physical meaning that the difference in
$A$
between two points in the plane is the amount of horizontal magnetic flux trapped under the surface
$z = h$
between those points, or more strictly vertical posts penetrating the thin layer of fluid at those points. The flux function may then be taken (in an appropriate gauge) to satisfy the advection–diffusion equation
\begin{equation} \partial _t A + \boldsymbol {u}\cdot {\boldsymbol {\nabla }} A = - \eta h \, \hat {\boldsymbol {z}}\cdot {\boldsymbol {\nabla }} \times \bigl [h^{-1} {\boldsymbol {\nabla }} \times (A\hat {\boldsymbol {z}}) \bigr ] , \end{equation}
whose curl is (3.34) with (3.35). This may be written as
\begin{equation} \partial _t A + ( \boldsymbol {u} + \eta h^{-1}{\boldsymbol {\nabla }} h ) \cdot {\boldsymbol {\nabla }} A = \eta \nabla ^2 A , \end{equation}
showing that the effect of the shallow-water geometry is to modify the advection velocity
$\boldsymbol {u}$
by a diffusion-dependent term. In the plane, the equation (3.45) for
$A$
is straightforwardly
\begin{equation} \partial _t A + \boldsymbol {u} \cdot {\boldsymbol {\nabla }} A = \eta ( \nabla ^2 A - h^{-1} \partial _x h \, \partial _x A - h^{-1} \partial _y h\, \partial _y A ) \end{equation}
in Cartesian coordinates, or
\begin{equation} \partial _t A + \boldsymbol {u} \cdot{\boldsymbol {\nabla }} A = \eta (\nabla ^2 A - h^{-1} \partial _r h\, \partial _r A - h^{-1} r^{-2} \partial _\theta h\, \partial _\theta A ) \end{equation}
in polar coordinates.
From the structure of (3.46), it is clear that the maximum value of
$A$
in a domain cannot increase in time, nor the minimum value decrease. Thus the flux between any two points is bounded by the difference between the maximum and minimum of
$A$
at time
$t=0$
. This precludes a growing magnetic eigenfunction in a steady flow
$\boldsymbol {u}$
, or one taking a Floquet form for a time-periodic flow
$\boldsymbol {u}$
. This straightforward antidynamo argument assumes suitable boundary conditions – for example, that
$A$
is constant and independent of time on any component of the boundary so that the normal magnetic field is zero there. A more formal antidynamo theorem, showing that
$A \to 0$
and
$\boldsymbol {b} \to 0$
in a suitable norm for general classes of flows, would be desirable and remains a topic for future study.
3.3. Magnetic field evolution with the correct magnetic diffusion term
Having shown in § 2.2 how it is possible to have kinematic exponential field growth under a flow driven by the forcing (2.11) with a Laplacian diffusion in the induction equation, it behoves us to consider the evolution of the magnetic field, under the same flow, but with the diffusion term (3.35). Figure 6 shows the long-term evolution of the magnetic energy, assuming Laplacian viscosity, for the same values of
$\hat \eta$
as shown in figure 3. The numerical method and resolution are the same as employed in § 2.2. The contrast between figure 3 and figure 6 is marked. With Laplacian diffusion for the magnetic field, the magnetic energy is exponentially growing; by contrast, with the diffusion term (3.35), the magnetic energy decays exponentially. As might be expected, the decay rate increases monotonically with
$\hat \eta$
. Snapshots of the long-term (decaying) forms of the flux function
$A$
and the
$z$
-component of the electric current are shown in figure 7.
Long-term kinematic evolution of
$\langle h {\boldsymbol {b}}^2 \rangle$
for the hydrodynamic flow resulting from the forcing (2.11) with
$A=1.5$
,
$F=\sqrt {2/3}$
,
$\hat \nu =0.1$
, with Laplacian viscosity and with the diffusion term (3.35) for the magnetic field. The different curves are for (a)
$\hat \eta ^{-1}=5$
, (b)
$\hat{\eta}^{-1}=10$
, (c)
$\hat{\eta}^{-1}=15$
, (d)
$\hat{\eta}^{-1}=20$
.
Snapshots of contours of (
$a$
)
$A$
and (
$b$
)
${\tilde {\kern-2pt \boldsymbol j}} \cdot {\boldsymbol {\hat z}}$
, for the kinematic field evolution driven by the stationary hydrodynamic flow resulting from the forcing (2.11) with
$A=1.5$
,
$F=\sqrt {2/3}$
,
$\hat \nu =0.1$
,
$\hat \eta =0.1$
, with Laplacian viscosity and with diffusion for the magnetic field given by (3.35). The plots are normalised such that
$\mathrm {max} |A|=1$
.
4. Spherical geometry
Many astrophysical applications involve flow on a sphere, and so here we consider briefly the form of the equations and the magnetic diffusion term in this geometry. We take the flow and field to be defined on a unit sphere
$S$
given by
$r=1$
in spherical polar coordinates
$(r,\theta ,\phi )$
. The fluid occupies a thin layer bounded by
$r=1$
and
$r = 1+ \varepsilon h(\theta ,\phi ,t)$
with
$\varepsilon \ll 1$
as usual. The flow and field are given by
$\boldsymbol {u}(\theta ,\phi ,t)$
and
$\boldsymbol {b}(\theta ,\phi ,t)$
, with the radial component and dependence on radius removed from consideration. We will derive the equations here using a general formulation, as we need to establish notation and appropriate spherical operators, but the reader may wish instead to read the discussion in Gilman & Dikpati (Reference Gilman and Dikpati2002), which gives the SWMHD system in the form of (4.4)–(4.6) with (4.1), or (4.10)–(4.14).
Here we first set up the equations for a flow and field on a general surface
$S$
embedded in ordinary three-dimensional space, following the approach of Il’in (Reference Il’in1991); see this paper and Gilbert et al. (Reference Gilbert, Riedinger and Thuburn2014) for more detail. We let
$\boldsymbol {n}$
be a unit vector field normal to the surface
$S$
, which is extended just off the surface in such a way that
$\boldsymbol {\nabla } \times {\boldsymbol {n}} = 0$
. In this section we will use
$\boldsymbol {\nabla }$
as the usual operator in the full three-dimensional space rather than
$\widetilde {\boldsymbol {\nabla }}$
as earlier, and use
$\boldsymbol {n}$
in preference to
$\tilde {{\boldsymbol {n}}}$
. Given a scalar field
$\chi$
and a vector field
$\boldsymbol {u}$
defined on the surface
$S$
(in other words vectors
$\boldsymbol {u}(\theta ,\phi )$
that are everywhere tangent to
$S$
), we set
\begin{equation} {\mathrm {curl}_{\mathrm s}}\, \chi = {\boldsymbol {\nabla }} \times (\chi {\boldsymbol {n}}) = - {\boldsymbol {n}}\times {\boldsymbol {\nabla }} \chi , \quad {\mathrm {curl}_{\mathrm v}}\, \boldsymbol {u} = {\boldsymbol {n}} \cdot {\boldsymbol {\nabla }} \times \boldsymbol {u} = - {\boldsymbol {\nabla }} \cdot ({\boldsymbol {n}}\times \boldsymbol {u}), \end{equation}
and we also write
${\mathrm {grad}}\, \chi$
and
${\mathrm {div}}\, \boldsymbol {u}$
for the gradient of
$\chi$
and the divergence of
$\boldsymbol {u}$
taken within the surface. Note that the layer thickness here is not being considered; the geometrical set-up is on the purely two-dimensional surface
$S$
. With these two operators, the Laplacian is defined on scalar functions by
\begin{equation} \nabla ^2 \chi = - {\mathrm {curl}_{\mathrm v}}\, {\mathrm {curl}_{\mathrm s}}\, \chi . \end{equation}
The key result of Il’in (Reference Il’in1991) we use is that the projection, say
$\boldsymbol {\pi }$
, of the
$\boldsymbol {u}\cdot {\boldsymbol {\nabla }} \boldsymbol {u}$
term on the surface
$S$
is given by
\begin{equation} \boldsymbol {\pi } \bigl ( \boldsymbol {u}\cdot {\boldsymbol {\nabla }} \boldsymbol {u} \bigr )= - \boldsymbol {u} \times {\boldsymbol {n}}\, {\mathrm {curl}_{\mathrm v}}\,\boldsymbol {u} + {\mathrm {grad}}\, \tfrac {1}{2} \boldsymbol {u}^2. \end{equation}
Within this framework, the equations for SWMHD on
$S$
take the form
\begin{align}\partial _t \boldsymbol {u} - \boldsymbol {u} \times {\boldsymbol {n}}\, {\mathrm {curl}_{\mathrm v}}\,\boldsymbol {u} +\boldsymbol {b} \times {\boldsymbol {n}}\, {\mathrm {curl}_{\mathrm v}}\,\boldsymbol {b} + {\mathrm {grad}}\, \tfrac {1}{2} (\boldsymbol {u}^2 - \boldsymbol {b}^2) + g\, {\mathrm {grad}}\, h = \boldsymbol {F}, \end{align}
\begin{align} \partial _t \boldsymbol {b} - {\mathrm {curl}_{\mathrm s}} ( {\boldsymbol {n}} \cdot \boldsymbol {u} \times \boldsymbol {b} ) - \boldsymbol {b}\, {\mathrm {div}}\, \boldsymbol {u} + \boldsymbol {u}\, {\mathrm {div}}\, \boldsymbol {b} = \boldsymbol {d}, \end{align}
\begin{align} \partial _t h + {\mathrm {div}} (h\boldsymbol {u}) = 0 , \quad {\mathrm {div}} (h \boldsymbol {b} ) = 0 , \end{align}
with the viscous diffusion term
$\boldsymbol {F}$
and magnetic diffusion term
$\boldsymbol {d}$
.
In spherical geometry, with
${\boldsymbol {n}} =\hat {\boldsymbol {r}}$
on the unit sphere and
$\boldsymbol {u} = u_\theta \hat {\boldsymbol {\theta }} + u_\phi \hat {\boldsymbol {\phi }}$
, we have
\begin{align} {\mathrm {grad}}\, \chi =\partial _\theta \chi \, \hat {\boldsymbol {\theta }} + s^{-1} \partial _\phi \chi \, \hat {\boldsymbol {\phi }}, \quad {\mathrm {div}}\, \boldsymbol {u} = s^{-1} \partial _\theta (su_\theta ) + s^{-1} \partial _\phi u_{\phi }, \end{align}
\begin{align} {\mathrm {curl}_{\mathrm s}}\, \chi = s^{-1} \partial _\phi \chi \,\hat {\boldsymbol {\theta }} - \partial _\theta \chi \, \hat {\boldsymbol {\phi }}, \quad {\mathrm {curl}_{\mathrm v}}\, \boldsymbol {u} = s^{-1} \partial _\theta (su_\phi ) - s^{-1} \partial _\phi u_{\theta }, \end{align}
\begin{align} \boldsymbol {\pi } \bigl ( \boldsymbol {u}\cdot {\boldsymbol {\nabla }} \boldsymbol {u} \bigr ) \!=\! \bigl [ ( u_\theta \partial _\theta \!+ \!s^{-1} u_\phi \partial _\phi ) u_{\theta } - s^{-1} c u_\phi u_\phi \bigr ] \hat {\boldsymbol {\theta }} + \bigl [ ( u_\theta \partial _\theta + s^{-1} u_\phi \partial _\phi ) u_{\phi } + s^{-1} c u_\theta u_\phi \bigr ] \hat {\boldsymbol {\phi }}, \end{align}
where we abbreviate
$s = \sin \theta$
,
$c = \cos \theta$
. We can use these expressions in (4.4)–(4.6) to write down the shallow-water equations as in Gilman & Dikpati (Reference Gilman and Dikpati2002), or expand out all the terms to obtain
\begin{align} \partial _t u_\theta + \boldsymbol {u} \cdot {\boldsymbol {\nabla }} u_{\theta } - s^{-1} c u_\phi u_\phi - \boldsymbol {b}\cdot {\boldsymbol {\nabla }} b_{\theta } + s^{-1} c b_\phi b_\phi + g \partial _\theta h = F_\theta , \end{align}
\begin{align}\partial _t u_\phi + \boldsymbol {u} \cdot {\boldsymbol {\nabla }} u_{\phi } + s^{-1} c u_\theta u_\phi - \boldsymbol {b}\cdot {\boldsymbol {\nabla }} b_{\phi } - s^{-1} c b_\theta b_\phi + s^{-1} g \partial _\phi h = F_\phi , \end{align}
\begin{align} \partial _t b_\theta + \boldsymbol {u} \cdot {\boldsymbol {\nabla }} b_{\theta } - \boldsymbol {b}\cdot {\boldsymbol {\nabla }} u_{\theta } = d_\theta , \end{align}
\begin{align}\partial _t b_\phi + \boldsymbol {u} \cdot {\boldsymbol {\nabla }} b_{\phi } + s^{-1} c u_\phi b_\theta - \boldsymbol {b}\cdot {\boldsymbol {\nabla }} u_{\phi } - s^{-1} c b_\phi u_\theta = d_\phi , \end{align}
\begin{align} \partial _t h + s^{-1} \partial _\theta ( shu_\theta ) + s^{-1} \partial _\phi (hu_\phi ) = 0 , \quad s^{-1} \partial _\theta ( shb_\theta ) + s^{-1} \partial _\phi (hb_\phi ) = 0 , \end{align}
with
$\boldsymbol {u} \cdot {\boldsymbol {\nabla }} = u_\theta \partial _\theta + s^{-1} u_\phi \partial _\phi$
and similarly for
$\boldsymbol {b}\cdot {\boldsymbol {\nabla }}$
.
We now consider the magnetic diffusion term
$\boldsymbol {d}$
; the viscous diffusion term
$\boldsymbol {F}$
is set out in Gilbert et al. (2014). The appropriate generalisation of (3.35) is
\begin{equation} \boldsymbol {d} = - h^{-1}\, {\mathrm {curl}_{\mathrm s}}( \eta h\, {\mathrm {curl}_{\mathrm v}}\, \boldsymbol {b}) . \end{equation}
After integration by parts, the magnetic energy equation, analogous to (3.39), is given by
\begin{align} & \frac {\mathrm d E_M}{{\mathrm d}t} \equiv \frac {\rm d}{{\rm d}t} \int \tfrac {1}{2} h \boldsymbol {b}^2 \, {{\mathrm d}S}\nonumber\\ &\qquad = \int \left [ \boldsymbol {b} \cdot {\mathrm {curl}_{\mathrm s}} ( h {\boldsymbol {n}}\cdot \boldsymbol {u}\times \boldsymbol {b}) + \tfrac {1}{2} \boldsymbol {b}^2\, {\mathrm {div}}(h\boldsymbol {u}) \right ] {{\mathrm d}S} - \int \eta h ({\mathrm {curl}_{\mathrm v}}\, \boldsymbol {b})^2 \, {{\mathrm d}S} , \end{align}
with the dissipative term correctly taking a negative semidefinite form.
For a vector potential defined on the surface by
\begin{equation} h\boldsymbol {b} = {\mathrm {curl}_{\mathrm s}}\, A, \end{equation}
the corresponding
$A$
equation is
\begin{equation} \partial _t A + \boldsymbol {u}\cdot {\boldsymbol {\nabla }} A = - \eta h\, {\mathrm {curl}_{\mathrm v}} (h^{-1}\, {\mathrm {curl}_{\mathrm s}}\, A) = \eta \bigl [ \nabla ^2 A + h^{-1} {\boldsymbol {n}} \cdot {\mathrm {grad}}\, h \times {\mathrm {curl}_{\mathrm s}}\, A\bigr ] \end{equation}
using the scalar Laplacian defined in (4.2). This amounts to
\begin{equation} \partial _t A + \boldsymbol {u} \cdot {\boldsymbol {\nabla }} A = \eta \bigl [ \nabla ^2 A - h^{-1} (\partial _\theta h\, \partial _\theta A + s^{-2} \partial _\phi h\, \partial _\phi A) \bigr ] , \end{equation}
where the Laplacian on the sphere is as usual given by
\begin{equation} \nabla ^2 \chi = \partial _\theta ^2 \chi + s^{-1} \partial _\theta \chi + s^{-2} \partial _\phi ^2 \chi . \end{equation}
For the components of diffusion of the magnetic field in spherical geometry, taking
$\eta$
constant, we add a term that is zero to
$\boldsymbol {d}$
in (4.15) to write
\begin{equation} \eta ^{-1} \boldsymbol {d} = {\mathrm {grad}} \bigl [h^{-1}\, {\mathrm {div}} (h\boldsymbol {b})\bigr ] - h^{-1}\, {\mathrm {curl}_{\mathrm s}}( h\, {\mathrm {curl}_{\mathrm v}}\, \boldsymbol {b}) , \end{equation}
which amounts to
\begin{align} \!\! \eta ^{-1} \boldsymbol {d} &\! = \!\left \{ \partial _\theta \bigl [ h^{-1} s^{-1} \partial _\theta (shb_\theta ) \!+\! h^{-1} s^{-1}\partial _\phi (hb_\phi )\bigr ] - h^{-1} s^{-1} \partial _\phi \bigl [ h s^{-1} \partial _\theta (sb_\phi ) - h s^{-1} \partial _\phi b_{\theta } \bigr ] \right \} \!\hat {\boldsymbol {\theta }} \notag \\ & + \left \{ s^{-1} \partial _\phi \bigl [ h^{-1} s^{-1} \partial _\theta (shb_\theta ) \!+\! h^{-1} s^{-1}\partial _\phi (hb_\phi )\bigr ] \!+\!h^{-1} \partial _\theta \bigl [ h s^{-1} \partial _\theta (sb_\phi )\! -\! h s^{-1} \partial _\phi b_{\theta }\bigr ] \right \}\! \hat {\boldsymbol {\phi }}, \end{align}
and then expand this to obtain
\begin{align} \eta ^{-1} d_\theta & = \nabla ^2 b_\theta - 2 s^{-2} c \partial _\phi b_{\phi } - s^{-2} b_\theta + h^{-1} \partial _\theta h\, \partial _\theta b_{\theta } + s^{-2} h^{-1} \partial _\phi h\, \partial _\phi b_{\theta } \notag \\ & + \partial _\theta (h^{-1} \partial _\theta h) b_\theta + s^{-1} \partial _\theta (h^{-1}\partial _\phi h) b_\phi - 2 s^{-2} c (h^{-1} \partial _\phi h) b_\phi , \end{align}
\begin{align} \eta ^{-1} d_\phi & = \nabla ^2 b_\phi + 2 s^{-2} c \partial _\phi b_{\theta } - s^{-2} b_\phi + h^{-1} \partial _\theta h\, \partial _\theta b_{\phi } + s^{-2} h^{-1} \partial _\phi h\, \partial _\phi b_{\phi } \notag \\ & + s^{-1} \partial _\phi ( h^{-1} \partial _\theta h) b_\theta + s^{-2} \partial _\phi (h^{-1} \partial _\phi h) b_\phi + s^{-1} c (h^{-1} \partial _\theta h ) b_\phi ; \end{align}
we observe numerous coupling terms between the magnetic and height fields.
5. Conclusions
The equations of SWMHD were introduced by Gilman (Reference Gilman2000) as a simplified system for modelling thin stratified fluid layers permeated by a magnetic field. They were derived for an ideal system, namely for an inviscid and perfectly conducting fluid. However, extending the system to allow for the dissipative processes of viscous diffusion and magnetic diffusion is valuable for two reasons. First, these processes exist in nature, will modify flows, waves and instabilities at appropriate length scales, and so may need to be quantified. Second, numerical models will generally need to incorporate dissipation, even if simulating turbulence or complex flows at scales much larger than some nominal dissipative scale.
The appropriate form to take for the magnetic diffusion term is not evident at the outset. Perhaps the most natural route is to place a term
$\boldsymbol {d} = \eta \nabla ^2 \boldsymbol {b}$
in the SWMHD induction equation in line with the full equations of three-dimensional magnetohydrodynamics (MHD), as adopted by Lillo et al. (Reference Lillo, Mininni and Gómez2005). In § 2, we explored the consequences of this, and showed that kinematic dynamo action – exponential growth of magnetic energy – is possible in a two-dimensional planar flow inspired by the Galloway & Proctor (Reference Galloway and Proctor1992) dynamo. However, given that the only processes present in the SWMHD induction equation are advection (or Lie-dragging, see Schutz Reference Schutz1980) of the magnetic field and magnetic diffusion, ensuring that the diffusion term represents the correct physics is crucial. As discussed in the introduction, there are two physical constraints that must be respected: the SWMHD solenoidal condition
$\boldsymbol {\nabla } \cdot (h\boldsymbol {b})=0$
in (1.10), and a negative semidefinite Ohmic dissipation term in (1.12). Unfortunately, the straightforward choice of a magnetic diffusion term
$\boldsymbol {d}=\eta \nabla ^2 \boldsymbol {b}$
violates (1.10), and generally does not respect (1.12) (Mak Reference Mak2013). In this way, the choice
$\boldsymbol {d}=\eta \nabla ^2 \boldsymbol {b}$
is both mathematically and physically inconsistent with the underlying system, and further analysis shows that the dynamo action of § 2 is illusory. This diffusion term redistributes magnetic energy in a way that is unphysical; analogously, an incorrect form of the viscous diffusion term can likewise give spurious sinks and sources of angular momentum (Gilbert et al. Reference Gilbert, Riedinger and Thuburn2014).
One approach to introducing magnetic diffusion in SWMHD is then to take an operator that is required only to satisfy the constraints (1.10) and (1.12). There is a wide possible choice here; for example, a term of the form of (1.13) with any value of
$p$
but with
$q=0$
satisfies these constraints. More satisfactory, though, is to derive systematically an operator with a particular choice of
$p$
from the underlying equations of three-dimensional MHD. In § 3, we showed how a physically consistent magnetic diffusion term can be obtained by an asymptotic reduction of the full three-dimensional induction equation, which results from integrating across the shallow fluid layer. The resulting SWMHD induction equation is
\begin{equation} \partial _{t} {\boldsymbol {b}} + {\boldsymbol {u}} \cdot {\boldsymbol {\nabla }} {\boldsymbol {b}} = {\boldsymbol {b}} \cdot {\boldsymbol {\nabla }} {\boldsymbol {u}} - h^{-1} {\boldsymbol {\nabla }} \times ( \eta h {\boldsymbol {\nabla }} \times \boldsymbol {b}) , \end{equation}
corresponding to the choice
$p=1$
and
$q=0$
in (1.13). As (5.1) is derived from the full three-dimensional equations, it should be consistent with other physics of SWMHD; it can also be used with a spatially varying magnetic diffusivity
$\eta (x,y)$
. With this form of the diffusion operator we derived a simple type of antidynamo theorem in § 3.2, which confirms that the dynamo action found in § 2 (and in Lillo et al. Reference Lillo, Mininni and Gómez2005) is unphysical. Further confirmation is provided by the numerical results in § 3.3. In hindsight, this is perhaps not surprising: while all three components of magnetic field are present in SWMHD, the vertical field is passive and not coupled back into the induction equation. Although there can be plenty of stretching of the horizontal components of the magnetic field in the thin layer, the resulting folding leads to fields with cancelling orientations, and so no net growth of magnetic flux. Lacking are the vertical dependence of the field and vertical motions that could constructively fold field lines, for example through the stretch–fold–shear mechanism (e.g. Bayly & Childress Reference Bayly and Childress1988).
Even if dynamo action is not involved, it is generally not permitted (for constant
$\eta$
) to use Laplacian diffusion in (5.1) because the SWMHD solenoidal constraint
$\boldsymbol {\nabla } \cdot (h {\boldsymbol {b}})=0$
will be violated as the flow evolves. However, if there is a regime in which the free surface perturbations are small, i.e.
$h(\boldsymbol {x},t)=H(1+\delta h'(\boldsymbol {x},t))$
for some constant
$H$
with
$\delta \ll 1$
, then an asymptotic reduction of (5.1) and
$\boldsymbol {\nabla } \cdot (h {\boldsymbol {b}})=0$
can be made as
$\delta \to 0$
. Then, for constant
$\eta$
, the leading-order dissipative term in (5.1) is simply
$\eta \nabla ^2 {\boldsymbol {b}}$
, and the leading-order solenoidal constraint is
$\boldsymbol {\nabla } \cdot {\boldsymbol {b}}=0$
, which together allow a consistent evolution. This now looks like two-dimensional MHD, although the coupled flow
${\boldsymbol {u}}(\boldsymbol {x},t)$
could still have shallow-water effects depending upon how the asymptotic reduction is made. For example, this would be the case for a diffusive extension of the equations of quasigeostrophic SWMHD introduced by Zeitlin (Reference Zeitlin2013), where the small parameter
$\delta$
is the Rossby number, and shallow-water effects are felt through the Rossby radius of deformation in the vorticity equation, with both
$\boldsymbol {\nabla } \cdot {\boldsymbol {u}}=0$
and
$\boldsymbol {\nabla } \cdot {\boldsymbol {b}}=0$
to leading order.
To conclude, we propose that the form (5.1) of the induction equation be used in future studies of SWMHD. Indeed, based on our analysis, (5.1), in its Cartesian form (3.42)–(3.43) and (3.47), has already been adopted in the recent hot Jupiter simulations of Hindle et al. (Reference Hindle, Bushby and Rogers2019, Reference Hindle, Bushby and Rogers2021). Since shallow-water systems are also used for global studies of magnetohydrodynamical waves and instabilities in spherical geometry (e.g. Gilman & Dikpati Reference Gilman and Dikpati2002; Dikpati et al. Reference Dikpati, Gilman and Rempel2003; Márquez Artavia et al. Reference Márquez Artavia, Jones and Tobias2017), we have set out the appropriate form of the magnetic diffusion term in (4.19) and (4.23)–(4.24) for spherical polar coordinates.
Funding.
A.D.G. is grateful to the Leverhulme Trust for a Research Fellowship when this work commenced, and more recently acknowledges support from the EPSRC grant EP/T023139/1. D.W.H. was supported in part by STFC grant ST/N000765/1. D.W.H. would also like to thank the Isaac Newton Institute for Mathematical Sciences (INI), Cambridge, for support and hospitality during the programme ‘Anti-diffusive dynamics: from sub-cellular to astrophysical scales’, where some of the work on this paper was undertaken. The INI programme was supported by EPSRC grant EP/R014604/1.
Declaration of interests.
The authors report no conflict of interest.
Data availability statement.
Data from numerical simulations was used in this study. The data could be reproduced from the details of the numerical simulations (the equations of motion, resolution and parameters) given in § 2.
Open access
