2.1 Hamilton's principle of stationary action

The process involved in constructing the action for fluid models has been well known since Lagrange (Reference Lagrange1789). Once the generalized coordinates $q_k(t)$ are chosen, where $k$ runs over all possible degrees of freedom, the action is determined via

(2.1)\begin{equation} S[q] = \int^{t_1}_{t_0}\textrm{d}t L\left(q,\dot{q},t\right), \end{equation}

with $L$ representing the Lagrangian. It must be noted that $S$ is a ‘functional’, i.e. its domain and range are functions and real numbers, respectively. Hamilton's principle states that that the equations of motion are the extrema of the action, i.e. we require $\delta S[q]/\delta q^k = 0$, where the functional derivative is defined as follows:

(2.2)\begin{equation} \delta S[q;\delta q]= \left.\frac{\textrm{d} S[q +\epsilon \delta q]}{\textrm{d} \epsilon}\right|_{\epsilon =0} =: \left\langle \frac{\delta S[q]}{\delta q^{i}},\delta q^{i}\right\rangle. \end{equation}

The continuum version is very similar to the discrete case since the discrete index $k$ is replaced by a continuous one, which we denote by $a$. The coordinate $q$ is a function of $a$ and $t$, and tracks the location of a fluid particle labelled by $a$. We also note the following important quantities which are used throughout the paper: the deformation matrix $\partial q^i/\partial a^j=:q^{\hspace {1pt}i}_{\,,\hspace {1pt}j}$ and the corresponding determinant, the Jacobian, $\mathcal {J}:= \det (q^{\hspace {1pt}i}_{\,,\hspace {1pt}j})$. The volume evolves in time via

(2.3)\begin{equation} \textrm{d}^3q=\mathcal{J} \,\textrm{d}^3a, \end{equation}

and the area is governed by

(2.4)\begin{equation} (\textrm{d}^2q)_i= \mathcal{J} a^{\hspace{1pt}j}_{\,,\hspace{1pt}i} (\textrm{d}^2a)_j, \end{equation}

where $\mathcal {J} a^{\hspace {1pt}j}_{\,,\hspace {1pt}i}$ is the transpose of the cofactor matrix of $q^{\hspace {1pt}j}_{\,,\hspace {1pt}i}$. The quantities and the relations introduced above can be used to generate a wide range of identities. One can find a detailed discussion of these, for example, in Serrin (Reference Serrin and Flügge1959), Morrison (Reference Morrison1998) and Bennett (Reference Bennett2006).

2.2 Two representations: the Lagrangian and the Eulerian points of view

The Lagrangian position $q$ evolves in time and is entirely characterized by its label $a$. But the fluid parcels are not solely determined by the position alone; they can also carry with them a certain density, entropy and magnetic field. As the fluid moves along its trajectory, these quantities are also transported along with it, and are consequently characterized only by the label $a$ as well. We will refer to these quantities as attributes. As the label $a$ is independent of time, these attributes serve as Lagrangian constants of motion. The subscript $0$ will be used to label the attributes, in order to distinguish them from their Eulerian counterparts.

Let us now consider the Eulerian picture. All Eulerian fields depend on the position $\boldsymbol {r}:=(x^1,x^2,x^3)$ and time $t$, which can both be measured in the laboratory. As a result, we shall refer to these fields as observables. Moving from the Eulerian to Lagrangian viewpoint and vice versa is accomplished with the Lagrange–Euler maps which we describe below in more detail.

The Eulerian velocity field $\boldsymbol {v}(\boldsymbol {r},t)$ is the velocity of the fluid element at a location $\boldsymbol {r}$ and time $t$. If we seek to preserve the equivalence of the Lagrangian and Eulerian pictures, this must also equal $\dot {\boldsymbol {q}}(\boldsymbol {a},t)$. As a result, it is evident that we require $\dot {\boldsymbol {q}}(\boldsymbol {a},t)=\boldsymbol {v}(\boldsymbol {r},t)$, where the dot indicates that the time derivative is obtained at fixed label $\boldsymbol {a}$. However, there is a discrepancy since the left-hand side is a function of $\boldsymbol {a}$ and $t$, while the right-hand side involves $\boldsymbol {r}$ and $t$. This conundrum is resolved by noting that the fluid element is at $\boldsymbol {r}$ in the Eulerian picture, and at $\boldsymbol {q}$ in the Lagrangian one. Hence, we note that $\boldsymbol {r}=\boldsymbol {q}(\boldsymbol {a},t)$, which implies that $\boldsymbol {a}=\boldsymbol {q}^{-1}(\boldsymbol {r},t)=:\boldsymbol {a}(\boldsymbol {r},t)$ upon inversion. As a result, our final Lagrange–Euler map for the velocity is

(2.5)\begin{equation} \boldsymbol{v}(\boldsymbol{r},t) =\left.\dot{\boldsymbol{q}}(\boldsymbol{a},t)\right|_{\boldsymbol{a}=\boldsymbol{a}(\boldsymbol{r},t)}. \end{equation}

Now we consider the attributes defined earlier, which we have noted are carried along by the fluid. The first attribute is the entropy of the fluid particle, which we shall label $s_0$. For ideal fluids, one expects the entropy to remain constant along the fluid trajectory. In other words, the Eulerian specific entropy $s(\boldsymbol {r},t)$ must also remain constant throughout, implying that $s=s_0$. Apart from entropy, the magnetic stream function $\psi$ for 2-D GVMHD (Andreussi et al. Reference Andreussi, Morrison and Pegoraro2013; Morrison et al. Reference Morrison, Lingam and Acevedo2014) also obeys this property.

Next, we can consider attributes which obey a conservation law similar to the density. The conservation law in this case is that of mass conservation. The attribute is denoted by $\rho _0(\boldsymbol {a})$ and the observable by $\rho (\boldsymbol {r},t)$. The statement of mass conservation in a given (infinitesimal) volume amounts to $\rho (\boldsymbol {r},t)\,\textrm {d}^3r=\rho _0(\boldsymbol {a})\,\textrm {d}^3a$. Using (2.3) we obtain $\rho _0=\rho \mathcal {J}$. As a result, we have found the Lagrange–Euler map for $\rho$. There exist other attribute-observable pairs in the literature, which also possess similar conservation laws, such as the entropy density.

In the case of magnetofluid models, it is often advantageous to introduce the magnetic field attribute $\boldsymbol {B}_0(\boldsymbol {a})$. In the case of ideal magnetofluid models, the conservation law of frozen-in magnetic flux is applicable. In algebraic terms, this amounts to $\boldsymbol {B}\boldsymbol {\cdot } \boldsymbol {d}^{\boldsymbol {2}}\boldsymbol {r}=\boldsymbol {B}_0\boldsymbol {\cdot } \boldsymbol {d}^{\boldsymbol {2}}{\boldsymbol {a}}$, and from (2.4) we obtain $\mathcal {J} B^i=q^{\hspace {1pt}i}_{\,,\hspace {1pt}j} \,B_0^j$.

In all of the above expressions, the picture is still incomplete since we need to remove the $\boldsymbol {a}$-dependence of the attributes. In a manner similar to that undertaken for the velocity, we evaluate the attributes at $\boldsymbol {a}=\boldsymbol {q}^{-1}(\boldsymbol {r},t)=:\boldsymbol {a}(\boldsymbol {r},t)$. This completes our prescription, and one can fully determine the observables once we are provided the attributes in conjunction with the Lagrangian coordinate $\boldsymbol {q}$.

We may also represent the Lagrange–Euler map in an integral form, which permits a more intuitive interpretation. We shall start with the assumption that the attribute-observable relations are found via appropriate conservation laws. We have stated before that one moves from the Lagrangian to the Eulerian picture by ‘plucking out’ the fluid element that happens to be at the Eulerian observation point $r$ at time $t$. Such a process is accomplished mathematically via the delta function $\delta (\boldsymbol {r}-\boldsymbol {q}(\boldsymbol {a},t))$. For instance, we see that the density can be treated as follows:

(2.6)\begin{align} \rho({\boldsymbol{r}},t)&=\int_D \textrm{d}^3a \, \rho_0(\boldsymbol{a}) \delta\left({\boldsymbol{r}}-{\boldsymbol{q}}\left(\boldsymbol{a},t\right)\right)\nonumber\\ &=\left. \frac{\rho_0}{\mathcal{J}}\right|_{\boldsymbol{a}=\boldsymbol{a}({\boldsymbol{r}},t)}. \end{align}

