2.1. Equations of compressible MHD

In standard notation, the governing equations for a perfect gas in the absence of viscous, thermal and magnetic diffusivities are given by

(2.1)\begin{gather} \frac{\partial \rho}{\partial t} + \boldsymbol{\nabla}\boldsymbol{\cdot} \left( \rho \boldsymbol{u} \right) = 0, \end{gather}

(2.2)\begin{gather}\rho \left( \frac{\partial \boldsymbol{u}}{\partial t} + \boldsymbol{u}\boldsymbol{\cdot} \boldsymbol{\nabla} \boldsymbol{u} \right) ={-}\boldsymbol{\nabla} p + \rho \boldsymbol{g} + \frac{1}{\mu_0} \left( \boldsymbol{\nabla}\times \boldsymbol{B} \right) \times \boldsymbol{B}, \end{gather}

(2.3)\begin{gather}\frac{\partial \boldsymbol{B}}{\partial t} + \left( \boldsymbol{u} \boldsymbol{\cdot} \boldsymbol{\nabla} \right) \boldsymbol{B} = \left( \boldsymbol{B} \boldsymbol{\cdot} \boldsymbol{\nabla} \right) \boldsymbol{u} - \boldsymbol{B} \left( \boldsymbol{\nabla}\boldsymbol{\cdot} \boldsymbol{u} \right), \end{gather}

(2.4)\begin{gather}\frac{\partial p}{\partial t} + \boldsymbol{u} \boldsymbol{\cdot} \boldsymbol{\nabla} p + \gamma p \boldsymbol{\nabla}\boldsymbol{\cdot} \boldsymbol{u} = 0, \end{gather}

(2.5)\begin{gather}p = \mathcal{R} \rho T, \end{gather}

where $\mathcal {R} = c_p - c_v$ is the specific gas constant; $c_p$ and $c_v$ are the specific heats at constant pressure and volume, respectively; $\gamma = c_p/c_v$ is the ratio of specific heats; $\mu _0$ is the magnetic permeability of free space; $\boldsymbol{g}$ is the gravitational acceleration. Throughout this paper, the geometry under consideration consists of a plane layer of fluid bounded by horizontal planes located at $z = 0$ and $z = d$, with the $z$-axis pointing downwards ($\boldsymbol{g} = g \hat {\boldsymbol{e}}_z$). In Cartesian coordinates, we write the fluid velocity as $\boldsymbol{u} = (u,v,w)$. Note that in the above system, which has no dissipation, the dynamics is described completely through (2.1)–(2.4), which involve only the thermodynamic variables $p$ and $\rho$. It is, however, helpful also to introduce the temperature through the perfect gas equation of state (2.5).

The governing equations (2.1)–(2.5) can be written in dimensionless form by scaling magnetic field, mass density, temperature and pressure with their values at the top of the layer ($z=0$): $B_r$, $\rho _r$, $T_r$ and $p_r = \mathcal {R} \rho _r T_r$, respectively. (The subscript ‘$r$’ is used throughout to denote representative values of quantities.) Furthermore, we scale lengths with the layer depth $d$, time with the acoustic time scale $d/\sqrt {\mathcal {R} T_r}$ and velocities with $\sqrt {\mathcal {R} T_r}$. Representative values for the square of the isothermal sound speed and the square of the Alfvén speed are, respectively, $c_{s,r}^2 = p_r/\rho _r$, $c_{A,r}^2 = B_r^2/(\mu _0 \rho _r)$. A representative pressure scale height in a hydrostatically balanced atmosphere is $H_r = p_r/(\rho _r g)$. This implies the following relation:

(2.6)\begin{equation} c_{s,r}^2 = \frac{p_r}{\rho_r} = \mathcal{R} T_r = g H_r. \end{equation}

The dimensionless equations of compressible MHD then take the form

(2.7)\begin{gather} \frac{\partial \rho}{\partial t} + \boldsymbol{\nabla}\boldsymbol{\cdot} \left( \rho \boldsymbol{u} \right) = 0, \end{gather}

(2.8)\begin{gather}\rho \left( \frac{\partial \boldsymbol{u}}{\partial t} + \boldsymbol{u}\boldsymbol{\cdot} \boldsymbol{\nabla} \boldsymbol{u} \right) ={-} \boldsymbol{\nabla} p + \lambda \rho \, \hat{\boldsymbol{e}}_z + M_A^2 \left( \boldsymbol{\nabla}\times \boldsymbol{B} \right) \times \boldsymbol{B}, \end{gather}

(2.9)\begin{gather}\frac{\partial \boldsymbol{B}}{\partial t} + \left( \boldsymbol{u} \boldsymbol{\cdot} \boldsymbol{\nabla} \right) \boldsymbol{B} = \left( \boldsymbol{B} \boldsymbol{\cdot} \boldsymbol{\nabla} \right) \boldsymbol{u} - \boldsymbol{B} \left( \boldsymbol{\nabla}\boldsymbol{\cdot} \boldsymbol{u} \right), \end{gather}

(2.10)\begin{gather}\frac{\partial p}{\partial t} + \boldsymbol{u} \boldsymbol{\cdot} \boldsymbol{\nabla} p + \gamma p \boldsymbol{\nabla}\boldsymbol{\cdot} \boldsymbol{u} = 0, \end{gather}

(2.11)\begin{gather}p = \rho T, \end{gather}

where

(2.12a,b)\begin{equation} \lambda = \frac{d}{H_r}, \quad M_A^2 = \frac{B_r^2/\mu_0}{p_r}. \end{equation}

Here, $\lambda$ is the ratio of the depth of the fluid to the hydrostatic pressure scale height; as such, it is a measure of atmospheric stratification, with small (large) $\lambda$ indicating weak (strong) stratification. The parameter $M_A$ is the Alfvén Mach number – the ratio of Alfvén speed to isothermal sound speed. Note that $M_A^2 = 2/\beta _p$, where $\beta _p$ – the plasma beta – is the ratio of gas pressure to magnetic pressure.

The equations of compressible MHD (2.7)–(2.11) form the most general set of equations describing MHD behaviour in a stratified, electrically conducting, diffusionless, perfect gas, allowing the description of phenomena that occur on a range of distinct time scales. Of particular note is that the equations describe sound waves (or fast magneto-acoustic waves), the time scale for which is often much shorter than that of other waves or instabilities. Furthermore, in scenarios in which the flows are highly subsonic, the high-frequency sound waves often play no dynamically significant role. It is therefore useful to have a system of equations that filter out sound waves and home in on the dynamics evolving on the slower time scale. Such slow-time dynamics is described by the anelastic approximation.

2.2. The anelastic approximation

The anelastic approximation posits an atmosphere that is almost neutrally buoyant – the reference state – and considers the evolution of perturbations on top of that state. Convective instability is governed by the well-known Schwarzschild criterion, which dictates that instability depends on the gradient of specific entropy $s = c_v \ln (p \rho ^{-\gamma })$: an unstable configuration is one in which specific entropy increases with depth, i.e. $\mathrm {d} s/\mathrm {d} z > 0$. The reference state for the anelastic approximation is thus taken to be one of hydrostatic balance in which entropy is nearly constant (also referred to as near-adiabatic stratification), but with significant variations individually in pressure and density. The small departure from adiabatic stratification induces convection if $\mathrm {d} s/\mathrm {d} z > 0$, or gravity waves if $\mathrm {d} s/\mathrm {d} z < 0$, each of which engenders small perturbations to the reference state.

By treating the departure from adiabatic stratification as a small parameter, one can conduct a formal asymptotic expansion of the compressible equations to extract equations – the anelastic equations – that describe how the perturbations evolve on top of the fixed thermodynamic background. We thus introduce the quantity $\Delta\!\!\nabla$, a measure of departure from adiabaticity, defined by

(2.13)\begin{equation} \Delta\!\!{\nabla} ={-}\frac{d}{H_p} \left( \frac{\mathrm{d} \ln \rho}{\mathrm{d} \ln p} - \left( \frac{\mathrm{d} \ln \rho}{\mathrm{d} \ln p} \right)_{ad} \right), \end{equation}

where

(2.14a,b)\begin{equation} H_p = \left( \frac{\mathrm{d} \ln p}{\mathrm{d} z} \right)^{{-}1} \quad \textrm{and} \quad \left( \frac{\mathrm{d} \ln \rho}{\mathrm{d} \ln p} \right)_{ad} = 1/\gamma \end{equation}

are, respectively, the pressure scale height and the adiabatic gradient. The stratification is superadiabatic (convectively unstable) in regions where $\Delta\!\!\nabla > 0$, and subadiabatic (convectively stable) where $\Delta\!\!\nabla < 0$. The departure from adiabaticity $\Delta\!\!\nabla$ is essentially the dimensionless gradient of specific entropy

(2.15)\begin{equation} \frac{\mathrm{d} s}{\mathrm{d} z} \equiv \frac{c_p}{d} \Delta\!\!{\nabla}. \end{equation}

Given that the fundamental requirement of the anelastic approximation is that the atmospheric stratification is close to adiabatic, this allows us to define a small parameter $\varepsilon$ by

(2.16)\begin{equation} \varepsilon = \left| \Delta\!\!{\nabla}_r \right| \ll 1. \end{equation}

Specific entropy is thus constant to the lowest (zeroth) order: $\mathrm {d} s /\mathrm {d} z = {O}(\varepsilon )$. Since the entropy gradient is the central quantity in the anelastic formulation, it is more convenient to use the entropy formulation of the conservation of internal energy (2.4),

(2.17)\begin{equation} \rho T \left( \frac{\partial s}{\partial t} + \boldsymbol{u}\boldsymbol{\cdot} \boldsymbol{\nabla} s \right) = 0.\end{equation}