Further below, we will also use a new variable, the canonical momentum density $\boldsymbol {M}^c=(M^c_1,M^c_2,M^c_3)$, which is related to its Lagrangian counterpart via

(2.7)\begin{align} \boldsymbol{M}^c(\boldsymbol{r},t)&=\int_D \textrm{d}^3a \, \boldsymbol{\varPi}(\boldsymbol{a},t) \delta\left({\boldsymbol{r}}-{\boldsymbol{q}}(\boldsymbol{a},t)\right)\nonumber\\ &= \left.\frac{\boldsymbol{\varPi}(\boldsymbol{a},t)}{\mathcal{J}} \right|_{\boldsymbol{a}=\boldsymbol{a}(\boldsymbol{r},t)}. \end{align}

For ideal MHD, the canonical momentum density is $\boldsymbol {\varPi }(a,t)=(\varPi _1,\varPi _2,\varPi _3)=\rho _0 \dot {\boldsymbol {q}}$. It is worth noting that $\boldsymbol {\varPi }(\boldsymbol {a},t)$ can be found from the Lagrangian through $\boldsymbol {\varPi }(\boldsymbol {a},t) = {\delta L}/{\delta \dot {\boldsymbol {q}} }$ and does not necessarily equal $\rho _0 \dot {\boldsymbol {q}}$ in general. One can also construct such integral relations for the entropy and the magnetic field. We refer the reader to Morrison et al. (Reference Morrison, Lingam and Acevedo2014) for a more detailed discussion along these lines.

3.1 The general action

The domain of integration $D$ is chosen to be a subset of $\mathbb {R}^{3}$. Central to our formulation is the Lagrangian coordinate $\boldsymbol {q}\colon D\rightarrow D$, which we shall assume to be a well-behaved function with the required smoothness, invertibility, etc. Next we need to specify our set of observables, or alternatively our set of attributes. For our models, we work with $\mathfrak {E}=\{\boldsymbol {v},\rho ,\sigma ,\boldsymbol {B}\}$ where $\sigma = \rho s$ is the entropy density. Finally, we shall impose the Eulerian closure principle, which is necessary for our model to be ‘Eulerianizable’. Mathematically, this principle amounts to the action being fully expressible in terms of the Eulerian observables. Physically, the principle states that our theory must be solely describable in terms of physically meaningful quantities, the observables, and must also give rise to equations of motion in terms of these observables. As a result, we require our action to be given via

(3.1)\begin{equation} S[\boldsymbol{q}]:= \int_{T} \textrm{d}t\int_D \textrm{d}^3a\,\mathcal{L}\left(\boldsymbol{q},\dot{\boldsymbol{q}},\partial \boldsymbol{q}/\partial \boldsymbol{a}\right) =: \bar{S}[\mathfrak{E}]. \end{equation}

As per the Eulerian closure principle, this amounts to finding an action $\bar {S} = \int _\mathcal {T} \textrm {d}t \int _D \textrm {d}^3r\,\bar {\mathcal {L}}$ in terms of the Eulerian observables. The presence of the bar indicates that the action and the Lagrangian density are expressed solely in terms of the observables.

3.2 Constructing the gyroviscous action

The first step in the process involves the construction of the kinetic energy, which must also satisfy the closure principle. Using the analogy with particle mechanics, we know that it equals

(3.2)\begin{equation} S_{\mathrm{kin}}:= \int_\mathcal{T} \textrm{d}t\int_D \textrm{d}^3a\, \frac{1}{2} \rho_0 |\dot{\boldsymbol{q}}|^2 = \int_\mathcal{T} \textrm{d}t\int_D \textrm{d}^3r\, \frac{1}{2} \rho |\boldsymbol{v}|^2, \end{equation}

where the last equality is obtained by using relations outlined in § 2.2.

The internal energy per unit mass is a function of the entropy density and the density, and in Eulerian terms it can be represented by $U(\rho ,\sigma )$. Using the inverse Lagrange–Euler maps, we can construct the Lagrangian internal energy density accordingly,

(3.3)\begin{align} S_{\mathrm{int}}&:= \int_\mathcal{T} \textrm{d}t\int_D \textrm{d}^3a\, \rho_0 U\left(\frac{\rho_0}{\mathcal{J}},\frac{\sigma_0}{\mathcal{J}}\right)\nonumber\\ &= \int_\mathcal{T} \textrm{d}t\int_D \textrm{d}^3r\, \rho U\left(\rho,\sigma\right). \end{align}

The next step is the construction of the magnetic energy, and we use the same process outlined for the internal energy, viz. we determine the Eulerian term and obtain the Lagrangian version consequently through the Lagrange–Euler map,

(3.4)\begin{align} S_\textrm{mag}&:= \int_\mathcal{T} \int_D \textrm{d}^3r\, \frac{1}{2} |\boldsymbol{B}|^2,\nonumber\\ &= \int_\mathcal{T} \textrm{d}t\int_D \textrm{d}^3a\, \frac{1}{2 \mathcal{J}} q^{\hspace{1pt}i}_{\,,\hspace{1pt}j} \, q^{\hspace{1pt}i}_{\,,\hspace{1pt}k} \, B_0^j B_0^k. \end{align}

The magnetic energy is actually $|\boldsymbol {B}|^2/8{\rm \pi}$ in CGS units but we drop the factor of $4{\rm \pi}$ henceforth by scaling it away through the adoption of Alfvénic units.

Now we are ready to construct the most important term which will be responsible for the gyroviscosity. The gyroviscous term is taken to be linear in $\dot {\boldsymbol {q}}$ and is given by

(3.5)\begin{equation} S_\textrm{gyro} = \int_\mathcal{T} \textrm{d}t\int_D \textrm{d}^3a\, \dot{\boldsymbol{q}} \boldsymbol{\cdot} \boldsymbol{\varPi}^{\star} = \int_\mathcal{T} \textrm{d}t \int_D \textrm{d}^3r\, \boldsymbol{v} \boldsymbol{\cdot} \boldsymbol{M}^\star. \end{equation}

In other words, we operate under the premise that $\boldsymbol {\varPi }^{\star }$ is solely a functional of $\boldsymbol {q}$ and $t$. As the Eulerian perspective is inherently endowed with physical variables (e.g. density and magnetic field), we will focus on the Eulerian equivalent of $\boldsymbol {\varPi }^{\star }$; from the Eulerian closure principle we obtain the relation

(3.6)\begin{equation} \boldsymbol{M}^{\star} = \left. \frac{\boldsymbol{\varPi}^\star(\boldsymbol{a},t)}{\mathcal{J}} \right|_{\boldsymbol{a}=\boldsymbol{a}(\boldsymbol{r},t)}. \end{equation}

The complete action functional is now given by

(3.7)\begin{equation} S = S_\textrm{kin} - S_\textrm{int} - S_\textrm{mag} + S_\textrm{gyro}. \end{equation}

The action of (3.7) is general, but not the most general second-order (in $\boldsymbol {v}$) action that satisfies the Eulerian closure principle. For example, the term $S_\textrm {kin}$ could be generalized by replacing its integrand with $\rho _0G |\dot {\boldsymbol {q}}|^2/2|_{\boldsymbol {a}} = \rho G(\rho , \sigma , \boldsymbol {B})|\boldsymbol {v}|^2/2$ and the integrand of $S_\textrm {int}$ could be replaced by $\rho _0 U|_{\boldsymbol {a}}= \rho U(\rho , \sigma , \boldsymbol {B})$, a form that was shown in Morrison (Reference Morrison1982) to allow for anisotropic pressure. Here both $G$ and $U$ could be arbitrary functionals (including derivatives) of their arguments. Similarly the term $S_\textrm {mag}$ could be generalized.

The Eulerian canonical momentum density is defined via (2.7), which can be computed by finding the Lagrangian canonical momentum using $\boldsymbol {\varPi }(\boldsymbol {a},t) = {\delta L}/{\delta \dot {\boldsymbol {q}} }$ and Eulerianizing it. Upon doing so, we arrive at the so-called gyromap, a device introduced in Morrison et al. (Reference Morrison, Caldas and Tasso1984) as follows:

(3.8)\begin{equation} \boldsymbol{M}^c := \rho \boldsymbol{v} + \boldsymbol{M}^\star =\boldsymbol{M} + \boldsymbol{M}^\star. \end{equation}

The benefit of employing the gyromap and its natural origin will be discussed in § 6 and further explicated in appendix B.

So far we have only required $\boldsymbol {M}^\star$ to satisfy the closure principle, i.e. that it be expressible in terms of the subset $\{\rho ,\sigma ,\boldsymbol {B}\}\subset \mathfrak {E}$, including all possible Eulerian derivatives. Given that $\boldsymbol {M}^\star$ is a momentum density, arising perhaps from underlying gyration of particles, a natural assumption is that it has the magnetization form

(3.9)\begin{equation} \boldsymbol{M}^\star = \boldsymbol{\nabla} \times \boldsymbol{L}^\star, \end{equation}

i.e. we assume that $\boldsymbol {M}^\star$ is divergence-free. Since we are interested in a gyroviscosity due to gyromotion, this is a physically reasonable assumption. However, one could replace (3.9) by a Helmholtz decomposition for a more general collisionless viscosity. The present choice is also motivated in part by the realization in Morrison et al. (Reference Morrison, Caldas and Tasso1984) and Morrison et al. (Reference Morrison, Lingam and Acevedo2014) that this choice is consistent with existing 2-D gyroviscous models. Because $\boldsymbol {M}^\star$ has the units of momentum density, from which we see that the quantity $\boldsymbol {J}^\star \propto (q/m) \boldsymbol {M}^\star$ resembles a current density. If one assumes that the fluid ‘particles’ possess a finite magnetic moment, it follows that the fluid must have a finite magnetization. In other words, one may identify $\boldsymbol {J}^\star$ with the magnetization current density, which is divergence-free (Jackson Reference Jackson1998) and the current through an area depends on flux through a bounding curve. Are other choices possible and do any of them conserve angular momentum? Perhaps an even simpler way of envisioning the ansatz for $\boldsymbol {M}^\star$ is that it must emerge from the gyration of particles. In pictorial terms, this gyration is reminiscent of the effect generated by the curl of a vector field, which motivates our choice of $\boldsymbol {M}^\star$. Further grounds for assuming this particular expression are described in Morrison et al. (Reference Morrison, Lingam and Acevedo2014, § 5). With this ansatz, evidently $\boldsymbol {\nabla } \boldsymbol {\cdot } \boldsymbol {M}^c = \boldsymbol {\nabla } \boldsymbol {\cdot } \boldsymbol {M}$, since the second term vanishes. Note that the right-hand side of this expression appears in the continuity equation, and we see that one could also replace it by the left-hand side if we operate with $\boldsymbol {M}^\star = \boldsymbol {\nabla } \times \boldsymbol {L}^\star$. Furthermore, dimensional analysis permits the identification of $\boldsymbol {L}^\star$ with the angular momentum density.

As we have reduced the question of determining $\boldsymbol {L}^\star$, we must ask ourselves as to whether any further simplifications are feasible. Once again, we can resort to physical intuition to gain an idea of what $\boldsymbol {L}^\star$ might look like. Without further special assumptions about the fluid, e.g. it having some intrinsic or extrinsic direction, the vectorial character of $\boldsymbol {L}^\star$ must come from $\boldsymbol {B}$ or from the set of gradients of the observables; these and their cross products are the only vectors available. Thus, for example, a general form for $\boldsymbol {L}^\star$ could be composed of a linear combination of these vectors with coefficients dependent on $\rho , \sigma$ and $|\boldsymbol {B}|$. If we assume $\boldsymbol {L}^\star$ constitutes an internal angular momentum density of some kind associated with particle gyration, then it is reasonable to posit that it would tend to align with the magnetic field $\boldsymbol {B}$. Moreover, in the limit of a large magnetic field, the corresponding gyroradii would become small, owing to which the fluid particle may not be significantly affected by gradients on these scales. Combining the preceding arguments leads to the generic form

(3.10)\begin{equation} \boldsymbol{L}^\star= \mathcal{F}\left(\rho,\sigma,|\boldsymbol{B}|\right) \boldsymbol{B}. \end{equation}

In § 4.2 we will argue for further specification of the properties of (3.10).

With the choice of (3.10), the gyroviscous term of the action, expressed in terms of the observables is given by

(3.11)\begin{align} S_\textrm{gyro} &= \int_\mathcal{T} \textrm{d}t\int_D \textrm{d}^3r\, \boldsymbol{v} \boldsymbol{\cdot} \boldsymbol{\nabla} \times \left[\boldsymbol{B} \mathcal{F}\left(\rho,\sigma,|\boldsymbol{B}|\right) \right] \nonumber\\ &= \int_\mathcal{T} \textrm{d}t\int_D \textrm{d}^3r\, \mathcal{F} \boldsymbol{B}\boldsymbol{\cdot} \boldsymbol{\nabla} \times {\boldsymbol{v}}, \end{align}

where the second equality follows from integrating by parts and neglecting the boundary term. We shall use the latter operation consistently throughout the rest of the paper. Now that we have constructed the gyroviscous term, we note that it is still generic since there is considerable freedom in the choice of $\mathcal {F}$.

4.1 The equations of motion

The equations for the density, entropy density and the magnetic field can be determined via the attributes/observables relations defined through the appropriate conservation laws and the Lagrange–Euler maps. The entropy density and the density obey similar laws, given by

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

(4.2)\begin{gather}\frac{\partial \sigma}{\partial t} + \boldsymbol{\nabla} \boldsymbol{\cdot} \left( \sigma \boldsymbol{v}\right) = 0. \end{gather}

The equation governing the magnetic field is

(4.3)\begin{equation} \frac{\partial \boldsymbol{B}}{\partial t} + \boldsymbol{B} \left(\boldsymbol{\nabla} \boldsymbol{\cdot} \boldsymbol{v}\right) - \left(\boldsymbol{B} \boldsymbol{\cdot} \boldsymbol{\nabla}\right) \boldsymbol{v} + \left(\boldsymbol{v} \boldsymbol{\cdot} \boldsymbol{\nabla}\right) \boldsymbol{B} = 0, \end{equation}

which can be recast into the more familiar induction equation if $\boldsymbol {\nabla } \boldsymbol {\cdot } \boldsymbol {B} = 0$ is satisfied. If the constraint is obeyed, then we obtain

(4.4)\begin{equation} \frac{\partial \boldsymbol{B}}{\partial t} = \boldsymbol{\nabla} \times \left(\boldsymbol{v} \times \boldsymbol{B}\right). \end{equation}

The dynamical equation for the momentum is derived from $\delta S = 0$, and is thus equal to

(4.5)\begin{align} &\frac{\partial}{\partial t}\left(\rho v^{k}\right) + \partial_j \left[\rho v^j v^k + \left(p+\frac{|\boldsymbol{B}|^2}{2}\right) \delta^{jk} - B^j B^k\right] \nonumber\\ &\quad -\partial_j \left[\boldsymbol{B}\boldsymbol{\cdot}\left(\boldsymbol{\nabla}\times \boldsymbol{v}\right)\left(\rho\frac{\partial{\mathcal{F}}}{\partial\rho}+\sigma\frac{\partial{\mathcal{F}}}{\partial\sigma}+|\boldsymbol{B}|\frac{\partial{\mathcal{F}}}{\partial|\boldsymbol{B}|}-\mathcal{F} \right)\delta^{jk}\right]\nonumber\\ &\quad+ \partial_j \left[\boldsymbol{B}\boldsymbol{\cdot}\left(\boldsymbol{\nabla}\times \boldsymbol{v}\right)\left(\frac{\partial{\mathcal{F}}}{\partial|\boldsymbol{B}|}\frac{B^{j}B^{k}}{|\boldsymbol{B}|}\right)\right] \nonumber\\ &\quad + \partial_j \left[\epsilon_{kji}B^{i}\left(\frac{\partial{\mathcal{F}}}{\partial\rho}\frac{\partial\rho}{\partial t}+\frac{\partial{\mathcal{F}}}{\partial\sigma}\frac{\partial\sigma}{\partial t}+ \frac{\partial{\mathcal{F}}}{\partial|\boldsymbol{B}|}\frac{\boldsymbol{B}}{|\boldsymbol{B}|}\boldsymbol{\cdot}\frac{\partial \boldsymbol{B}}{\partial t}\right)\right] + \partial_j \left[\epsilon_{kji}{\mathcal{F}}\frac{\partial B^{i}}{\partial t}\right] \nonumber\\ &\quad -\partial_{j}\left[\epsilon_{lji}{\mathcal{F}}B^{l}\left(\partial_{k}v^{i}\right)-\epsilon_{kli}{\mathcal{F}}B^{j}\left(\partial_{l}v^{i}\right)\right] = 0, \end{align}