The assumed small departure from adiabaticity sets the scale for all thermodynamic perturbations, allowing a decomposition of the thermodynamic variables into $z$-dependent steady reference states, denoted by overbars, and dynamically induced time-dependent fluctuations, denoted by asterisks:

(2.18)\begin{equation} \left.\begin{gathered} p = p_r (\bar{p} + \varepsilon p^*),\quad \rho = \rho_r (\bar{\rho} +\varepsilon \rho^*), \\ T = T_r (\bar{T} + \varepsilon T^*),\quad s = s_r + c_p \varepsilon(\bar{s} + s^*), \end{gathered}\right\} \end{equation}

where $\bar p$, $p^*$, etc. are dimensionless quantities. The scaling for velocity follows from the physical picture that fluid motions arise owing to buoyancy variations; thus, balancing inertia against buoyancy perturbations implies that $u_r^2 \sim \varepsilon g d$. This, in turn, dictates that the time scale is long (slow evolution), with $t \sim d/u_r \sim \varepsilon ^{-1/2} \sqrt {d/g}$. Note that the ordering of velocity means that the flow Mach number (the ratio of flow speed to sound speed) is small, i.e. from (2.6), $M = u_r/c_{s,r} \sim \varepsilon ^{1/2}$, provided that $\lambda$ is ${O}(1)$. In the standard formulation of the anelastic equations, it is assumed that the Lorentz force does not contribute to the hydrostatic balance at lowest order; this implies that the magnetic field must also scale with $\varepsilon$. Balancing the Lorentz force with pressure perturbations yields $B_r^2/\mu _0 \sim \varepsilon p_r$: thus the Alfvén Mach number is small, i.e. $M_A = c_{A,r}/c_{s,r} \sim \varepsilon ^{1/2}$. Based on these considerations, it is useful to introduce scaled parameters, defined by $B_r = \varepsilon ^{1/2} \tilde {B}_r$, $c_{A,r} = \varepsilon ^{1/2} \tilde {c}_{A,r}$, $M_A = \varepsilon ^{1/2} \tilde {M}_A$, where tilde variables and parameters are ${O}(1)$. In terms of the ordering with $\varepsilon$, the natural scalings of length, time, velocity and magnetic field are therefore as follows:

(2.19a–d)\begin{equation} \boldsymbol{x}= d \boldsymbol{x}' , \quad t = \varepsilon^{{-}1/2} \left(\frac{d}{c_{s,r}} \right) t^*, \quad \boldsymbol{u} = \varepsilon^{1/2} c_{s,r} \boldsymbol{u}^* , \quad \boldsymbol{B} = \varepsilon^{1/2} \tilde{B}_r \boldsymbol{B}^*, \end{equation}

where, as earlier, we use asterisks to denote (dimensionless) variables that have been scaled with a power of $\varepsilon$. The formulation proceeds by substituting expressions (2.19a–d) into the compressible MHD equations (2.7)–(2.11) (dropping $'$ superscripts on length variables) and equating terms at successive powers of $\varepsilon$. At ${O}(\varepsilon ^0)$ we obtain non-trivial expressions only from the $z$-component of (2.8) and from (2.11); these define the reference state by

(2.20)\begin{gather} \frac{\mathrm{d} \bar{p}}{\mathrm{d} z} =\lambda \bar{\rho} , \end{gather}

(2.21)\begin{gather}\bar{p} = \bar{\rho} \bar{T}. \end{gather}

At ${O}(\varepsilon ^1)$ we obtain the following evolution equations governing the perturbations:

(2.22)\begin{gather} \boldsymbol{\nabla}\boldsymbol{\cdot} \left( \bar{\rho} \boldsymbol{u}^* \right) = 0, \end{gather}

(2.23)\begin{gather}\bar{\rho} \left( \frac{\partial \boldsymbol{u}^*}{\partial t^*} + \boldsymbol{u}^*\boldsymbol{\cdot} \boldsymbol{\nabla} \boldsymbol{u}^* \right) ={-} \boldsymbol{\nabla} p^* + \lambda \rho^* \hat{\boldsymbol{e}}_z + \tilde{M}_A^2 \left( \boldsymbol{\nabla}\times \boldsymbol{B}^* \right) \times \boldsymbol{B}^*, \end{gather}

(2.24)\begin{gather}\frac{\partial \boldsymbol{B}^*}{\partial t^*} + \left( \boldsymbol{u}^* \boldsymbol{\cdot} \boldsymbol{\nabla} \right) \boldsymbol{B}^* = \left( \boldsymbol{B}^* \boldsymbol{\cdot} \boldsymbol{\nabla} \right) \boldsymbol{u}^* - \boldsymbol{B}^* \left( \boldsymbol{\nabla}\boldsymbol{\cdot} \boldsymbol{u}^* \right), \end{gather}

(2.25)\begin{gather}\frac{\partial s^*}{\partial t^*} + \boldsymbol{u}^* \boldsymbol{\cdot} \boldsymbol{\nabla} s^* + w^* \frac{\mathrm{d} \bar{s}}{\mathrm{d} z}= 0, \end{gather}

(2.26)\begin{gather}\frac{p^*}{\bar{p}} = \frac{\rho^*}{\bar{\rho}} + \frac{T^*}{\bar{T}}, \end{gather}

(2.27)\begin{gather}s^* = \frac{T^*}{\bar{T}} - \frac{(\gamma - 1)}{\gamma} \frac{p^*}{\bar{p}}, \end{gather}

where $\boldsymbol{u}^*=(u^*, v^*, w^*)$. The reference state entropy gradient, of ${O}(\varepsilon )$, is

(2.28)\begin{equation} \varepsilon \frac{\mathrm{d} \bar{s}}{\mathrm{d} z} = \frac{1}{\gamma}\frac{\mathrm{d} }{\mathrm{d} z} \left( \ln \bar{p} \bar{\rho}^{-\gamma} \right). \end{equation}

Equations (2.22)–(2.27) constitute the equations of anelastic MHD in dimensionless form.

A further simplification of the anelastic momentum equation may be realised by subsuming the reference state density into the pressure gradient term and also by using (2.26), (2.27) to eliminate the density perturbation in the buoyancy term in favour of $s^*$ and $p^*$. Thus we obtain

(2.29)\begin{equation} -\frac{1}{\bar{\rho}} \boldsymbol{\nabla} p^* + \frac{\lambda \rho^*}{\bar{\rho}} \hat{\boldsymbol{e}}_z ={-} \boldsymbol{\nabla} \left( \frac{p^*}{\bar{\rho}} \right) - \lambda s^* \hat{\boldsymbol{e}}_z - \frac{p^*}{\bar{\rho}} \left[ \frac{1}{\bar{\rho}} \frac{\mathrm{d} \bar{\rho}}{\mathrm{d} z} - \frac{\lambda \bar{\rho} }{\gamma \bar{p}} \right] \hat{\boldsymbol{e}}_z. \end{equation}

For a near-adiabatic reference state,

(2.30)\begin{equation} \frac{\mathrm{d} }{\mathrm{d} z} \left( \frac{1}{\gamma}\ln \bar{p}\bar{\rho}^{-\gamma} \right) = \frac{1}{\gamma \bar{p}}\frac{\mathrm{d} \bar{p}}{\mathrm{d} z} - \frac{1}{\bar{\rho}}\frac{\mathrm{d} \bar{\rho}}{\mathrm{d} z} = {O}(\varepsilon) . \end{equation}

On using the equation of hydrostatic balance (2.20), it therefore follows that

(2.31)\begin{equation} \frac{1}{\bar{\rho}}\frac{\mathrm{d} \bar{\rho}}{\mathrm{d} z} = \frac{ \lambda \bar{\rho} }{\gamma \bar{p}} + {O}(\varepsilon) . \end{equation}

The term in the square brackets on the right-hand side of (2.29) is thus formally smaller than the other two terms; the momentum equation (2.23) therefore becomes

(2.32)\begin{equation} \frac{\partial \boldsymbol{u}^*}{\partial t^*} + \boldsymbol{u}^*\boldsymbol{\cdot} \boldsymbol{\nabla} \boldsymbol{u}^* ={-} \boldsymbol{\nabla} \left( \frac{p^*}{\bar{\rho}} \right) - \lambda s^* \hat{\boldsymbol{e}}_z + \frac{\tilde{M}_A^2 }{\bar{\rho}} \left( \boldsymbol{\nabla}\times \boldsymbol{B}^* \right) \times \boldsymbol{B}^* . \end{equation}

This simplification of the momentum equation was demonstrated by Lantz (Reference Lantz1992) and Braginsky & Roberts (Reference Braginsky and Roberts1995). It is important to note that it requires no additional approximation, but follows immediately from the original assumption of a near-adiabatic reference state; i.e. (2.23) and (2.32) are asymptotically equivalent leading-order expressions.

In the course of the following analysis it will sometimes be advantageous to work with the governing equations in their dimensional form, which we therefore state here for completeness. Equations (2.22), (2.24), (2.25) and (2.26) are unchanged. In dimensional form, the two versions of the momentum equation (2.23) and (2.32), and the thermodynamic relation (2.27) read as

(2.33)\begin{gather} \bar{\rho} \left( \frac{\partial \boldsymbol{u}^*}{\partial t^*} + \boldsymbol{u}^*\boldsymbol{\cdot} \boldsymbol{\nabla} \boldsymbol{u}^* \right) ={-} \boldsymbol{\nabla} p^* + \rho^* g \hat{\boldsymbol{e}}_z + \frac{1}{\mu_0} \left( \boldsymbol{\nabla}\times \boldsymbol{B}^* \right) \times \boldsymbol{B}^*, \end{gather}

(2.34)\begin{gather}\frac{\partial \boldsymbol{u}^*}{\partial t^*} + \boldsymbol{u}^*\boldsymbol{\cdot} \boldsymbol{\nabla} \boldsymbol{u}^* ={-} \boldsymbol{\nabla} \left( \frac{p^*}{\bar{\rho}} \right) - \frac{s^*}{c_p} g \hat{\boldsymbol{e}}_z + \frac{1}{\mu_0 \bar{\rho}} \left( \boldsymbol{\nabla}\times \boldsymbol{B}^* \right) \times \boldsymbol{B}^* . \end{gather}

(2.35)\begin{gather}s^* = c_p \frac{T^*}{\bar{T}} - (c_p - c_v) \frac{p^*}{\bar{p}} . \end{gather}

The dimensional equations governing the reference state are

(2.36a–c)\begin{equation} \frac{\mathrm{d} \bar{p}}{\mathrm{d} z} = \bar{\rho}g, \quad \bar{p} = \mathcal{R} \bar{\rho} \bar{T},\quad \varepsilon \frac{\mathrm{d} \bar{s}}{\mathrm{d} z} = c_v \frac{\mathrm{d} }{\mathrm{d} z} \ln{\bar{p} \bar{\rho}^{-\gamma}}. \end{equation}

2.3. The subtle ordering for magnetic buoyancy

As can be seen from the formulation above, the magnetic field is incorporated into the anelastic approximation in a straightforward way, with the only requirement being that the field is sufficiently weak. On the other hand, including the effects of magnetic buoyancy in the Boussinesq approximation is a subtle procedure, which, inter alia, involves consideration of the length scale characteristic of magnetic buoyancy perturbations (Spiegel & Weiss Reference Spiegel and Weiss1982; Corfield Reference Corfield1984; Bowker et al. Reference Bowker, Hughes and Kersalé2014). Specifically, the requirement that the length scale of motions in the direction of the imposed horizontal magnetic field be long in comparison with the transverse scale has to be built into the approximation. In the anelastic system, however, no such special measures are taken to ensure that magnetic buoyancy is included consistently. The question thus arises as to whether the anelastic equations do indeed retain the effects of magnetic buoyancy, and, if so, why are no special measures required governing the length scale of the perturbation? Here, we seek to elucidate this matter by examining the fundamental ordering of the physical quantities necessary for magnetic buoyancy instability, thus clarifying the relationship between the descriptions of magnetic buoyancy in the full (compressible) and reduced (anelastic, Boussinesq) systems. By striking the appropriate balances between terms in the governing equations, we derive the general scalings necessary to account for the effects of magnetic buoyancy in the reduced equations. The resulting scalings allow us to identify, and distinguish between, different reduced systems (or regimes) of the full compressible MHD equations.

As above, we consider an atmosphere of depth $d$ and pressure scale height $H_p = p_r/(\rho _r g)$. In the basic state, we assume an imposed horizontal magnetic field, stratified in the vertical direction with scale height $H_B$. When the system is perturbed, and motions ensue, the presence of the imposed horizontal field introduces a length scale in the direction of the field that may be distinct from both $d$ and $H_B$: we denote this length scale by $L_B$. For any variable $f$, we define $f_r$ and $\delta f$ to be, respectively, representative values of $f$ in equilibrium and of the magnitude of fluctuations of $f$. For vector fields, it is necessary to distinguish between components aligned with, and those perpendicular to, the imposed magnetic field. We denote the magnitudes of the components of the fluctuations parallel and perpendicular to the imposed field by subscripts $\parallel$ and $\perp$, respectively.

By considering the magnitudes of the fluctuating quantities in the dimensional compressible equations (2.1)–(2.5), we will establish the relation that must be obeyed between the various length scales of the problem if magnetic buoyancy is to be of significance. We will pursue, in a slightly more general fashion, the line of argument expounded by Bowker et al. (Reference Bowker, Hughes and Kersalé2014). First, from the momentum equation (2.2), balancing the inertia term with those of total pressure fluctuations, buoyancy and magnetic tension gives

(2.37)\begin{equation} \rho_r \delta u_\perp^2 \sim \delta \varPi \sim \delta\rho g d \sim \frac{B_r \delta B_\perp}{\mu_0} \left( \frac{d}{L_B} \right) . \end{equation}

With our focus on buoyancy-driven instabilities, the time scale is determined by the balance of vertical acceleration and buoyancy, thus giving the scaling

(2.38)\begin{equation} \frac{\partial }{\partial t}\sim \frac{\delta u_\perp}{d}. \end{equation}

From the balance between inertia and buoyancy in (2.37), we may express the kinetic energy of the transverse flow in terms of the density variation as

(2.39)\begin{equation} \delta u_\perp^2 \sim \frac{\delta \rho}{\rho_r} gd \sim \frac{\delta \rho}{\rho_r} \left( \frac{d}{H_p} \right) c_s^2, \end{equation}

where $c_s^2 = p_r/\rho _r$. The scaling for total pressure fluctuations, $\delta \varPi = \delta p + \delta p_m$, results from balancing the gradient of total pressure fluctuations with the inertia terms, and hence with the buoyancy, thus,

(2.40)\begin{equation} \delta \varPi \sim \rho_r \delta u_\perp^2 \sim \frac{\delta \rho}{\rho_r} \frac{d}{H_p} p_r. \end{equation}

Balancing advection and stretching terms in the parallel and perpendicular components of the induction equation (2.3) results in the following relations:

(2.41a–c)\begin{equation} \delta B_\parallel{\sim} \left( \frac{d}{H_B} \right) B_r, \quad \delta B_\perp{\sim} \left( \frac{d}{L_B} \right) B_r, \quad \delta u_\parallel{\sim} \left( \frac{L_B}{H_B} \right) \delta u_\perp . \end{equation}

From the balance between buoyancy and magnetic tension in (2.37), together with (2.41b), the density variation may be expressed in terms of the representative field strength as

(2.42)\begin{equation} \frac{\delta \rho}{\rho_r} \sim \left(\frac{B_r^2}{\mu_0 p_r} \right) \left( \frac{H_p}{d} \right) \left( \frac{d}{L_B} \right)^2 . \end{equation}

It should be noted that (2.42) expresses the contribution to the density variation arising from the magnetic field, which is of particular interest here. There will of course also be a contribution arising from entropy variations, which is present even in the absence of field. From (2.41a) and (2.42), the magnetic pressure variation may be expressed as

(2.43)\begin{equation} \delta p_m \sim \frac{B_r \delta B_\parallel}{\mu_0} \sim \frac{\delta \rho}{\rho_r} \left( \frac{d}{H_p} \right) \left( \frac{d}{H_B} \right) \left( \frac{L_B}{d} \right)^2 p_r. \end{equation}

The above orderings are quite general, insofar as they arise simply from balancing terms in the governing equations, but with no assumptions having been made about the magnitude of the density fluctuation $\delta \rho$ relative to the background $\rho _r$, nor of the relative magnitudes of the length scales $d$, $H_p$, $H_B$ and $L_B$.

Note that both the anelastic and Boussinesq approximations assume that the size of thermodynamic fluctuations is small compared with their background values: ${\delta \rho /\rho _r \ll 1}$. The anelastic approximation is valid for stratified atmospheres whose vertical extent can span many pressure scale heights: $d \geqslant H_p$. The Boussinesq approximation is more restrictive as it applies only to weakly stratified atmospheres – ones whose depth is much smaller than the pressure scale height: $d \ll H_p$.

For magnetic pressure fluctuations to be influential – a necessary condition for magnetic buoyancy instability – they must be comparable in magnitude to density fluctuations: $\delta p_m / p_r \sim \delta \rho /\rho _r$. From (2.43), we deduce that this requirement imposes the following important relation between length scales:

(2.44)\begin{equation} L_B^2 \sim H_B H_p . \end{equation}

To establish the consistency of the above argument, it is instructive to substitute for $L_B^2$ from (2.44) into (2.42), which gives the density variation due to the magnetic field as

(2.45)\begin{equation} \frac{\delta \rho}{\rho_r} \sim \left( \frac{d}{H_B}\right) \left(\frac{B_r^2}{\mu_0 p_r} \right) . \end{equation}

Thus, for a fixed $H_B$, the field will not influence the density if $B_r$ is sufficiently weak; conversely, the influence of the field is accentuated by very small values of $H_B$.

From the relation (2.44), we are able to identify four distinct regimes in which magnetic buoyancy is faithfully described – two in the anelastic approximation and two in the Boussinesq approximation – together with a further Boussinesq regime in which the effects of magnetic buoyancy are explicitly excluded. These distinct regimes, or reduced systems, are explained in more detail in the following subsection. Furthermore, each reduced system has its own distinctive set of governing equations, which are provided in the Appendix.

2.4. Distinguished regimes of the anelastic and Boussinesq approximations

In order to categorise the five possible regimes, we define the following ordering parameters:

(2.46a–c)\begin{equation} \varepsilon_1 = \frac{d}{H_p}, \quad \varepsilon_2 = \frac{\delta \rho}{\rho_r}, \quad \varepsilon_B = \frac{d}{H_B}. \end{equation}

The definitions of $\varepsilon _1$ and $\varepsilon _2$ are in accordance with the notation used in Corfield (Reference Corfield1984) and Bowker et al. (Reference Bowker, Hughes and Kersalé2014). Note also that $\varepsilon _2$ is equivalent to the parameter $\varepsilon$ used to derive the anelastic equations in § 2.2. For all the various approximations, thermodynamic fluctuations are considered small: $\varepsilon _2 \ll 1$.

It is important to consider how the entropy gradient comes into this ordering. The background entropy gradient must be of the size of thermodynamic fluctuations (at most); otherwise, the assumption that $\delta \rho /\rho _r \ll 1$ would not be justified. Indeed, this is an explicit requirement of the anelastic approximation. Likewise, the same ordering of the entropy gradient ($\mathrm {d} s/\mathrm {d}z \sim \varepsilon _2$) pertains to the Boussinesq approximation (Spiegel & Weiss Reference Spiegel and Weiss1982).

2.4.1. Standard anelastic approximation

The derivation of the anelastic equations in § 2.2 involved no special consideration of the magnetic scale height $H_B$, de facto treating it as on a par with the layer depth: ${d \sim H_B}$. Likewise, making no distinction between the parallel and perpendicular directions amounts to an implicit assumption that $L_B \sim d$. With $d/H_p = {O}(1)$, it follows that all the length scales stand on an equal footing, i.e.

(2.47)\begin{equation} d \sim H_p \sim H_B \sim L_B. \end{equation}

Crucially, such an ordering satisfies (2.44), the essential length scale relation for magnetic buoyancy. Thus, possibly fortuitously, no special measures need to be taken to incorporate magnetic buoyancy into the standard anelastic approximation. Total pressure variations (2.40) and magnetic pressure variations (2.43) are of the same order of magnitude, with

(2.48)\begin{equation} \frac{\delta \varPi}{p_r} \sim \frac{\delta p_m}{p_r} = {O}(\varepsilon_2). \end{equation}

With the ordering of length scales (2.47), the relation (2.42) requires that the magnetic field strength, expressed through the parameter $M_A^2 = B_r^2/\mu _0 p_r$, is ${O}(\varepsilon _2)$.

2.4.2. Weak field-gradient anelastic approximation

For the case of weak magnetic stratification, i.e. $H_B \gg d$, it turns out that there is another possible regime within the anelastic framework. In this case, the parallel length scale $L_B$ is necessarily larger than $d$ in order to satisfy the crucial relation for magnetic buoyancy (2.44). More precisely,

(2.49a–c)\begin{equation} d \sim H_p, \quad d \ll H_B, \quad L_B \sim \sqrt{ d H_B}. \end{equation}

As for the ordering (2.47), variations in total and magnetic pressure obey (2.48). However, with a weak field gradient ($H_B \gg d$), the magnetic field needs to be stronger than in § 2.4.1 in order to compensate for the weak gradient; thus here, again using (2.42), $M_A^2 = {O}(\varepsilon _2/\varepsilon _B)$, with the requirement that $\varepsilon _2/\varepsilon _B \ll 1$, so that the Alfvén speed remains much slower that the sound speed.

2.4.3. Standard magneto-Boussinesq approximation

The length scales in the standard magneto-Boussinesq equations of Spiegel & Weiss (Reference Spiegel and Weiss1982) obey the ordering

(2.50a,b)\begin{equation} d \ll H_p, \quad H_B \sim L_B \sim H_p. \end{equation}

The magnitudes of total and magnetic pressure fluctuations are given by

(2.51a,b)\begin{equation} \frac{\delta \varPi}{p_r} = {O}(\varepsilon_1 \varepsilon_2), \quad \frac{\delta p_m}{p_r} = {O}(\varepsilon_2). \end{equation}

Thus,

(2.52)\begin{equation} \frac{\delta p}{p_r} ={-} \frac{\delta p_m}{p_r} + {O}(\varepsilon_1 \varepsilon_2), \end{equation}

thereby guaranteeing that magnetic pressure variations enter the equation of state to produce density variations. The field strength is given by $M_A^2 = {O}(\varepsilon _2/\varepsilon _1)$, subject to $\varepsilon _2 \ll \varepsilon _1$, so that the Alfvén speed is much smaller than the sound speed.

2.4.4. Strong field-gradient magneto-Boussinesq approximation

The magneto-Boussinesq orderings, which lead to (2.51a,b) and (2.52), can also be obeyed if the field gradient is much stronger than the pressure gradient (Bowker et al. Reference Bowker, Hughes and Kersalé2014). The precise ordering required is

(2.53a–c)\begin{equation} d \ll H_p, \quad H_B \sim d, \quad L_B \sim \sqrt{d H_p}. \end{equation}

Here, the field is weaker than in § 2.4.3, with $M_A^2 = {O}(\varepsilon _2)$.

2.4.5. Boussinesq magnetoconvection

Under the Boussinesq orderings, it is also possible to include the effects of magnetic tension, but neglect the dynamical influence of magnetic pressure. This gives the well-studied system of Boussinesq magnetoconvection (e.g. Weiss & Proctor Reference Weiss and Proctor2014), in which the length scales obey the following ordering:

(2.54a–c)\begin{equation} d \ll H_p, \quad L_B \sim d, \quad H_B \gtrsim d. \end{equation}

As expected, the necessary scaling for the inclusion of magnetic buoyancy, (2.44), is not now satisfied. The key feature of the Boussinesq magnetoconvection ordering is that the variations of both total and magnetic pressure are smaller than density variations (${\delta \rho /\rho _r = {O}(\varepsilon _2)}$), with

(2.55)\begin{equation} \frac{\delta \varPi}{p_r} \sim \frac{\delta p_m}{p_r} = {O}(\varepsilon_1 \varepsilon_2). \end{equation}

It thus follows that thermodynamic pressure fluctuations $\delta p/p_r$ are, at most, ${O}(\varepsilon _1 \varepsilon _2)$; hence, to leading order, density variations arise on account only of temperature variations. In this regime the magnetic field is also very weak, with $M_A^2 = {O} (\varepsilon _1 \varepsilon _2)$.

2.4.6. Summarising the five systems

We have described how, within each of the anelastic and magneto-Boussinesq approximations, there are two distinct orderings of the length scales of the problem that allow for the consistent inclusion of the effects of magnetic buoyancy. In addition to that of the standard anelastic system, there is a further distinct ordering with a weaker field gradient, a stronger field and a longer characteristic length scale of the perturbations. Similarly, in addition to the standard magneto-Boussinesq system, there is a distinct ordering with a stronger field gradient, a weaker field and a shorter characteristic length scale of the perturbations. Within the Boussinesq approximation, there is also the ordering that leads to the system of Boussinesq magnetoconvection, in which magnetic buoyancy is excluded.

It is important to note also that the relative magnitudes of $d$, $H_B$ and $L_B$ determine the relative sizes of the perpendicular and parallel components of the velocity and magnetic fluctuations. Naturally, the ordering of length scales also has implications for the scalings of spatial derivatives parallel and perpendicular to the field: $\boldsymbol{\nabla } = \boldsymbol{\nabla }_\perp + (d/L_B) \boldsymbol{\nabla }_\parallel$. The scalings of all quantities for each regime are summarised in table 1. On applying the orderings in table 1 to the governing compressible equations, we obtain five distinct reduced systems; these are described in detail in the Appendix. As a final point, we should note why there are three reductions possible under the Boussinesq approximation, but only two under the anelastic approximation. Within the assumptions of the Boussinesq approximation, it is possible to obtain an asymptotically consistent set of equations that includes the effects of magnetic tension but not of magnetic buoyancy – namely, the equations of Boussinesq magnetoconvection. However, within the anelastic approximation, such a reduction is not possible; the influence of the magnetic field is felt both via tension and buoyancy.

Table 1.

Summary of orderings in different regimes.

3.1. Compressible MHD

In the absence of fluid motion, the (dimensionless) compressible equations admit a steady $z$-dependent basic-state solution. The basic-state variables (denoted by subscript ‘$b$’) satisfy the well-known relations for a magnetohydrostatic perfect gas,

(3.1a,b)\begin{equation} \frac{\mathrm{d} p_b}{\mathrm{d} z} = \lambda \rho_b - \frac{1}{2} M_A^2 \frac{\mathrm{d} B_b^2}{\mathrm{d} z}, \quad p_b= \rho_b T_b . \end{equation}

Elimination of $p_b$ between the two equations in (3.1a,b) results in the ordinary differential equation (ODE) for $\rho _b$,

(3.2)\begin{equation} \frac{\mathrm{d} \rho_b}{\mathrm{d} z} + \rho_b \varGamma (z) ={-} \frac{1}{2} M_A^2 \frac{1}{T_b} \frac{\mathrm{d} B_b^2}{\mathrm{d} z}, \quad \text{where} \ \varGamma (z) = \frac{\mathrm{d} \ln T_b}{\mathrm{d} z} - \frac{\lambda}{T_b}. \end{equation}

On applying the boundary condition $\rho _b(z=0) = 1$, (3.2) can be integrated to give

(3.3)\begin{equation} \rho_b(z) = \exp \left[ -\int _0^z \varGamma (z) \, \mathrm{d} z \right] \left\lbrace 1 - \frac{1}{2} M_A^2 \int _0^z \frac{1}{T_b} \frac{\mathrm{d} B_b^2}{\mathrm{d} z} \exp \left[ \int _0^z \varGamma (z) \, \mathrm{d} z \right] \, \mathrm{d} z \right\rbrace. \end{equation}

We consider perturbations to the basic state, expressing velocity, magnetic field and thermodynamic variables in the perturbed state as $\boldsymbol{\delta u}$, $\boldsymbol{B}_b + \boldsymbol{\delta b}$, $p_b + \delta p$, etc., respectively. On assuming that the perturbations are small, and hence that nonlinear terms may be neglected, we can express all disturbances in the following form:

(3.4)\begin{equation} \delta p = \hat{p}(z) \exp({\rm i} \omega t + {\rm i}k_x x +{\rm i} k_y y), \quad \text{etc.}, \end{equation}

where $\omega$ is the (complex) oscillation frequency and $k_x$ and $k_y$ are the wavenumbers in the $x$ and $y$ horizontal directions. Substituting into the governing equations (2.1)–(2.4) and retaining only the lowest-order terms in the perturbations leads to the following linear system for the perturbations:

(3.5)\begin{gather} {\rm i} \omega \hat{R} ={-} {\rm i}k_x \hat{u} - {\rm i}k_y \hat{v} - \left( \frac{\mathrm{d} }{\mathrm{d} z} + \frac{1}{H_\rho(z)} \right) \hat{w} , \end{gather}

(3.6)\begin{gather}{\rm i} \omega \hat{u} ={-} {\rm i}k_x \hat{P} + M_A^2 \frac{B_b(z)^2}{\rho_b(z)} \left( {\rm i}k_y \hat{F}_x - {\rm i}k_x \hat{F}_y \right) , \end{gather}

(3.7)\begin{gather}{\rm i} \omega \hat{v} ={-} {\rm i}k_y \hat{P} + M_A^2 \frac{B_b(z)^2}{\rho_b(z) H_B(z)} \hat{F}_z , \end{gather}

(3.8)\begin{gather}{\rm i} \omega \hat{w} ={-} \left( \frac{\mathrm{d} }{\mathrm{d} z} +\frac{1}{H_\rho(z)} \right)\hat{P} + \lambda \hat{R} + M_A^2 \frac{B_b(z)^2}{\rho_b(z)} \left\lbrace- \left( \frac{\mathrm{d} }{\mathrm{d} z} + \frac{2}{H_B(z)} \right) \hat{F}_y + {\rm i} k_y \hat{F}_z \right\rbrace , \end{gather}

(3.9)\begin{gather}{\rm i}\omega \hat{F}_x = {\rm i} k_y \hat{u} , \end{gather}

(3.10)\begin{gather}{\rm i}\omega \hat{F}_y ={-} {\rm i}k_x \hat{u} - \left( \frac{\mathrm{d} }{\mathrm{d} z} + \frac{1}{H_B(z)} \right) \hat{w} , \end{gather}

(3.11)\begin{gather}{\rm i}\omega \hat{F}_z = {\rm i}k_y \hat{w} , \end{gather}

(3.12)\begin{gather}{\rm i} \omega \hat{P} ={-} \frac{\gamma p_b(z)}{\rho_b(z)}\left\lbrace {\rm i}k_x \hat{u} + {\rm i}k_y \hat{v} + \left( \frac{\mathrm{d} }{\mathrm{d} z} + \frac{1}{\gamma H_p(z)} \right) \hat{w} \right\rbrace , \end{gather}

where $\hat {R} = \hat {\rho }/\rho _b$, $\hat {P} = \hat {p}/\rho _b$, $\hat {F}_x = \hat {b}_x/B_b$, $\hat {F}_y = \hat {b}_y/B_b$, $\hat {F}_z = \hat {b}_z/B_b$ and

(3.13a–c)\begin{equation} \frac{1}{H_\rho} = \frac{1}{\rho_b} \frac{\mathrm{d} \rho_b}{\mathrm{d} z}, \quad \frac{1}{H_p} = \frac{1}{p_b} \frac{\mathrm{d} p_b}{\mathrm{d} z}, \quad \frac{1}{H_B} = \frac{1}{B_b} \frac{\mathrm{d} B_b}{\mathrm{d} z}. \end{equation}

Equations (3.5)–(3.12), together with the boundary conditions $\hat {w} = 0$ at $z = 0, 1$, constitute an eigenvalue problem for the frequency $\omega$. In general, this set of equations requires a numerical solution, which we obtain using Chebyshev differentiation matrices (Trefethen Reference Trefethen2000).

3.2. Anelastic MHD

In the absence of motion, the (dimensional) anelastic equations admit a steady $z$-dependent basic-state solution (denoted by index ‘0’), which satisfies

(3.14)\begin{gather} \frac{\mathrm{d} }{\mathrm{d} z} \left( p_0 + \frac{\tilde{M}_A^2}{2} B_0^2 \right) = \lambda \rho_0 , \end{gather}

(3.15a,b)\begin{gather} \frac{p_0}{\bar{p}} = \frac{\rho_0}{\bar{\rho}} + \frac{T_0}{\bar{T}} , \quad s_0 = \frac{T_0}{\bar{T}} - \frac{(\gamma-1)}{\gamma} \frac{p_0}{\bar{p}} . \end{gather}

Note that, since we are considering an ideal system with no diffusivity, we have three equations for the four variables $\rho _0$, $p_0$, $T_0$, $s_0$. If the thermal and magnetic diffusivity were non-zero, the balance between thermal diffusion and Ohmic heating in the energy equation would provide a relation between $T_0$ and $B_0$. In the ideal case, however, we have the freedom to set $T_0 = 0$. In other words, the imposed field $B_0$ influences the density and pressure stratification (through (3.14)), but does not affect the (steady state) temperature distribution. The basic-state entropy gradient can thus be expressed solely in terms of $\bar {p}$ and $B_0$ as

(3.16)\begin{equation} \frac{\mathrm{d} s_0}{\mathrm{d} z} = \frac{(\gamma-1)}{\gamma} \tilde{M}_A^2 \frac{B_0}{\bar{p}} \frac{\mathrm{d} B_0}{\mathrm{d} z}. \end{equation}

We express velocity, magnetic field and thermodynamic variables in the perturbed state as $\boldsymbol{\delta u}^*$, $\boldsymbol{B}_0 + \boldsymbol{\delta b}^*$, $s_0 + \delta s^*$, etc., respectively. As for the full compressible equations, we now consider small perturbations to the basic state, allowing us to express all disturbances as

(3.17)\begin{equation} \delta s^* = \hat{s}^*(z) \exp({\rm i} \tilde{\omega} t + {\rm i}k_x x +{\rm i}k_y y), \quad \text{etc.}, \end{equation}

where the tilde notation on the eigenvalue reflects the different time scaling between the compressible and anelastic systems (cf. (3.4))). Substituting expressions (3.17) into the governing equations (2.22), (2.32), (2.24), (2.25), and retaining only the lowest-order terms in the perturbations, leads to the following set of (dimensionless) perturbation equations:

(3.18)\begin{gather} 0={\rm i} k_x \hat{u}^* + {\rm i} k_y \hat{v}^* + \left( \frac{\mathrm{d} }{\mathrm{d} z} + \frac{1}{ \bar{H}_\rho(z)} \right) \hat{w}^* , \end{gather}

(3.19)\begin{gather}{\rm i} \tilde{\omega} \hat{u}^* ={-} {\rm i} k_x \hat{P}^* + \tilde{M}_A^2 \frac{B_0(z)^2}{\bar{\rho}(z)} \left( {\rm i} k_y \hat{F}_x^* - {\rm i} k_x \hat{F}_y^* \right) , \end{gather}

(3.20)\begin{gather}{\rm i} \tilde{\omega} \hat{v}^* ={-} {\rm i} k_y \hat{P}^* + \tilde{M}_A^2 \frac{B_0(z)^2}{\bar{\rho}(z) H_B(z)} \hat{F}_z^* , \end{gather}

(3.21)\begin{gather}{\rm i} \tilde{\omega} \hat{w}^* ={-} \frac{\mathrm{d} \hat{P}^*}{\mathrm{d} z} - \lambda \hat{s}^* +\tilde{M}_A^2 \frac{B_0(z)^2}{\bar{\rho}(z)}\left\lbrace - \left( \frac{\mathrm{d} }{\mathrm{d} z} +\frac{2}{H_B(z)} \right) \hat{F}_y^* + {\rm i} k_y \hat{F}_z^* \right\rbrace , \end{gather}