where repeated indices indicate summation (as per the Einstein convention), and we have employed the standard relationship between the internal energy and the scalar pressure $p$. We note that (4.5) can be obtained in two different ways from the action. The first is to follow the conventional variation with respect to $q$ and obtain it accordingly. The second method involves the use of the procedure outlined in Frieman & Rotenberg (Reference Frieman and Rotenberg1960) and Newcomb (Reference Newcomb1962) and is described in appendix A. For our model, (4.1)–(4.3) and (4.5) constitute the complete set of dynamical equations.

Before discussing the ansatz in more detail, a few observations regarding (4.5) are in order. The second term occurring in the first line of this equation represents the ideal MHD momentum flux (enclosed in square brackets), which is seen from the absence of $\mathcal {F}$ in it. The second and third lines contain terms that are purely symmetric under the interchange $k \leftrightarrow j$. The fourth line contains terms that are wholly antisymmetric under $k \leftrightarrow j$. The fifth (and final) line contains terms that are neither purely symmetric nor purely antisymmetric. As a result, we see that the entire momentum flux tensor is not symmetric, as opposed to the ideal MHD tensor, or the 2-D gyroviscous tensor for the specific model considered in Morrison et al. (Reference Morrison, Lingam and Acevedo2014). Note that we refer to the terms from line two onwards as gyroviscous because they are expressed in terms of the velocity shear, akin to viscous hydrodynamics. The gyroviscous tensor thus obtained above can be compared against the general expression(s) presented in Ramos (Reference Ramos2005b). Furthermore, these effects arise from charged particle gyration – the latter aspect is explored below.

4.2 The origin of the gyroviscous ansatz

In § 3.2, we briefly outlined the process involved in constructing a generic gyroviscous term. Now, we shall draw upon further physics to select a specific choice for the ansatz.

First, let us suppose that we start out with the notion of an internal angular momentum $\boldsymbol {L}^\star$. In order to understand where this angular momentum originates, we recall an identity from electromagnetism which relates the angular momentum to the magnetic moment via the gyromagnetic ratio, $(2m)/e$. If we consider a two-species model of ions and electrons, then the ions will play the dominant role, owing to their higher mass. Hence, we know that $\boldsymbol {L}^\star = ({2m}/{e}) \boldsymbol {\mu }$. The magnetic moment $\boldsymbol{\mu}$ is typically an adiabatic invariant in plasmas, and its magnitude is given by $|\boldsymbol {\mu }|={m \boldsymbol {v}_\perp ^2}/{2|\boldsymbol {B}|}$, which is proportional to $P_\perp /|\boldsymbol {B}|$ where $P_\perp$ denotes the perpendicular component of the (anisotropic) pressure. But, the magnetic moment is a vector and the most natural way to construct a vector is through the unit vector of the magnetic field. Putting these results together, we find that a natural ansatz (albeit a specific one) for $\boldsymbol {L}^\star$ is given by

(4.6)\begin{equation} \boldsymbol{L}^\star = \alpha \frac{m}{2e} \frac{P_\perp}{|\boldsymbol{B}|^2} \boldsymbol{B}, \end{equation}

where $\alpha$ is a dimensionless proportionality constant, which can be arbitrarily specified; in the ensuing analysis, we set $\alpha = 1$ for simplicity. By comparison with the more general ansatz outlined in § 3.2, we find that they are identical when $\mathcal {F} = \alpha ({m}/{2e}) P_\perp /|\boldsymbol {B}|^2$.

The function $P_\perp$ is a function of $\sigma$, $\rho$ and $|\boldsymbol {B}|$. For a more detailed discussion of the anisotropic pressure, we refer the reader to Kimura & Morrison (Reference Kimura and Morrison2014). It is defined as

(4.7)\begin{equation} P_\perp = \rho^2 \frac{\partial U}{\partial \rho} + \rho |\boldsymbol{B}| \frac{\partial U}{\partial |\boldsymbol{B}|}, \end{equation}

an expression that first appeared in Morrison (Reference Morrison1982), where $U$ is the internal energy that is a function of $\rho$ and $\sigma$, but also of the magnetic field; see also Hazeltine, Mahajan & Morrison (Reference Hazeltine, Mahajan and Morrison2013). If we wish to forgo anisotropy, then we assume that $U$ is independent of $B$, and hence the second term in the above term vanishes. This assumption was used in deriving the equation of motion (4.5) since the internal energy introduced in (3.3) had no $\boldsymbol {B}$-dependence. Such an assumption also leads to the pressure tensor becoming isotropic, given by the first term of (4.7) alone.

In summary, the ansatz constructed was chosen such that the gyroviscosity (and consequently the momentum transport) arises via the gyration of charged particles, thereby lending the term its name. The fact that momentum transport could take place via such gyrations was first noted by Chapman & Cowling (Reference Chapman and Cowling1970) and Kaufman (Reference Kaufman1960) in the 1950s and 1960s. This principle was applied to incompressible gyrofluids in Newcomb (Reference Newcomb1972, Reference Newcomb1973, Reference Newcomb1983) and compressible gyrofluids in Morrison (Reference Morrison, Eliasson and Shukla2009) and Morrison et al. (Reference Morrison, Lingam and Acevedo2014), who showed that this specific ansatz yielded results that were fully compatible with the 2-D version of the Braginskii tensor (Braginskii Reference Braginskii1965).

Lastly, we note that substituting (4.6) in (3.8) after employing $\boldsymbol {M}^\star = \boldsymbol {\nabla } \times \boldsymbol {L}^\star$ will yield a number of extra terms with the same dimensions as $\boldsymbol {M} = \rho \boldsymbol {v}$. Hence, if one divides the expression throughout by $\rho$, the contributions arising from $\boldsymbol {M}^\star$ have the dimensions of velocity and possess physical interpretations. The first term, which is proportional to $(\boldsymbol {B} \times \boldsymbol {\nabla } P_\perp )/|\boldsymbol {B}|^2$, amounts to the diamagnetic drift velocity. The second term, which is proportional to $P/|\boldsymbol {B}|^3 (B \times \boldsymbol {\nabla } |\boldsymbol {B}|)$ is analogous to the $\boldsymbol {\nabla } |\boldsymbol {B}|$ drift velocity for charged particles. This correspondence has been pointed out in Morrison et al. (Reference Morrison, Caldas and Tasso1984, § 6).

5.1 Constructing a hybrid conserved angular momentum

When we perform the constrained variation of our action, we recover

(5.1)\begin{equation} \left.\begin{array}{c} \dfrac{\partial M^c_i}{\partial t} + \partial_j T_{ij} = 0, \\ T_{ij} = M^c_i v_j - \dfrac{\partial \mathcal{L}}{\partial \left(\partial_j v_k\right)} \left(\partial_i v_k\right) + \dfrac{\partial \mathcal{L}}{\partial B^i} B_j + \delta_{ij} \left[\mathcal{L} - \dfrac{\partial \mathcal{L}}{\partial B^k} B^k - \dfrac{\partial \mathcal{L}}{\partial \rho} \rho \right]. \end{array}\right\} \end{equation}

Additional details can be found in Holm et al. (Reference Holm, Marsden and Ratiu1998, equations (7.6)–(7.8)) and Lingam & Morrison (Reference Lingam and Morrison2014, § 3). Note that the Lagrangian density $\mathcal {L}$ in the above expression refers to the one present in (3.7). A rather unusual fact emerges if one inspects the above energy-momentum tensor: when one considers ideal MHD, or even Hall and extended MHD, the tensor $T_{ij}$ is symmetric. In turn, this ensures that the angular momentum $\boldsymbol {M}={\boldsymbol {r}} \times \rho \boldsymbol {v}$ is conserved. However, this is evidently not the case for the above energy-momentum tensor.