(3.22)\begin{gather}{\rm i} \tilde{\omega} \hat{F}_x^* = {\rm i} k_y \hat{u}^* , \end{gather}

(3.23)\begin{gather}{\rm i} \tilde{\omega} \hat{F}_y^* = {\rm i} k_y \hat{v}^* + \left( \frac{1}{\bar{H}_\rho(z)} - \frac{1}{H_B(z)} \right) \hat{w}^* , \end{gather}

(3.24)\begin{gather}{\rm i} \tilde{\omega} \hat{F}_z^* = {\rm i} k_y \hat{w}^* , \end{gather}

(3.25)\begin{gather}{\rm i} \tilde{\omega} \hat{s}^* ={-} \hat{w}^* \left( \frac{\mathrm{d} \bar{s}}{\mathrm{d} z} + \frac{\mathrm{d} s_0}{\mathrm{d} z} \right) , \end{gather}

where $\hat {P}^* =\hat {p}^*/\bar {\rho }$, $\hat {F}_x^* = \hat {b}_x^*/B_0$, $\hat {F}_y^* = \hat {b}_y^*/B_0$, $\hat {F}_z^* = \hat {b}^*_z/B_0$ and

(3.26a,b)\begin{equation} \frac{1}{\bar{H}_\rho} = \frac{1}{\bar{\rho}} \frac{\mathrm{d} \bar{\rho}}{\mathrm{d} z}, \quad \frac{1}{H_B} = \frac{1}{B_0} \frac{\mathrm{d} B_0}{\mathrm{d} z} . \end{equation}

As for the compressible system, the general solution of the eigenvalue problem defined by (3.18)–(3.25) requires numerical treatment.

3.3. Comparing solutions to the compressible and anelastic systems

Before proceeding to the analysis of the solutions of the compressible and anelastic systems, it is important to explain the relations between the parameters and variables pertaining to the two systems. The two systems differ in terms of time and velocity scales; by design, the anelastic system governs the dynamics of slow motions evolving on a long time scale. Thus, comparing the velocities and eigenvalues of the compressible system with those of the anelastic system requires the following scaling: $| \boldsymbol{u}| = \varepsilon ^{1/2} |\boldsymbol{u}^*|$, $\omega = \varepsilon ^{1/2} \tilde {\omega }$. Since we have used the same scaling for the magnetic field ($B_r$) for the two systems, no re-scaling is needed when comparing magnetic field perturbations. Note, however, that $M_A^2= \varepsilon \tilde {M}_A^2$.

In the compressible system, the basic-state density and pressure correspond to the sum of the reference and basic states of the anelastic system; thus, $\rho _b = \bar {\rho } + \varepsilon \rho _0$ and $p_b = \bar {p} + \varepsilon p_0$. For the temperature, however, $T_b = \bar {T}$ since, in ideal MHD, we have the freedom to impose one of the equilibrium thermodynamic quantities. An important difference between the two systems is that, in the anelastic case, the background reference state is strictly adiabatic (i.e. the static solution profiles $\bar {\rho }$, $\bar {p}$, $\bar {T}$ are such that the departure from adiabaticity $\Delta\!\!\nabla$ is strictly zero); the asymptotically small entropy gradient enters only into the prognostic equations, but not the background state. The basic state entropy gradient in the compressible system is given by

(3.27)\begin{equation} \frac{\mathrm{d} s_b}{\mathrm{d} z} = \frac{\mathrm{d} }{\mathrm{d} z} \left( \frac{1}{\gamma} \ln \left( \frac{p_b}{\rho_b^{\gamma}} \right) \right) = \left( \frac{\mathrm{d} T_b}{\mathrm{d} z} - \frac{\left( \gamma - 1 \right)}{\gamma} \lambda \right) \frac{1}{T_b} + \frac{\left( \gamma - 1 \right)}{\gamma} M_A^2 \left( \frac{B_b}{p_b} \right) \frac{\mathrm{d} B_b}{\mathrm{d} z}; \end{equation}

this is related to the reference and basic-state entropy gradients of the anelastic system as follows:

(3.28)\begin{equation} \frac{\mathrm{d} s_b}{\mathrm{d} z} = \varepsilon \frac{\mathrm{d} \bar{s}}{\mathrm{d} z} + \varepsilon \frac{\mathrm{d} s_0}{\mathrm{d} z} + {O}(\varepsilon^2). \end{equation}

4.1. Compressible case

The equilibrium state is given by

(4.1a–d)\begin{equation} \rho_b = \rho_r \exp\left( \frac{z}{H} \right), \quad p_b = p_r \exp\left( \frac{z}{H} \right), \quad T_b = T_r, \quad B_b = B_r \exp\left( \frac{z}{2H} \right), \end{equation}

where the (magnetohydrostatic) scale height $H$ satisfies

(4.2)\begin{equation} \frac{1}{H} = \frac{1}{\rho_b} \frac{\mathrm{d} \rho_b}{\mathrm{d} z} = \frac{1}{p_b} \frac{\mathrm{d} p_b}{\mathrm{d} z} = \frac{2}{B_b} \frac{\mathrm{d} B_b}{\mathrm{d} z} = \frac{g}{c_s^2 + \frac{1}{2} c_A^2} . \end{equation}

The linearised equations of compressible MHD (the dimensional versions of (3.5)–(3.11)) can be combined to give the following single second-order ODE for $\hat {w}$:

(4.3)\begin{align} &\left( c_A^2 \omega^2 + \gamma c_s^2(\omega^2 - c_A^2 k_y^2) \right) \left( \omega^2 - c_A^2 k_y^2 \right) \left( \frac{\mathrm{d}^2 \hat{w}}{\mathrm{d} z ^2} + \frac{1}{H} \frac{\mathrm{d} \hat{w}}{\mathrm{d} z} \right)\nonumber\\ &\hskip1.8pc\quad + \left\lbrace \vphantom{\left( \gamma c_s^2 c_A^2 k_y^2 - g\left( g-\frac{\gamma c_s^2}{H} \right) \right)}\omega^6 - \left( \gamma c_s^2 k_H^2 + c_A^2 (k_x^2 +2k_y^2) \right) \omega^4 \right. \nonumber\\ &\hskip3.5pc\quad + \left( c_A^2 k_y^2 k_H^2 (2\gamma c_s^2 +c_A^2) -g\left( g - \frac{\gamma c_s^2}{H} \right)k_H^2 + \frac{g}{H}c_A^2 k_x^2 \right) \omega^2 \nonumber\\ &\hskip10pc\quad \left.- c_A^2 k_y^2 k_H^2 \left( \gamma c_s^2 c_A^2 k_y^2 - g\left( g-\frac{\gamma c_s^2}{H} \right) \right) \right\rbrace \hat{w} = 0, \end{align}

where $k_H^2 = k_x^2 + k_y^2$. Subject to the impermeable boundary conditions $\hat {w} = 0$ on $z = 0$ and $z=d$, the ODE (4.3) admits a solution $\hat {w} \propto \exp \left ( -z/2H \right ) \sin \left ( k_z z \right )$, with $k_z = n {\rm \pi}/d$ (where mode $n = 1$ is the most readily destabilised). Substituting the ansatz for $\hat w$ into (4.3) yields the dispersion relation