This fact is not unusual because a number of hydrodynamic models are known to possess asymmetric energy-momentum tensors. In particular, if the constituent ‘particles’ (which may be fluid parcels) have an internal degree of freedom (i.e. spin), the energy-momentum tensor of the fluid will manifest a non-symmetric component (Papapetrou Reference Papapetrou1949; Snider & Lewchuk Reference Snider and Lewchuk1967; Olmsted & Snider Reference Olmsted and Snider1976; Dewar Reference Dewar1977; Evans Reference Evans1979; Kopczyński Reference Kopczyński1990; Lingam Reference Lingam2015a). Examples of hydrodynamic models with asymmetric energy-momentum tensors include ferrohydrodynamics (Rosensweig Reference Rosensweig1985; Billig Reference Billig2005) and nematics (de Gennes & Prost Reference de Gennes and Prost1993). Although many core plasma models are characterized by symmetric energy-momentum tensors (Pfirsch & Morrison Reference Pfirsch and Morrison1985; Similon Reference Similon1985), other plasma models feature asymmetric energy-momentum tensors (e.g. Brizard Reference Brizard2010a). In consequence, not all components of the angular momentum will be conserved, although the toroidal component is conserved in such models (Scott & Smirnov Reference Scott and Smirnov2010).

To resolve this, we will adopt the procedure delineated in McLennan (Reference McLennan1966). We begin with the observation that the first expression in (5.1) remains invariant under the transformations $M^c_i \rightarrow M^c_i + \partial _j {\varSigma }_{ij}$ and $T_{ij} \rightarrow T_{ij} - \partial {\varSigma }_{ij}/\partial t$. Let us suppose that we choose $\partial {\varSigma }_{ij}/\partial t$ to be the antisymmetric part of $T_{ij}$, thereby ensuring that $T_{ij} - \partial {\varSigma }_{ij}/\partial t$ is purely symmetric. Hence, by utilizing this choice of ${\varSigma }_{ij}$, we find that

(5.2)\begin{equation} \frac{\partial {\varSigma}_{ij}}{\partial t} = T^A_{ij} = \frac{1}{2} \left(T_{ij} - T_{ji}\right) = \epsilon_{ijk} \tau_k, \end{equation}

where $\boldsymbol {\tau }$ has the units of torque density and is given by

(5.3)\begin{equation} \tau_k = \frac{1}{2}\left[ \epsilon_{kab} M^c_a v_b + \frac{m}{2e} \frac{P_\perp}{|\boldsymbol{B}|^2} \left(B_k \partial_l v_l - B_l \partial_k v_l\right) \right]. \end{equation}

The first term in the above expression is ${{\boldsymbol {M}}}^c\times {\boldsymbol {v}}$, which can also be expressed as ${{\boldsymbol{M}}}^\star \times {\boldsymbol{v}}$ since $\rho {{\boldsymbol {v}} \times {\boldsymbol {v}}} = 0$. The second and third terms are proportional to $(\boldsymbol {\nabla } \boldsymbol {\cdot } {{\boldsymbol {v}}}) {\boldsymbol {B}}$ and $(\boldsymbol {\nabla } {\boldsymbol {v}}) \boldsymbol {\cdot } {\boldsymbol {B}}$, respectively. Since we know that $\tau$ behaves as a torque density, let us define a dynamical variable $\mathcal {S}$ such that $\partial \mathcal {S}_k/ \partial t = \tau _k$; this constitutes a relation that mirrors the conventional torque-angular momentum relation in classical mechanics. Using this in (5.2), we find that $\varSigma _{ij} = \epsilon _{ijk} \mathcal {S}_k$. With these ingredients, we can now construct a symmetric momentum conservation law as follows:

(5.4)\begin{equation} \frac{\partial M^\textrm{tot}_i}{\partial t} + \partial_j T^S_{ij} = 0, \end{equation}

with $T^S_{ij}$ representing the symmetric energy-momentum tensor and $M^\textrm {tot}_i = M^c_i + \epsilon _{ijk} \partial _j \mathcal {S}_k$. As the resultant energy-momentum tensor is symmetric, it follows that the corresponding angular momentum ${{\boldsymbol {r}}} \times {{\boldsymbol {M}}}^\textrm {tot}$ is conserved.

The ramifications of $\mathcal {S}$ are manifold. It can be interpreted as an intrinsic angular momentum density generated from the torque density (5.3). This is consistent with prior works (Papapetrou Reference Papapetrou1949; Snider & Lewchuk Reference Snider and Lewchuk1967; Olmsted & Snider Reference Olmsted and Snider1976; Dewar Reference Dewar1977; Evans Reference Evans1979; Kopczyński Reference Kopczyński1990) that outlined the connections between intrinsic angular momentum and a non-symmetric energy-momentum tensor. A second justification arises from ${{\boldsymbol {M}}}^\textrm {tot} = {\boldsymbol {M}}^c + \boldsymbol {\nabla } \times \mathcal {S}$, implying by dimensional analysis that $\mathcal {S}$ has the dimensions of angular momentum density. If we define ${\boldsymbol {M}}^\textrm {int} = \boldsymbol {\nabla } \times \mathcal {S}$, we see that $\boldsymbol {\nabla } \boldsymbol {\cdot } {\boldsymbol {M}}^{\textrm {int}} = 0$. The kinship between ${\boldsymbol {M}}^\star$ and ${\boldsymbol {M}}^{\textrm {int}}$ is obvious as they are both generated via an internal angular momentum mechanism and are divergence-free.

Let us now summarize our results. We defined a dynamical variable $\mathcal {S}$ such that it obeys $\partial \mathcal {S}_i/\partial t = \tau _i$ where $\tau$ is given by (5.3), and it emerges from the antisymmetric part of the original energy-momentum tensor. We also find that the new momentum ${\boldsymbol {M}}^{\textrm {tot}} = {\boldsymbol {M}} + ({\boldsymbol {M}}^\star + {\boldsymbol {M}}^{\textrm {int}})$ yields a symmetric momentum tensor (which is the symmetric part of the old one). Using the expressions for ${\boldsymbol {M}}^\star$ and ${\boldsymbol {M}}^{\textrm {int}}$, we have

(5.5)\begin{equation} {\boldsymbol{M}}^\star + {\boldsymbol{M}}^{\textrm{int}} = \boldsymbol{\nabla} \times \left({\boldsymbol{L}}^\star + \mathcal{S}\right). \end{equation}

Hence, we can define a composite intrinsic angular momentum ${\boldsymbol {J}} = {\boldsymbol {L}}^\star + \mathcal {S}$, akin to the total angular momentum in quantum mechanics (Weinberg Reference Weinberg2015). The introduction of ${\boldsymbol {J}}$ yields ${\boldsymbol {M}}^{\textrm {tot}} = {\boldsymbol {M}} + \boldsymbol {\nabla } \times {\boldsymbol {J}}$, which is simple in form and has an immediate physical interpretation. The angular momentum corresponding to ${\boldsymbol {M}}^{\textrm {tot}}$ is conserved, and is given by ${\boldsymbol {r}} \times {\boldsymbol {M}}^{\textrm {tot}}$. Hence, the total angular momentum defined below is an invariant,

(5.6)\begin{equation} \int \textrm{d}^3r\left[{\boldsymbol{r}} \times {\boldsymbol{M}} + {\boldsymbol{r}}\times \left(\boldsymbol{\nabla} \times {\boldsymbol{J}}\right)\right]. \end{equation}

Before proceeding further, some major aspects concerning the 2-D GVMHD model described in Morrison et al. (Reference Morrison, Caldas and Tasso1984) and Morrison et al. (Reference Morrison, Lingam and Acevedo2014) merit further explication. To begin with, we can rewrite (5.1) as follows:

(5.7)\begin{equation} \frac{\partial M_i}{\partial t} + \partial_j \tilde{T}_{ij} + \left[ \frac{\partial M^\star_i}{\partial t} + \partial_j \left(M^\star_i v_j \right) \right] = 0, \end{equation}

where we have introduced the new energy-momentum tensor

(5.8)\begin{equation} \tilde{T}_{ij} = \rho v_i v_j - \frac{\partial \mathcal{L}}{\partial \left(\partial_j v_k\right)} \left(\partial_i v_k\right) + \frac{\partial \mathcal{L}}{\partial B_i} B_j + \delta_{ij} \left[\mathcal{L} - \frac{\partial \mathcal{L}}{\partial B_k} B_k - \frac{\partial \mathcal{L}}{\partial \rho}. \rho \right]. \end{equation}