(4.4)\begin{gather} \omega^6 - \left( (\gamma c_s^2 \!+\! c_A^2)k'^2 \!+\! c_A^2 k_y^2 \right) \omega^4 \!+\! \left( c_A^2 k_y^2 k'^2 (2\gamma c_s^2 \!+\!c_A^2) \!-\!g\left( g - \frac{\gamma c_s^2}{H} \right)k_H^2 + \frac{g}{H}c_A^2 k_x^2 \right) \omega^2\nonumber\\ \hskip10pc\quad - c_A^2 k_y^2 \left( \gamma c_s^2 c_A^2 k_y^2 k'^2 - g\left( g-\frac{\gamma c_s^2}{H} \right) k_H^2 \right) = 0, \end{gather}

where $k'^2 = k_x^2 + k_y^2 + k_z^2 +(4H^2)^{-1}$. Expression (4.4), which is a cubic in $\omega ^2$, describes three different wave modes: the acoustic and gravity modes found in the non-magnetic case, modified owing to the presence of the magnetic field, together with the slow magneto-acoustic wave. The former two modes are always stable (for sub-adiabatic atmospheres), whilst the latter can be destabilised through magnetic buoyancy. Yu (Reference Yu1965) showed that stability requires that the coefficient of $\omega ^0$ in (4.4) be negative. For long wavelength modes ($k_y \to 0$), this results in the following stability condition:

(4.5)\begin{equation} \gamma - 1 > \frac{1}{2} \frac{c_A^2}{c_s^2}. \end{equation}

Note that, according to Newcomb's energy criterion (see (1.4)), for the particular atmosphere considered here, two-dimensional interchange modes ($k_y \equiv 0$) are never unstable since the gradient of $B/\rho$ is negative. For two-dimensional undular modes ($k_x \equiv 0$), the stability criterion is

(4.6)\begin{equation} \gamma - 1 > \frac{c_s^2 c_A^2 + c_A^4}{4 c_s^4 + 3 c_s^2 c_A^2}. \end{equation}

This is precisely the condition given by Parker (Reference Parker1966) (his equation (7) in the absence of cosmic ray pressure).

4.2. Anelastic case

The background state in the anelastic case is described by

(4.7a–d)\begin{equation} \bar{\rho} =\rho_r \exp\left( \frac{z}{\bar{H}} \right), \quad \bar{p} = p_r \exp\left( \frac{z}{\bar{H}} \right), \quad \bar{T} = T_r, \quad B_0 = \tilde{B}_r \exp\left( \frac{z}{2\bar{H}} \right) , \end{equation}

where the (hydrostatic) scale height $\bar {H}$ satisfies

(4.8)\begin{equation} \frac{1}{\bar{H}} = \frac{1}{ \bar{\rho}} \frac{\mathrm{d} \bar{\rho}}{\mathrm{d} z} = \frac{1}{\bar{p}} \frac{\mathrm{d} \bar{p}}{\mathrm{d} z} = \frac{2}{B_0} \frac{\mathrm{d} B_0}{\mathrm{d} z} = \frac{g}{c_s^2}. \end{equation}

For such a reference state, the departure from adiabaticity, defined by (2.13), is

(4.9)\begin{equation} \Delta\!\!{\nabla} ={-} \frac{d}{\bar{H}} \frac{(\gamma - 1)}{\gamma} . \end{equation}

For this particular atmosphere, the stratification is near adiabatic, i.e. $\varepsilon = \left | \Delta\!\!\nabla \right | \ll 1$, when $\gamma$ is close to 1. The (sub-adiabatic) reference state entropy gradient is given by

(4.10)\begin{equation} \frac{\mathrm{d} \bar{s}}{\mathrm{d} z} ={-} \frac{c_p}{d} . \end{equation}

An important consequence of the restriction $\gamma \cong 1$ is that the basic-state entropy gradient $\mathrm {d} s_0/\mathrm {d}z$, defined by (3.16), is ${O}(\varepsilon )$ smaller than the reference state entropy gradient $\mathrm {d} {\bar s}/\mathrm {d}z$, and hence does not enter the linearised entropy equation (3.25).

In a similar fashion to the compressible case, the linearised anelastic equations (the dimensional versions of (3.18)–(3.24)) can be combined to yield a single second-order ODE for $\hat {w}$, with a solution of the form $\hat {w} \propto \exp \left ( -z/2\bar {H} \right ) \sin \left ( k_z z \right )$. The resulting dispersion relation for the anelastic system is

(4.11)\begin{align} &- \tilde{\omega}^4 \overline{k'}^2 + \left( 2\tilde{c}_A^2 k_y^2 \overline{k'}^2 + \frac{g }{d} k_H^2 + \frac{\tilde{c}_A^2}{2\bar{H}^2}(2k_x^2 - k_H^2) \right) \tilde{\omega}^2 \nonumber\\ &\qquad\hskip10pc- \tilde{c}_A^2 k_y^2 \left( \tilde{c}_A^2 k_y^2 \overline{k'}^2 + \frac{g }{d} k_H^2 - \frac{\tilde{c}_A^2}{2\bar{H}^2} k_H^2 \right) =0, \end{align}

where $\overline {k'}^2 = k_x^2 +k_y^2 +k_z^2 + (4\bar {H}^2)^{-1}$. The dispersion relation (4.11), which is a quadratic in $\omega ^2$, thus supports two waves – the fast mode associated with the sound wave in the compressible system has been filtered out. The remaining two waves are the gravity wave (modified by the magnetic field), and the slow magneto-acoustic mode.

Stability requires that the $\omega ^0$ term is negative. Hence, for long wavelength modes ($k_y \to 0$), the stability criterion reduces to

(4.12)\begin{equation} \frac{d}{\bar{H}}\frac{\tilde{c}_A^2}{c_s^2} < 2. \end{equation}

This is the anelastic version of the condition given by Yu (Reference Yu1965) for the compressible case. In fact, the anelastic condition (4.12) can be recovered from the compressible condition (4.5) by setting $\gamma = 1+ \varepsilon (\bar {H}/d)$ and $c_A^2 = \varepsilon \tilde {c}_A^2$. Similarly, the stability criterion for two-dimensional undular modes (4.6) reduces in the anelastic case to

(4.13)\begin{equation} \frac{d}{\bar{H}}\frac{\tilde{c}_A^2}{c_s^2} < 4. \end{equation}

4.3. Comparison of the two systems

Although the derivation of the dispersion relations (4.4) and (4.11), individually, is performed most clearly in dimensional variables, it is more instructive, when comparing the two, to work with dimensionless quantities. On scaling length with the layer depth $d$, and time with $d/c_{s}$ (i.e. $\omega \to \left ( c_{s}/d \right ) \omega$, $k \to d^{-1} k$), the compressible dispersion relation (4.4) becomes

(4.14)\begin{align} &\omega^6 - \left( \gamma k'^2 + M_A^2 (k'^2 +k_y^2) \right) \omega^4\nonumber\\ &\quad + \left( M_A^2 k_y^2 k'^2 (2\gamma + M_A^2) + (\gamma-1) \lambda^2 k_H^2 + \frac{ \lambda M_A^2}{2H} \left( 2k_x^2 - \gamma k_H^2 \right) \right) \omega^2\nonumber\\ & \hskip10pc\quad - M_A^2 k_y^2 \left( \gamma M_A^2 k_y^2 k'^2 + (\gamma-1) \lambda^2 k_H^2 - \frac{\gamma \lambda M_A^2 }{2H} k_H^2 \right) = 0, \end{align}

with dimensionless scale height $H = \lambda ^{-1} \left ( 1+ \frac {1}{2}{M}_A^2 \right )$.

Similarly, the anelastic dispersion relation (4.11) in dimensionless form is

(4.15)\begin{align} &- \overline{k'}^2 \tilde{\omega}^4 + \left( 2\tilde{M}_A^2 k_y^2 \overline{k'}^2 + \lambda k_H^2 + \frac{\lambda \tilde{M}_A^2}{2 \bar{H}} \left( 2k_x^2 - k_H^2 \right) \right) \tilde{\omega}^2 \nonumber\\ &\hskip10pc\quad \ - \tilde{M}_A^2 k_y^2 \left( \tilde{M}_A^2 k_y^2 \overline{k'}^2 + \lambda k_H^2 - \frac{\lambda \tilde{M}_A^2}{2 \bar{H}} k_H^2 \right) = 0, \end{align}

where $\bar {H} = \lambda ^{-1}$. To compare directly the dispersion relation for the compressible system (4.14) with that for the anelastic system (4.15), the following additional rescaling is necessary: $\omega = \varepsilon ^{1/2} \tilde {\omega }$, $M_A = \varepsilon ^{1/2} \tilde {M}_A$. Furthermore, the restriction that the stratification of the atmosphere be close to adiabatic pins down the value of $\gamma$ (cf. (4.9)) as

(4.16)\begin{equation} \gamma = \frac{\lambda}{\lambda - \varepsilon} = 1+ \frac{\varepsilon}{\lambda} + {O}(\varepsilon^2). \end{equation}

On substitution for $\omega$, $M_A$ and $\gamma$, (4.14) may be expressed as the following asymptotic expression:

(4.17)\begin{align} - \overline{k'}^2 \tilde{\omega}^4 &+ \left( 2\tilde{M}_A^2 k_y^2 \overline{k'}^2 + \lambda k_H^2 + \frac{\lambda^2 \tilde{M}_A^2}{2 } \left( 2k_x^2 - k_H^2 \right) \right) \tilde{\omega}^2 - \tilde{M}_A^2 k_y^2 \left( \tilde{M}_A^2 k_y^2 \overline{k'}^2 + \lambda k_H^2 - \frac{\lambda^2 \tilde{M}_A^2}{2 } k_H^2 \right)\nonumber\\ &+ \varepsilon \left\lbrace \tilde{\omega}^6 - \left( \frac{\overline{k'}^2}{\lambda} + \tilde{M}_A^2 \left( k^2 +k_y^2 \right) \right) \tilde{\omega}^4 \right. \nonumber\\ &\qquad\left. + \left( \tilde{M}_A^2 k_y^2 \overline{k'}^2 \left( \frac{2}{\lambda} + \tilde{M}_A^2 \right) - \frac{\lambda \tilde{M}_A^2}{2 } \left( k_H^2 + \frac{\lambda \tilde{M}_A^2}{2} k_H^2 \right) \right) \tilde{\omega}^2 \right. \nonumber\\ &\qquad\quad \left. - \tilde{M}_A^2 k_y^2 \left( \tilde{M}_A^2 k_y^2 \frac{\overline{k'}^2}{\lambda} - \frac{ \lambda \tilde{M}_A^2 }{2} \left( k_H^2 - \frac{\lambda \tilde{M}_A^2}{2} k_x^2 \right) \right) \right\rbrace + {O}\left( \varepsilon^2 \right) =0, \end{align}

where $k^2 = k_x^2 + k_y^2 + k_z^2$. Setting $\varepsilon = 0$ in (4.17) recovers the anelastic dispersion relation (4.15). By inspecting the expression (4.17), we can assess under what conditions do the dispersion relations resulting from the compressible and anelastic systems produce significantly different results. In particular, we are interested in seeing if there are circumstances under which the asymptotic ordering can be broken, with terms that are formally ${O}\left ( \varepsilon \right )$ being ‘promoted’ to ${O}\left ( 1 \right )$.

For any of the ${O}\left ( \varepsilon \right )$ terms (in curly brackets) to come into play at ${O}\left ( 1 \right )$, $\tilde {M}_A^2$ would have to be as large as ${O}\left ( \varepsilon ^{-1} \right )$ – a value far beyond the regime of validity of the anelastic approximation; with this ordering, the Alfvén wave speed would be comparable to the speed of sound, and the magnetic field would be strong enough to upset the assumed hydrostatic balance at leading order. At first glance one might also intuit that some of the ${O}\left ( \varepsilon \right )$ terms could contribute to the dominant balance if $\lambda \sim \varepsilon$. However, the derivation of (4.17) involves expanding $\gamma$ according to (4.16) which is valid only for $\lambda \gg \varepsilon$. The validity of the expression (4.17) requires that $\tilde {M}_A^2$ be of order unity and $\varepsilon /\lambda \ll 1$ – under these restrictions it is not possible to upset the ${O}(\varepsilon ^0)$ balance.

Figure 1 compares (squared) frequencies of the slow magneto-acoustic mode for a range of wavenumbers $k_x$, $k_y$ and increasing values of $\tilde {M}_A^2$. Positive (negative) values correspond to stability (instability). As for all diffusionless magnetic buoyancy instabilities (see, for example, the energy principle analysis of Hughes & Cattaneo Reference Hughes and Cattaneo1987), instability is favoured for large $k_x$ and small $k_y$. As expected, qualitative and quantitative agreement between the two systems is very good when $\tilde {M}_A^2 = 10$. The agreement between the two systems worsens as $\tilde {M}_A^2$ is increased. When $\tilde {M}_A^2 = 50$, quantitative agreement is less good, but the overall features of the contour plot are still captured by the anelastic system. When $\tilde {M}_A^2$ is as large as $\varepsilon ^{-1}$, in this case $\tilde {M}_A^2 = 100$, the agreement between compressible and anelastic systems breaks down and there are large differences between the solutions of the two systems.

Figure 1.

Contours of the (squared) frequency $\tilde {\omega }^2$ of the slow magneto-acoustic mode for the compressible system, with $\varepsilon = 0.01$ (solid lines), and the anelastic system (dashed lines); $\lambda =1$, $k_z = {\rm \pi}$. In (a) $\tilde {M}_A^2=10$; in (b) $\tilde {M}_A^2=50$; in (c) $\tilde {M}_A^2=100$.

5.1. Isothermal atmosphere with linear magnetic stratification

Here, as in § 4, we again consider an isothermal atmosphere, but now with the magnetic field linearly stratified with depth. The equilibrium temperature and magnetic field distributions are given by

(5.1a,b)\begin{equation} T_b = \bar{T} = 1, \quad B_b = B_0 = 1 + \zeta z . \end{equation}

From (3.3), the basic-state density in the compressible system is then given by

(5.2)\begin{equation} \rho_b(z) = p_b(z) = \exp( \lambda z) \left\lbrace 1 + \frac{ \zeta M_A^2 }{ \lambda } \left[ \zeta z \exp(-\lambda z) - \left( 1 + \frac{\zeta}{\lambda} \right) \left( 1- \exp(-\lambda z) \right) \right] \right\rbrace . \end{equation}