The first key point worth highlighting here is that Morrison et al. (Reference Morrison, Caldas and Tasso1984, Reference Morrison, Lingam and Acevedo2014) adopted: (i) a specific equation of state (EOS) for $P_\perp$ wherein $P_\perp /|\boldsymbol {B}|$ was a Lie-dragged scalar density; and (ii) the choice $\boldsymbol {B} = B_z \hat {z}$ for the magnetic field. These two conditions collectively ensured that $\boldsymbol {L}^\star$ had only one component and that the components of $\boldsymbol {M}^\star$ behaved as scalar densities that underwent Lie-dragging; in other words, the term inside the square brackets of (5.7) vanishes identically for the 2-D GVMHD model.

The second essential point is that 2-D GVMHD did not include any variables that were Lie-dragged as vector densities of rank unity. In contrast, the magnetic field in 3-D MHD and GVMHD plays this role (Morrison Reference Morrison1982; Lingam & Morrison Reference Lingam and Morrison2014),Footnote ⁴ but $B_z$ in 2-D GVMHD is a Lie-dragged scalar density as seen from Morrison et al. (Reference Morrison, Caldas and Tasso1984, equation (3)); to put it differently, $B_z$ in 2-D GVMHD is advected the same way as the plasma density $\rho$. Thus, the terms in (5.8) involving $B_i$ are rendered irrelevant because they were derived under the assumption that the magnetic field is a Lie-dragged vector density. Hence, these two facts collectively ensure that the only potential source of asymmetry in the energy-momentum tensor of 2-D GVMHD is the second term on the right-hand side of (5.8). When one utilizes the particular EOS for this model in conjunction with $M_z = 0$ and $\boldsymbol {B} = B_z \hat {z}$, it can be shown (Morrison et al. Reference Morrison, Caldas and Tasso1984, Reference Morrison, Lingam and Acevedo2014) that the gyroviscous term of 2-D GVMHD yields the contribution

(5.9a,b)\begin{equation} {T}_{ij}^{(\mathrm{2D-GV})} = \frac{m}{2e} \frac{P_\perp}{B_z} \mathcal{N}_{jlik} \partial_k v_l, \quad \quad \mathcal{N}_{jlik} = \delta_{jk} \epsilon_{li} - \delta_{li} \epsilon_{jk}, \end{equation}

to the energy-momentum tensor, which turns out to be fully symmetric.

The above discussion serves to illustrate how and why the energy-momentum tensor of the simplified 2-D GVMHD model of Morrison et al. (Reference Morrison, Caldas and Tasso1984, Reference Morrison, Lingam and Acevedo2014) is symmetric in nature. However, in order to achieve this symmetry, a number of restrictions on the EOS as well as the magnetic field and momentum density had to be imposed. When all of these constraints are relaxed, which is the case for 3-D GVMHD, one finds that an asymmetric energy-momentum tensor is obtained.

5.2 An illustration of the formalism

We have already noted earlier that the kinetic angular momentum ${r} \times {\boldsymbol {M}}$ is not conserved. However, we have seen that the angular momentum described in (5.6) is conserved. Together, these imply that the rate of loss (or gain) of the kinetic angular momentum ${r} \times {\boldsymbol {M}}$ is precisely equal to the rate of gain (or loss) of the intrinsic angular momentum ${\boldsymbol {J}}$. Let us recall that $\mathcal {S}$ comprises a part of ${\boldsymbol{J}}$, and we know that $\partial \mathcal {S}_i/\partial t = \tau _i$ where $\boldsymbol {\tau }$ is given by (5.3). The first term in (5.3) reduces to ${\boldsymbol {M}}^\star \times {\boldsymbol {v}}$, as noted earlier. It is worth mentioning that the additional two terms are quite different, but exhibit a similar scaling. Hence, we shall use only the first term in our subsequent analysis. The total torque (denoted by $\tilde {\mathcal {T}}$) is found by integrating this term over the volume, and thus gives rise to the scaling

(5.10)\begin{equation} \tilde{\mathcal{T}} = \int {\boldsymbol{M}}^\star \times \boldsymbol{v}\,\textrm{d}^3r \sim \frac{m}{e} \left(\frac{P_\perp}{|\boldsymbol{B}|}\right) \varOmega R^3, \end{equation}

where we have dropped the numerical factors and used a characteristic velocity of $\varOmega R$, with $R$ denoting the radius of the (spherical) object. It is evident that the scaling will be entirely determined by the EOS that is adopted.

Next, let us evaluate the spin-down rate, by using the relation $\tilde {\mathcal {T}} = I \dot {\varOmega }$, from classical mechanics. The moment of inertia, dropping all numerical factors, is approximately $M R^2 \sim \rho R^5$. Using this in (5.10), we find that

(5.11)\begin{equation} \dot{\varOmega} \sim \frac{m}{e} \left( \frac{P_\perp}{\rho |\boldsymbol{B}|} \right) \varOmega R^{-2}. \end{equation}

The above relation indicates that $\dot {\varOmega } \propto \varOmega$ (holding other quantities fixed). The EOS depends only on $\rho$, $s$ and $|\boldsymbol {B}|$ and hence we can conclude that the relation $\dot {\varOmega } \propto \varOmega$ is likely to be independent of the choice of the EOS. If we treat $\rho$ and $R$ to be independent variables, i.e. by choosing $M$ to be the dependent variable, one can also conclude that $\dot {\varOmega } \propto R^{-2}$ will be independent of the EOS. The characteristic time $t_c = \varOmega /\dot {\varOmega }$, is expected to be independent of $\varOmega$ and is given by

(5.12)\begin{equation} t_c \sim \frac{e}{m} \left(\frac{\rho |\boldsymbol{B}|}{P_\perp}\right) R^2, \end{equation}

and we see that it is proportional to $R^2$, when the other parameters are held constant. The Chew–Goldberger–Low EOS for $P_\perp$ (Chew, Goldberger & Low Reference Chew, Goldberger and Low1956) is of particular interest since the characteristic time $t_c$ and the rate $\dot {\varOmega }$ are both independent of the density and the magnetic field, thereby demonstrating an unexpected universality. The resulting spin-down corresponds to the dissipation of kinetic angular momentum, which must imply that there is a corresponding increase in the intrinsic angular momentum ${\boldsymbol{J}}$ (which comprises the other fluid variables).

The spin rates of low-mass stars are found to slow down by approximately two orders of magnitude over a span of $10^9$ years (Scholz Reference Scholz and Stempels2009). Modelling stellar spin-down is important for a multitude of reasons, including the fact that the older stars (with lower rotation rates) display lower activity in general, which has numerous ramifications for planetary habitability (Lingam & Loeb Reference Lingam and Loeb2018, Reference Lingam and Loeb2019). We can estimate the characteristic time by choosing solar parameters (i.e. a solar-type star) for an order-of-magnitude calculation. In particular, we substitute $|\boldsymbol {B}| \sim 10^{-4}$ T, $R \sim 7 \times 10^8\ \textrm {m}$ and $T \sim 5.8 \times 10^3$ K (Priest Reference Priest2014) in (5.12), which yields $t_c \sim 3 \times 10^6$ years. The two leading candidates invoked to explain stellar spin-down, star-disk and stellar wind braking, operate on time scales of ${\sim }10^6\text {--}10^7$ years and ${\sim }10^8$ years, respectively (Bouvier et al. Reference Bouvier, Matt, Mohanty, Scholz, Stassun, Zanni, Beuther, Klessen, Dullemond and Henning2014, § 4.1). Hence, we see that our semiquantitative estimate is comparable to these two time scales, and may therefore constitute a viable mechanism for governing angular momentum evolution of solar-mass stars.

The issue of angular momentum losses in protostars is another closely related topic (Bodenheimer Reference Bodenheimer1995; Matt & Pudritz Reference Matt and Pudritz2005; Hartmann, Herczeg & Calvet Reference Hartmann, Herczeg and Calvet2016) which might also be resolvable through the same mechanism. We emphasize that the heuristic treatment in this subsection has primarily relied on simple scaling arguments, and a complete picture can only emerge through the synthesis of rigorous analytical models and numerical simulations. We note that this only represents the tip of the iceberg – other potential applications include pulsar braking, transport in accretion discs and associated phenomena. In the realm of fusion, we note that the formalism developed herein may prove to be useful in explaining intrinsic rotation observed in tokamaks (Gürcan et al. Reference Gürcan, Diamond, Hahm and Singh2007; de Grassie Reference de Grassie2009; Diamond et al. Reference Diamond, Kosuga, Gürcan, McDevitt, Hahm, Fedorczak, Rice, Wang, Ku and Kwon2013; Rice Reference Rice2016).