In the anelastic system, the non-magnetised reference state is described by

(5.3a,b)\begin{equation} \bar{\rho}(z) = \bar{p}(z) = \exp( \lambda z), \quad \frac{\mathrm{d} \bar{s}}{\mathrm{d} z} ={-}1. \end{equation}

Recall that the departure from adiabaticity $\varepsilon$ for an isothermally stratified atmosphere is governed by (4.9). Thus, in the compressible system, a subadiabatic value of $\gamma$ is chosen according to $\gamma = \lambda /(\lambda - \varepsilon )$. In the anelastic system, we use the value of $\gamma$ that gives exactly adiabatic stratification, i.e. $\gamma _{ad} = 1$. As a consequence, the basic state entropy gradient does not enter the picture (i.e. $\mathrm {d} s_0/ \mathrm {d}z = 0$).

Figure 2 shows the growth rate of the magnetic buoyancy instability as a function of $\tilde {M}_A^2$ for various magnetic field gradients $\zeta$ and various degrees of departure of adiabaticity $\varepsilon$, for fixed $k_x = 4$, $k_y = 0.01$, $\lambda = 1$. Since our interest lies in comparison of the anelastic and compressible systems, rather than in an extensive investigation of the instability, we have concentrated on fixed values of the wavenumbers, but recognising the fact that the instability mechanism is favoured when $k_y \ll k_x$.

Figure 2.

Growth rates $(\tilde \sigma = -{\rm Im}(\tilde{\omega}))$ of magnetic buoyancy instability in compressible ($\varepsilon =10^{-3}\ \text{and} \varepsilon = 5 \times 10^{-3}$) and anelastic systems as a function of $\tilde {M}_A^2$, for various field gradients ($\zeta =5, 10, 20$). In (a) $\lambda =1$; in (b) $\lambda =5$. Stars denote the positions where $G_{max} = 0.7/\varepsilon$.

For a given value of $\zeta$, the agreement between the compressible and anelastic results improves as $\varepsilon$ is reduced. This is expected and reassuring. For a given $\varepsilon$ and $\tilde {M}_A^2$, the discrepancy between the two systems becomes larger as $\zeta$ is increased. In the case of the compressible system, the basic state can be substantially altered by the presence of the strong field gradient (large $\zeta$); i.e. magnetic stratification can provide a significant contribution to the background entropy gradient. The basic-state entropy gradient (3.27) can be thought of as a sum of atmospheric gradient $\beta ^A$ and magnetic gradient $\beta ^M$ contributions, $\mathrm {d} s_b/ \mathrm {d}z = \beta ^A + \beta ^M$. For this particular nearly adiabatic basic state, the atmospheric contribution is $\left | \beta ^A \right | = \varepsilon$, and the magnetic part is, to leading order, $\beta ^M = \varepsilon ^2 G/\lambda$, where $G = \tilde {M}_A^2 \zeta B_0/\bar {p}$. In the anelastic system, the corresponding part of the entropy gradient due to magnetic stratification does not enter the governing equations ($\mathrm {d} s_0/ \mathrm {d}z = 0$). This difference between the two systems does not result in a significant disagreement between them provided that $G = {O}(1)$; in this case, the background entropy gradient is, to leading order (${O}(\varepsilon )$), given by the atmospheric gradient $\beta ^A$, and the magnetic component $\beta ^M$ is ${O}(\varepsilon ^2)$. However, we can expect breakdown of the agreement between the two systems when the magnetic part of the entropy gradient – which does not enter the picture in the anelastic system – becomes large enough to be influential in the compressible system, i.e. when $\beta ^M = {O}(\varepsilon )$. This requires $G={O}(\varepsilon ^{-1})$. The maximum value of $G$ in the interval $0 \leqslant z \leqslant 1$ is

(5.4)\begin{equation} G_{max} = \max_{z \in [0,1]} G = \left\{\begin{array}{@{}ll} \dfrac{ \zeta^2 \tilde{M}_A^2}{\lambda \exp(1-\lambda/\zeta)}, & \text{for } \zeta \geqslant \lambda, \\ \zeta \tilde{M}_A^2 , & \text{for } \zeta < \lambda. \end{array}\right. \end{equation}

Thus we can expect breakdown of the agreement between the two systems when $G_{max}$ becomes as large as $1/\varepsilon$; i.e. the validity of the anelastic approximation can be broken by a large enough field gradient $\zeta = {O}(\varepsilon ^{-1/2})$, even if the field is weak $\tilde {M}_A^2 = {O}(1)$. The argument above is confirmed in figure 2, where we mark values of $\tilde {M}_A^2$ for which $G_{max} = 0.7/\varepsilon$.

5.2. Linear temperature and magnetic stratification

Now we consider the case where the equilibrium temperature and magnetic field distributions are given by

(5.5a,b)\begin{equation} T_b =\bar{T}= 1 + \theta z, \quad B_b = B_0 = 1 + \zeta z. \end{equation}

The basic-state density in the compressible system is (from (3.3))

(5.6)\begin{equation} \rho_b(z) = \left( 1+\theta z \right)^m \left\lbrace 1 - \zeta M_A^2 \frac{\left( \zeta + (m-1)\theta \right)\left[ (1+\theta z)^m - 1 \right] - m\zeta \theta z }{ m(m-1)\theta^2 (1+\theta z)^m }\right\rbrace, \end{equation}

where $m = (\lambda /\theta - 1)$ is the polytropic index. The basic-state pressure is $p_b = \rho _b T_b$. In the anelastic system, the non-magnetised reference state density and pressure take the form of polytropes,

(5.7a,b)\begin{equation} \bar{\rho} = \left( 1+\theta z \right)^{m}, \quad \bar{p} = \left( 1+\theta z \right)^{m+1}. \end{equation}

The departure from adiabaticity for such an atmospheric stratification is given by

(5.8a,b)\begin{equation} \varepsilon = \left| \Delta\!\!{\nabla}_r \right|, \quad \Delta\!\!{\nabla}_r = \frac{\theta (m+1)}{\gamma} \left( 1 - \frac{m \gamma}{m+1} \right) . \end{equation}

For a given $\gamma$, the polytropic index corresponding to adiabatic stratification is $m_{ad} = 1/(\gamma - 1)$. For the isothermal atmospheres discussed in §§ 4 and 5.1, the constraint of near adiabaticity required $\gamma$ to be close to unity; here, the value of $\gamma$ is unrestricted, so we adopt the standard value of $\gamma =5/3$. In the compressible system, with temperature stratification $\theta$, $\varepsilon$ is adjusted by selecting the appropriate polytropic index $m$, according to (5.8a,b). In the anelastic case we take the value of $m$ that gives exactly adiabatic stratification; with $\gamma =5/3$, $m_{ad} = 1.5$. The sub-adiabatic reference state entropy gradient is given by

(5.9)\begin{equation} \frac{\mathrm{d} \bar{s}}{\mathrm{d} z} ={-} \frac{1}{1+ \theta z}. \end{equation}

Note that, in contrast to the isothermal atmospheres, where the basic state entropy gradient does not enter the anelastic equations, here, it does come into play.

Figure 3 compares the growth rates of the magnetic buoyancy instability between the compressible and anelastic systems, again with $k_x = 4$, $k_y = 0.01$. The agreement between the two systems is better for larger stratification $\theta$.

Figure 3.

Growth rates $(\tilde \sigma = -{\rm Im}(\tilde{\omega}))$ of magnetic buoyancy instability in compressible ($\varepsilon =10^{-3}$ and $\varepsilon = 5 \times 10^{-3}$) and anelastic systems as a function of $\tilde {M}_A^2$ for the case of equilibria with linear temperature and magnetic field stratification ($\zeta =5, 10, 20$). In (a) $\theta =5$; in (b) $\theta =20$. Stars denote the positions where $H_{max} = 0.7/\varepsilon$.

The anelastic system is formally valid when the magnetic contribution to the entropy gradient $\mathrm {d} s_0/ \mathrm {d}z$ is order unity. We have seen above how, through a combination of large field strength and strong field gradient, the magnetic field can substantially affect the background entropy gradient. Here, it is also true that for sufficiently large $\tilde {M}_A^2$ and $\zeta$, the basic-state entropy gradient can become as large as $\mathrm {d} s_0/ \mathrm {d}z = {O}(\varepsilon ^{-1})$. In this case, the formal asymptotic ordering is broken, and the influence of magnetic field is strong enough to modify the assumed hydrostatic equilibrium. The maximum value of $\mathrm {d} s_0/ \mathrm {d}z$ in the interval $0 \leqslant z \leqslant 1$ is

(5.10)\begin{equation} H_{max} = \max_{z \in [0,1]} \frac{\mathrm{d} s_0}{\mathrm{d} z} = \left\{\begin{array}{@{}ll} \dfrac{(\gamma-1)}{\gamma} \dfrac{ \zeta (\zeta - \theta) \tilde{M}_A^2}{ (\gamma(1-\theta/\zeta) )^{m_{ad}+1}}, & \text{for } \zeta \geqslant (m_{ad}+1)\theta, \\ \dfrac{(\gamma-1)}{\gamma} \zeta \tilde{M}_A^2 , & \text{for } \zeta < (m_{ad}+1)\theta. \end{array}\right. \end{equation}

Thus we anticipate breakdown of the asymptotic ordering, and hence a breakdown of the agreement between the results of the compressible and anelastic systems, when $H_{max}$ becomes as large as $1/\varepsilon$. For illustration, in figure 3, we mark the values of $\tilde {M}_A^2$ for which $H_{max} = 0.7/\varepsilon$.