6.1 The Lagrangian view point and the Lagrange–Euler map

First, note that the Hamiltonian can be obtained from the Lagrangian via a Legendre transform, akin to the usual process in particle mechanics. The Hamiltonian is given by

(6.1)\begin{equation} H[\boldsymbol{q},\boldsymbol{\varPi}] = \int_D \textrm{d}^3a\,\dot{\boldsymbol{q}}\boldsymbol{\cdot} \boldsymbol{\varPi} - L, \end{equation}

where

(6.2)\begin{equation} L[\dot{\boldsymbol{q}},\boldsymbol{q}]= \int_D \textrm{d}^3a\,\mathcal{L}\left(\boldsymbol{q},\dot{\boldsymbol{q}},\partial \boldsymbol{q}/\partial \boldsymbol{a}\right) , \end{equation}

with $L$ defined so that the action of (3.7) is given by $S=\int _T \,\textrm {d}t\, L$. Consequently, the canonical momentum is given by

(6.3)\begin{equation} \boldsymbol{\varPi} = \frac{\delta L}{\delta \dot{\boldsymbol{q}}}=\rho_0\dot{\boldsymbol{q}} + \boldsymbol{\varPi}^{\star}, \end{equation}

and we see that we have a field theory counterpart to the finite-dimensional case for particle motion in a magnetic field where the kinetic momentum differs from the canonical momentum, here with the role of the vector potential being played by $\boldsymbol {\varPi }^{\star }$. Thus (6.1) gives the Hamiltonian

(6.4)\begin{equation} H= \int_D \textrm{d}^3a\left( \frac{|\boldsymbol{\varPi}-\boldsymbol{\varPi}^{\star}|^2}{2\rho_0} +\rho_0 U\left(\frac{\rho_0}{\mathcal{J}},\frac{\sigma_0}{\mathcal{J}}\right) + \frac{1}{2 \mathcal{J}} q^{\hspace{1pt}i}_{\,,\hspace{1pt}j} \, q^{\hspace{1pt}i}_{\,,\hspace{1pt}k} \, B_0^j B_0^k \right). \end{equation}

This Hamiltonian (6.4) together with the canonical Poisson bracket,

(6.5)\begin{equation} \left\{F,G\right\} = \int_D \textrm{d}^3a \left(\frac{\partial F}{\partial \boldsymbol{q}} \cdot \frac{\partial G}{\partial \boldsymbol{\varPi}} - \frac{\partial G}{\partial \boldsymbol{q}} \cdot \frac{\partial F}{\partial \boldsymbol{\varPi}}\right), \end{equation}

generates the Hamiltonian equations of motion in Lagrangian variables for our class of 3-D GVMHD models as follows:

(6.6)\begin{equation} \dot{\boldsymbol{q}}= \{\boldsymbol{q},H\}= \frac{\delta H}{\delta \boldsymbol{\varPi}}\quad \mathrm{and} \quad \dot{\boldsymbol{\varPi}}=\{\boldsymbol{\varPi},H\}= - \frac{\delta H}{\delta \boldsymbol{q}}, \end{equation}

equations equivalent to the Euler–Lagrange equations obtained via $\delta S = 0$.

Now, one can use the Lagrange–Euler maps to convert both the Hamiltonian and the bracket into Eulerian variables. The procedure is described in the next section. We will see that the origin of the gyromap lies in (6.3) and how this expression relates to different choices of Eulerian variables. The bracket obtained in terms of any of these choices is endowed with Lie algebraic properties (Morrison Reference Morrison1998), most importantly the Jacobi identity, but it does not possess the canonical form of (6.5) because the Eulerian variables are not a set of canonical variables. As a result, one refers to the Hamiltonian and the bracket as being noncanonical in nature, and indeed one version is identical to that originally given in Morrison & Greene (Reference Morrison and Greene1980).

As the Lagrange–Euler maps are not one-to-one, the noncanonical brackets are degenerate in general, which gives rise to the existence of invariants – the Casimirs. The theory of Casimir invariants has been studied quite extensively (Morrison Reference Morrison1998, Reference Morrison2005; Holm Reference Holm2008), but there are still unresolved subtleties regarding their incompleteness, see for example Yoshida, Morrison & Dobarro (Reference Yoshida, Morrison and Dobarro2014) and Yoshida & Morrison (Reference Yoshida and Morrison2014, Reference Yoshida and Morrison2016).

The Casimirs also possess several advantages of their own, such as variational principles for Eulerian equilibria of the form

(6.7)\begin{equation} \delta F := \delta (H + \lambda C) = 0, \end{equation}

where $C$ represents any combination of all the known Casimirs. This procedure is known as the energy–Casimir method. Once the equilibria are known, the following symmetric operator can be constructed:

(6.8)\begin{equation} \varLambda_{jk}:= \frac{\delta^2 F}{\delta \psi^j \delta \psi^k}, \end{equation}

where $F$ is defined in (6.7) and the $\psi$ denote the Eulerian (noncanonical) variables. The energy–Casimir method states that the positive-definiteness of this operator is a sufficient condition for stability, although there are mathematical intricacies involved (Holm et al. Reference Holm, Marsden, Ratiu and Weinstein1985; Rein Reference Rein1994; Batt, Morrison & Rein Reference Batt, Morrison and Rein1995; Morrison Reference Morrison1998; Yoshida et al. Reference Yoshida, Ohsaki, Ito and Mahajan2003). Thus, the Eulerian noncanonical Hamiltonian description we obtain allows for implementation of such energy principles, although we will not pursue this application here.

6.2 The gyro-bracket

We shall choose our new set of observables to be the Eulerian variables $\{\boldsymbol {M}^{c},\rho ,\sigma ,\boldsymbol {B}\}$, where $\boldsymbol {M}^c$ was defined in (3.8). The reason for this choice will soon become obvious. Recall that the Lagrange–Euler maps can be expressed in an integral form, as they were for the density $\rho$ and canonical momentum density $\boldsymbol {M}^c$ in (2.6) and (2.7), respectively. The remaining Eulerian variables are given by

(6.9)\begin{gather} \sigma=\int_D \textrm{d}^{3}a\delta\left(\boldsymbol{r}-\boldsymbol{q}\right)\sigma_{0}(\boldsymbol{a}), \end{gather}

(6.10)\begin{gather}B^{j}=\int_D \textrm{d}^{3}a\delta\left(\boldsymbol{r}-\boldsymbol{q}\right)q_{,k}^{j}B_{0}^{k}(\boldsymbol{a}). \end{gather}

We use these expressions to obtain the noncanonical bracket from the canonical counterpart by the functional chain rule. Any functional of the Eulerian observables can be expressed in terms of $\boldsymbol {\varPi }$ and $\boldsymbol {q}$; hence to delineate, we denote functionals of $\boldsymbol {\varPi }$ and $\boldsymbol {q}$ by $\bar {F}$ and those in terms of the observables by $F$, and note symbolically that $\bar {F}={F}\circ \mathfrak {E}$; consequently,

(6.11)\begin{equation} \int_D \textrm{d}^{3}a\left[\frac{\delta\bar{F}}{\delta\boldsymbol{\varPi}}\cdot\delta\boldsymbol{\varPi}+ \frac{\delta\bar{F}}{\delta \boldsymbol{q}}\cdot\delta \boldsymbol{q}\right]=\int_D \textrm{d}^{3}r \left[\frac{\delta F}{\delta \boldsymbol{M}^{c}}\cdot\delta \boldsymbol{M}^{c}+\frac{\delta F}{ \delta\boldsymbol{B}}\cdot\delta \boldsymbol{B}+\frac{\delta F}{\delta\rho}\delta\rho+ \frac{\delta F}{\delta\sigma}\delta\sigma\right]. \end{equation}

From (2.6), we can conclude that

(6.12)\begin{equation} \delta\rho=-\int_D \textrm{d}^{3}a\rho_{0}\boldsymbol{\nabla}\delta\left(\boldsymbol{r}-\boldsymbol{q}\right)\boldsymbol{\cdot}\delta \boldsymbol{q}, \end{equation}

and similar identities can be found for (2.6), (2.7), (6.9) and (6.10) as well. We substitute these identities into (6.11) and carry out integrations by parts, followed by a subsequent change in the order of integration. This results in terms that are dotted with $\delta \boldsymbol {q}$ and terms dotted with $\delta \boldsymbol {\varPi }$ on both the left- and right-hand sides of the expression. As $\delta \boldsymbol {q}$ and $\delta \boldsymbol {\varPi }$ are independent, these terms must balance and thereby we obtain relationships between the Eulerian and Lagrangian functional derivatives. The algebra involved is complicated, but quite straightforward, and we refer the reader to Morrison (Reference Morrison, Eliasson and Shukla2009) for a more pedagogical version. The final bracket that we obtain is found to be

(6.13)\begin{align} \left\{F,G\right\} &= -\int \textrm{d}^{3}r \left[M_{i}^{c}\left(\frac{\delta F}{\delta M_{j}^{c}} \partial_{j}\frac{\delta G}{\delta M_{i}^{c}}-\frac{\delta G}{\delta M_{j}^{c}}\partial_{j} \frac{\delta F}{\delta M_{i}^{c}}\right) \right. \nonumber\\ &\quad + \rho\left(\frac{\delta F}{\delta M_{j}^{c}}\partial_{j}\frac{\delta G}{\delta\rho}- \frac{\delta G}{\delta M_{j}^{c}}\partial_{j}\frac{\delta F}{\delta\rho}\right) + \sigma\left(\frac{\delta F}{\delta M_{j}^{c}}\partial_{j}\frac{\delta G}{\delta\sigma}- \frac{\delta G}{\delta M_{j}^{c}}\partial_{j}\frac{\delta F}{\delta\sigma}\right) \nonumber\\ &\left.\quad +\, B^{i}\left(\frac{\delta F}{\delta M_{j}^{c}}\partial_{j}\frac{\delta G}{\delta B^{i}}- \frac{\delta G}{\delta M_{j}^{c}}\partial_{j}\frac{\delta F}{\delta B^{i}}\right) + B^{i}\left(\frac{\delta G}{\delta B^{j}}\partial_{i}\frac{\delta F}{\delta M_{j}^{c}}- \frac{\delta F}{\delta B^{j}}\partial_{i}\frac{\delta G}{\delta M_{j}^{c}}\right)\right]. \end{align}

By inspection, one notices that the bracket derived above is exactly the same as the 3-D ideal MHD bracket (Morrison & Greene Reference Morrison and Greene1980); however, here the canonical momentum $\boldsymbol {M}^c$ replaces the kinetic momentum $\boldsymbol {M} = \rho \boldsymbol {v}$.

Since the bracket (6.13) uses $\boldsymbol {M}^{c}$ as one of its observables, we must express our Hamiltonian in terms of this observable (and the others) as well. Because of the closure principle we know this is possible; indeed, (6.4) in Eulerian variables becomes

(6.14)\begin{equation} H=\int \textrm{d}^{3}r\left[\frac{\left|\boldsymbol{M}^{c}-\boldsymbol{M}^\star\right|^{2}}{2\rho}+\rho U\left(\rho,\sigma,| \boldsymbol{B}|\right)+\frac{|\boldsymbol{B}|^{2}}{2}\right]. \end{equation}

The Hamiltonian of (6.14) with the bracket of (6.13) generates our class of 3-D GVMHD models in the form

(6.15)\begin{equation} \frac{\partial \mathfrak{E}}{\partial t}=\{\mathfrak{E},H\}. \end{equation}

On account of the fact that $\boldsymbol {M} = \boldsymbol {M}^c - \boldsymbol {M}^\star$, the energy has the form identical to that of ideal MHD. This is analogous to the fact that the kinetic energy for a charged particle in a magnetic field is identical to that for a free particle. Hence, there is a choice: one can either work with the standard ideal MHD bracket and the more complicated Hamiltonian of (6.14) in terms of the canonical momentum $\boldsymbol {M}^c$, or work with a complicated bracket written in terms of the variable $\boldsymbol {M}$, the conventional variable of magnetofluid theories, and the simpler ideal MHD Hamiltonian. To obtain the bracket in terms of $\boldsymbol {M}$ we can use the gyromap (3.8), $\boldsymbol {M} = \boldsymbol {M}^c - \boldsymbol {M}^\star$, in another chain rule calculation to transform from $\boldsymbol {M}^c$ to the variable $\boldsymbol {M}$. This is worked out in appendix B for the case $\boldsymbol{M}^\star = \boldsymbol {\nabla } \times (\mathcal {F} \boldsymbol {B})$, giving rise to a complicated Poisson bracket.

Given that the noncanonical Poisson bracket in terms of $\boldsymbol {M}^c$ is the same as that of ideal MHD, it possesses the same Casimir invariants as the ideal MHD case if we replace $\boldsymbol {M}$ with $\boldsymbol {M}^c$. This use of the gyromap to obtain Casimirs, which first appeared in Morrison et al. (Reference Morrison, Caldas and Tasso1984) and subsequently in other cases (Hazeltine et al. Reference Hazeltine, Hsu and Morrison1987; Izacard et al. Reference Izacard, Chandre, Tassi and Ciraolo2011; Lingam & Morrison Reference Lingam and Morrison2014; Morrison et al. Reference Morrison, Lingam and Acevedo2014), differs from most of the prior studies that have sought to derive Casimirs and other conserved invariants via the HAP approach using a variety of methods, see for example Morrison (Reference Morrison1982, Reference Morrison1998), Padhye & Morrison (Reference Padhye and Morrison1996a,Reference Padhye and Morrisonb), Hameiri (Reference Hameiri2004) and Webb et al. (Reference Webb, Dasgupta, McKenzie, Hu and Zank2014a, Reference Webb, Dasgupta, McKenzie, Hu and Zankb) for a comprehensive discussion of the same.

So, for our present general gyroviscous models, the gyromap tells us that the $\boldsymbol {M}$-independent Casimirs of ideal MHD will be unchanged, an example being the magnetic helicity $\int \textrm {d}^3r\, \boldsymbol {A}\boldsymbol {\cdot } \boldsymbol {B}$. On the other hand, the cross-helicity and other $\boldsymbol {M}$-dependent invariants are modified by the replacement $\boldsymbol {M}\rightarrow \boldsymbol {M}^c$. Thus, the new cross-helicity Casimir is given by

(6.16)\begin{align} \int \textrm{d}^3r \frac{\boldsymbol{M}^c \boldsymbol{\cdot} \boldsymbol{B}}{\rho} &= \int \textrm{d}^3r \frac{(\boldsymbol{M} + \boldsymbol{M}^{\star}) \boldsymbol{\cdot}\boldsymbol{B}}{\rho} \end{align}

(6.17)\begin{align} &=\int \textrm{d}^3r \left[\boldsymbol{v} \boldsymbol{\cdot} \boldsymbol{B} + \left(\frac{\mathcal{F}}{\rho}\right) \boldsymbol{B} \boldsymbol{\cdot} \boldsymbol{\nabla} \times \boldsymbol{B} \right]. \end{align}

Equation (6.16) is conserved for any choice of $\boldsymbol {M}^{\star }$ that satisfies the closure principle with a provision similar to that for conservation of the usual helicity of MHD, viz. that the flow be barotropic. In (6.17) we have inserted the special case of $\boldsymbol {M}^{\star }=\boldsymbol {\nabla }\times \boldsymbol {L}^{\star }$ with (3.10). The second term of (6.17) is proportional to the current helicity density, which is encountered regularly in the context of MHD (Moffatt Reference Moffatt1978; Krause & Raedler Reference Krause and Raedler1980; Brandenburg & Subramanian Reference Brandenburg and Subramanian2005; Rincon Reference Rincon2019) and Hall MHD (Mininni, Gómez & Mahajan Reference Mininni, Gómez and Mahajan2003; Lingam & Mahajan Reference Lingam and Mahajan2015; Lingam & Bhattacharjee Reference Lingam and Bhattacharjee2016a,Reference Lingam and Bhattacharjeeb; Mahajan & Lingam Reference Mahajan and Lingam2015, Reference Mahajan and Lingam2020) turbulence and dynamo theory.