1. Introduction
Flows in magnetic confinement devices are known to have a strong influence in suppressing turbulence and affecting neoclassical (Hinton & Wong Reference Hinton and Wong1985; Helander Reference Helander1998) and turbulent transport (Burrell Reference Burrell1997), as well as stability (Chu et al. Reference Chu, Greene, Jensen, Miller, Bondeson, Johnson and Mauel1995). In particular, sheared flows like zonal flows can reduce turbulent transport by shearing the turbulent eddies (Lin et al. Reference Lin, Hahm, Lee, Tang and White1998; Terry Reference Terry2000; Fujisawa Reference Fujisawa2008). Study of the interaction of plasma flows and turbulence is critical (Stroth, Manz & Ramisch Reference Stroth, Manz and Ramisch2011; Groebner, Burrell & Seraydarian Reference Groebner, Burrell and Seraydarian1990) in understanding the edge region of tokamaks and the Low-High (L-H) transition. It has been experimentally verified that the radial electric field plays a key role in the transition (Burrell Reference Burrell1999). Flows can also have applications in healing magnetic islands in stellarators (Hegna Reference Hegna2012). Flows are also an essential object of study in space and astrophysical plasmas (Istomin & Pariev Reference Istomin and Pariev1996; Goedbloed, Goedbloed & Poedts Reference Goedbloed, Goedbloed and Poedts2004; Beskin Reference Beskin2009; Davis & Tchekhovskoy Reference Davis and Tchekhovskoy2020).
If the flow speed is assumed of the order of the sound speed, ideal magnetohydrodynamics (MHD) is a good description of plasmas since gyroradius corrections do not enter (Morozov & Solov'ev Reference Morozov and Solov'ev1980; Goedbloed et al. Reference Goedbloed, Goedbloed and Poedts2004; Freidberg Reference Freidberg2014). The theory of steady axisymmetric plasma flows is well developed (Hameiri Reference Hameiri1983; Hassam Reference Hassam1996; Beskin Reference Beskin1997; Tasso & Throumoulopoulos Reference Tasso and Throumoulopoulos1998; Guazzotto & Betti Reference Guazzotto and Betti2005; McClements & Hole Reference McClements and Hole2010; Abel et al. Reference Abel, Plunk, Wang, Barnes, Cowley, Dorland and Schekochihin2013). Toroidal angular momentum conservation, which follows from Noether's theorem, facilitates the reduction of axisymmetric flows to the study of a modified Grad–Shafranov (GS) equation which includes the effect of steady plasma flow (Hameiri Reference Hameiri1983; Hassam Reference Hassam1996; Beskin Reference Beskin1997; Goedbloed & Lifschitz Reference Goedbloed and Lifschitz1997; Tasso & Throumoulopoulos Reference Tasso and Throumoulopoulos1998). Similarly, for plasma systems with helical symmetry, a helical GS with large flows can be obtained (Villata & Tsinganos Reference Villata and Tsinganos1993; Andreussi, Morrison & Pegoraro Reference Andreussi, Morrison and Pegoraro2012). Global solutions of such flow-modified GS equations provide useful models of plasma jets in space (Villata & Ferrari Reference Villata and Ferrari1994; Beskin Reference Beskin1997; Bogoyavlenskij Reference Bogoyavlenskij2000a,Reference Bogoyavlenskijb).
While the standard GS equation for static equilibrium is always elliptic, the flow-modified GS equation can be elliptic, parabolic or hyperbolic, depending on the flow's Mach number. Therefore, even with symmetry, complicated flows such as transonic flows (Morawetz Reference Morawetz1985) that modify the characteristics in a given domain can appear (Lifschitz & Goedbloed Reference Lifschitz and Goedbloed1997).
There is considerable interest in understanding the effects of symmetry breaking three-dimensional (3-D) perturbations and stability of 3-D MHD flows in astrophysical plasmas (Istomin & Pariev Reference Istomin and Pariev1996; Igumenshchev, Narayan & Abramowicz Reference Igumenshchev, Narayan and Abramowicz2003; Igumenshchev Reference Igumenshchev2008; McKinney & Blandford Reference McKinney and Blandford2009). In particular, 2-D and 3-D simulations of magnetically arrested disks show very different behaviours (Igumenshchev Reference Igumenshchev2008; White, Stone & Quataert Reference White, Stone and Quataert2019). However, high-resolution 3-D simulations are required to carefully capture the effects of symmetry breaking (White et al. Reference White, Stone and Quataert2019).
In the absence of symmetry, reduction to a GS-like equation fails in general. Non-symmetric toroidal MHD equilibria with flows share the same mathematical difficulty as the static non-symmetric MHD equilibrium, namely, the appearance of singularities on rational surfaces.
Although large steady shear flows are desirable in various confinement geometries, steady flows show a strong alignment with a symmetry direction when present (Spong Reference Spong2005; Helander Reference Helander2007; Helander & Simakov Reference Helander and Simakov2008; Helander Reference Helander2014). In axisymmetric tokamaks, the persistent flow is in the toroidal symmetry direction (Tasso & Throumoulopoulos Reference Tasso and Throumoulopoulos1998; Throumoulopoulos, Weitzner & Tasso Reference Throumoulopoulos, Weitzner and Tasso2006). While the toroidal angular momentum in a tokamak is typically conserved over several inverse collision frequency (Helander et al. Reference Helander, Beidler, Bird, Drevlak, Feng, Hatzky, Jenko, Kleiber, Proll and Turkin2012), poloidal flows are damped due to neoclassical effects. Axisymmetry makes tokamaks intrinsically ambipolar, i.e. the radial fluxes of ions and electrons are equal and independent of the radial electric field in the leading-order gyroradius expansion, leading to flow in the symmetry direction being unconstrained.
In the absence of intrinsic ambipolarity, the rotation is damped to the value required for ambipolar radial transport and typically leads to flows of the order of diamagnetic flows (Helander & Simakov Reference Helander and Simakov2008; Helander Reference Helander2014). Plasma rotations in tokamaks and stellarators, therefore, differ significantly because stellarators do not, in general, have a symmetry direction and are not intrinsically ambipolar (Helander et al. Reference Helander, Beidler, Bird, Drevlak, Feng, Hatzky, Jenko, Kleiber, Proll and Turkin2012). The only class of stellarators that have been shown to support large flows are quasisymmetric (Tessarotto et al. Reference Tessarotto, Johnson, White and Zheng1996; Spong Reference Spong2005; Canik et al. Reference Canik, Anderson, Anderson, Likin, Talmadge and Zhai2007; Helander Reference Helander2007; Helander & Simakov Reference Helander and Simakov2008). In a quasisymmetric stellarator, continuous symmetry comes from the symmetry of the magnitude of the magnetic field. As a result of the symmetry, such stellarators are also approximately intrinsically ambipolar (Helander & Simakov Reference Helander and Simakov2008). However, it is difficult, if not impossible, to obtain exact global quasisymmetry in a volume (Garren & Boozer Reference Garren and Boozer1991; Landreman & Sengupta Reference Landreman and Sengupta2019). The breakdown of quasisymmetry leads to flow damping (Spong Reference Spong2005; Helander & Simakov Reference Helander and Simakov2008; Simakov & Helander Reference Simakov and Helander2011). In the literature, it has been further argued (Tessarotto et al. Reference Tessarotto, Johnson, White and Zheng1996; Simakov & Helander Reference Simakov and Helander2011; Sugama et al. Reference Sugama, Watanabe, Nunami and Nishimura2011) that quasisymmetry is not enough to accommodate the large centrifugal forces from sonic flows, and one needs axial/helical symmetries.
As a result of the difficulties imposed by the 3-D geometry, a lack of symmetry that leads to magnetic resonances on rational surfaces and stringent restrictions on plasma flows from neoclassical physics, a complete mathematical description of the characteristics of a generic 3-D ideal MHD equilibrium, let alone the stability of such systems, is not available in the standard literature. The analytic description of non-symmetric 3-D flows has employed various approximations such as incompressibility (Kamchatnov Reference Kamchatnov1982) or large aspect ratio (Kovrizhnykh & Shchepetov Reference Kovrizhnykh and Shchepetov1989, Reference Kovrizhnykh and Shchepetov1980).
Several variational formalisms for 3-D plasma flow problems exist (Greene & Karlson Reference Greene and Karlson1969; Vladimirov & Moffatt Reference Vladimirov and Moffatt1995; Hameiri Reference Hameiri1998; Ilgisonis & Pastukhov Reference Ilgisonis and Pastukhov2000). A non-trivial result from the variational formalism is the proof that, under certain conditions, ideal MHD allows more than one kind of cross-helicity (Hameiri Reference Hameiri1981; Vladimirov & Moffatt Reference Vladimirov and Moffatt1995; Ilgisonis & Pastukhov Reference Ilgisonis and Pastukhov2000). We shall call the additional symmetry vector associated with the conserved quantity Hameiri's vector because it was first discovered in Hameiri (Reference Hameiri1981) and later derived independently in Ilgisonis & Pastukhov (Reference Ilgisonis and Pastukhov2000) and Vladimirov & Moffatt (Reference Vladimirov and Moffatt1995) (see discussion in Hameiri (Reference Hameiri1998) and Ilgisonis & Pastukhov (Reference Ilgisonis and Pastukhov2000)). A major disadvantage of the flow-modified variational formalism is that, unlike the static limit, the energy functional for the flow problem does not possess a minimum unless the equilibrium flow is small enough or is parallel to the magnetic field (Hameiri Reference Hameiri1998).
In this work, instead of addressing the flow damping problem in a 3-D system, we take a step back and analyse the ideal MHD equilibrium with steady flows in a non-symmetric geometry. Our goal here is to construct a class of perturbative solutions to non-symmetric ideal MHD with steady flows and obtain consistency conditions for steady flows when the magnetic fields close on themselves. The flows that we consider are larger than diamagnetic flows but are not sonic. We show that the same procedure that allowed us to obtain the non-symmetric ideal MHD static equilibrium in a low shear stellarator can be generalized to MHD with such flows. In particular, we highlight the crucial role of closed field lines (Grad Reference Grad1971; Weitzner Reference Weitzner2014, Reference Weitzner2016) in avoiding magnetic resonances and going to arbitrarily high order in the perturbative analysis. We also show that for a special class of flows that possess the additional symmetry, Hameiri's vector, one can obtain a generalized GS equation (GGS) similar to the quasisymmetric generalized GS equation first derived in Burby, Kallinikos & MacKay (Reference Burby, Kallinikos and MacKay2020a) (see also related works by Burby, Kallinikos & MacKay Reference Burby, Kallinikos and MacKay2020b; Constantin, Drivas & Ginsberg Reference Constantin, Drivas and Ginsberg2021). The GGS reduces to the standard flow-modified GS equation in symmetric geometries. We explicitly obtain the necessary conditions that need to be satisfied on rational surfaces. If these conditions are not satisfied current singularities can develop on rational surfaces (Grad Reference Grad1967; Solov'ev & Shafranov Reference Solov'ev and Shafranov1970; Boozer Reference Boozer1981).
The outline of the paper is as follows. First, in § 2, we study the ideal MHD equilibrium with flows and identify the characteristic surfaces. Then, in § 3 we study the perturbation of a class of exact solutions and obtain the dispersion relation. Then, we outline in § 4 how the perturbation theory can be carried out to all orders for nearly parallel flows. Next, we utilize Hameiri's vector to derive the GGS and the necessary conditions on rational surfaces in § 5. Finally, we discuss the implications of our results in § 6.
2. Preliminaries
We employ the standard model of ideal MHD to represent the plasma in steady flow state but we add the simplifying assumption that the entropy is a constant in the entire plasma. The system is
Here, $\rho$, $p$, $\boldsymbol{u}$, $\boldsymbol{B}$, $\boldsymbol{J}$ denotes density, pressure, flow velocity, magnetic field and current respectively. It is often convenient to replace (2.1d) by the equivalent expression
While (2.1d) and (2.1d*) are equivalent, the form (2.1d) makes explicit the need to select specific Electromotive force (EMF), equivalently periods associated with the potential $\varPhi$. As noted we add the assumption
to close the system. For concreteness, we shall assume an adiabatic equation of state of the form $p\sim \rho ^\gamma$, where $\gamma$ is the ratio of specific heats. If we introduce
then an equivalent form for (2.1b) is
We start with a formal mathematical exploration of the system (2.1a)–(2.2). It is clear there are eight scalar first-order differential equations for the eight scalar unknowns $\rho,\varPhi,\boldsymbol {u},\boldsymbol {B}$. If we were to use (2.1d*) instead of (2.1d), there would be one fewer unknown but also a hidden solvability condition for (2.1d*) that adds an unknown and leads back to (2.1d). We now wish to determine the ‘type’ of the system, elliptic, hyperbolic, parabolic, mixed or whatever may be. A closely related issue is the description of the standing waves in the linearized system. We linearize about a constant state for $\rho, \boldsymbol {u}$ and $\boldsymbol {B}$ and look for waves propagating as $\exp {(\textrm {i} \boldsymbol {k} \boldsymbol {\cdot } \boldsymbol {x})}$. Such solutions are special cases of simple waves in fluid dynamics (Courant & Friedrichs Reference Courant and Friedrichs1999). With $c^2\equiv \textrm {d}p/\textrm {d}\rho$, we obtain the following linear system
It follows easily that, given the first-order perturbations $\partial \rho$, $\delta \boldsymbol {u}$ and $\delta \boldsymbol {B}$, they satisfy the system (2.4a)–(2.4d).
From (2.5c) we find
while (2.5b) together with (2.5a) yields
Taking components of (2.5) in $\boldsymbol {k}$ and $\boldsymbol {B}$ directions, we obtain
A simple calculation shows that the condition for a solution is the ‘dispersion relation’
We observe that this condition is satisfied identically for all vectors $\boldsymbol {k}$ perpendicular to $\boldsymbol {u}$ and $\boldsymbol {B}$. The expression also shows connections to the fast–slow wave combination found in ideal MHS wave propagation studies. Note, however, that no Alfvén wave appears. Moreover, the magnitude of $\boldsymbol {k}$ does not appear in the system; only its direction occurs. We study the dispersion relation (2.7) and its solutions in much more detail in Appendix A. In particular, we show that real physical solutions of (2.7) exist.
An alternate interpretation of the relation is also significant. We might pose the question: What are the characteristic surfaces for the steady flow system? Characteristic surfaces are defined by the following property. Suppose one specifies $\rho$, $\boldsymbol {u}$ and $\boldsymbol {B}$ on a given surface. Can one then determine the normal derivative of these quantities on that surface? If so, the surface is not characteristic. On a characteristic surface, the flow variables are not arbitrary; they must be related. These surfaces define important properties of any solution. Suppose $\chi =0$ is such a surface and define $\boldsymbol {k}$ as the gradient of $\chi$. Hence, $\boldsymbol {k}$ is normal to the surface. If one further defines $\delta \rho, \delta \boldsymbol {u}$ and $\delta \boldsymbol {B}$ as the gradients of these variables dotted into the normal to the surface, then these quantities satisfy the system (2.7), except that there would be inhomogeneous terms added to the equations, terms given solely in terms of the variables on the given surface. With this interpretation of the system, the relation is the condition that one cannot give the unknown functions arbitrarily on the surface in question. Such peculiar surfaces, the characteristic surface, are thus defined by (2.7).
A physical problem without such real surfaces is typically elliptic. When they exist, there are usually hyperbolic properties of the system as well. Their appearance indicates the complexity of steady MHD flow. Flows with both elliptic and hyperbolic characteristics can occur in ideal MHD equilibria. The mixed nature (neither fully elliptic nor hyperbolic) of MHD leads to difficult unresolved mathematical issues.
3. A simple linearized problem
In order to gain some insight into the nature of steady flow states we examine a simple linearized problem. We start from an exact solution of the system (2.1a)–(2.1d)
We linearize about the state and introduce the perturbed variables with the structure
We next give the linearized system of equations for the many amplitudes in (3.2). The two divergent relations are
while $\varPi ^{(1)}$ is specified as
where $c^2$ is the sound speed $\textrm {d}p/\textrm {d}\rho$ at the corresponding value of $\rho (x)$. Ohm's law reads
Finally the three momentum balance equations are
We introduce the notation
We see that $\rho ^{(1)}$ occurs only in (3.3b) and (3.4). It is convenient to eliminate $\rho ^{(1)}$ from the problem and we replace (24) by
Provided the series retains the structure in the angles given in (3.2) at all orders $m=n=0$ is never possible in the two equations: Ohm's law (3.5b) and (3.5c). Thus, we may replace (3.5b), (3.5c) by the following
We may carry out similar reductions of pressure balance and from (3.5b) and (3.5c) we obtain
The $x$ component of (3.6a) has the simpler form
At this point we have five equations (3.3a), (3.5a), (3.4*) and (3.6b*), (3.6c*) in which we may eliminate $u^{(1)}$ and $\boldsymbol {B}^{(1)}_x$ by means of (3.5b*), (3.5c*) and with the unknowns $(\varPhi ^{(1)}, \varPhi ^{(1)}_{,x},v^{(1)},w^{(1)},B^{(1)}_y,B^{(1)}_z)$ and $\varPi ^{(1)}$. Hence if we consider $\varPhi ^{(1)}$ and $\varPi ^{(1)}$ as given then we may solve for $\varPhi ^{(1)}_{,x}$. Thus, we obtain a first-order homogeneous linear differential equation for $\varPhi ^{(1)}$ of the form
However, the construction fails with the vanishing of the determinant of the system of five equations for $(\varPhi ^{(1)}_{,x},v^{(1)},w^{(1)},B^{(1)}_y,B^{(1)}_z)$ in terms of $\varPhi ^{(1)}$ and $\varPi ^{(1)}$. The zeros of this determinant constitute the resonant singularities of the system. We then obtain the determinant of this system by replacing $B^{(1)}_{,x}$ and $u^{(1)}$ by $\varPhi ^{(1)}$ and ignoring the terms proportional to $\varPhi ^{(1)}$ or $\varPi ^{(1)}$. We find
and after straightforward reduction we obtain
For the particular background state, resonances are associated with Alfvén waves, the first factor, and a particular form of the fast–slow wave, the second factor.
It is clear that for a general equilibrium there will be resonances when $\bar {v}=\pm \bar {f}/\sqrt {\rho }$ or $\pm \bar {f} c/\sqrt {c^2\rho + (f^2+g^2)}$. For wide ranges of the ratio of the mode indices, $m/n$ resonant surfaces will appear. Just as it was possible to find formal expansions in mode amplitudes of ideal MHD equilibrium, it may be of interest to explore the possibility of formal expansions of flows with no resonant singularities.
In a more mathematically precise language we can state the main idea as follows. One starts with a homogeneous linear system of five equations in seven unknowns of the form $A X = 0$, where $X$ is a 7-component vector and $A$ is a $5\times 7$ rectangular matrix. Next, one writes $X = (x,y)$, where $x$ denotes the first 5 components of X, and $y$ denotes the last two components, $\varPhi ^{(1)}$ and $\varPi ^{(1)}$. Then the value of $y$ is fixed, which leads to the inhomogeneous linear system $ax = b$, where $a$ is the restriction of $A$ to the subspace $y=0$ and $b = A(0,y)$. When $\text {det}(a)$ is non-zero, one therefore gets expressions for each component of $x$ that are linear in $y$. When $\text {det}(a)=0$, there is no solution for $x$ in general unless $b$ is in the image of $a$. The quantity $\varDelta$ in (3.10) is just that $\text {det}(a)$.
Such a possibility depends on arranging the state so that it is in resonance everywhere or nowhere. Note that, when $\varDelta = \text {det}(a)=0$, there is still a solution for $x$ provided $y$ can be adjusted so that $A(0,y)$ is in the image of $a$. Such a situation is possible for parallel flows, i.e. flows where $\boldsymbol {B}$ and $\boldsymbol {u}$ are parallel, or
In this case $\bar {v}=\lambda \bar {f}$, so that (3.11) becomes
If $\lambda (x)$ is chosen so that the second or the third term vanishes, then $\varDelta$ is identically zero for every mode. Such flows may well be of some interest, but they seem extremely pathological and do not appear to be realizable. However, the case $\bar {f}=0$ is of more interest and similar to the equilibrium case. In particular, we assume the second factor never vanishes and $\bar {f}=0$ for $m=M, n=N$ where $M$ and $N$ are relatively prime. Thus, we assume
where $\mu (x)$ is an arbitrary function, while $\boldsymbol {u}$ is given by (3.12), (3.14). We expand in amplitude about this state and describe expansions to all orders in the amplitude.
4. Nearly parallel flows
In the last section we showed that perturbations of steady flows are typically subject to the appearance of singularities on particular surfaces. However, for parallel flows given by (3.12), (3.14) singular surfaces need not appear. In this section we examine the flows with these properties and we lay out the argument that one may construct a formal series solution to all orders without the appearance of any singularities when the lowest-order system is a parallel flow. We do not give the full details of the proof of expansions to all orders but we consider the conclusion valid. The process closely parallels the development in Weitzner (Reference Weitzner2014), where such a formalism is fully developed.
The expansion parameter of the series is the amplitude of the flow and field components. We assume that the lowest-order field and flow state is given by (3.12), (3.14). In first order we add fields and flows which is a sum of terms given by (3.2) subject to the condition that none of them is resonant. We may then continue to construct an expansion order by order. The series we construct are assumed to have the same structure of angular dependence at all orders. The linearized system used to determine the characteristics of the system closely parallels (2.5).
In this expansion at each order we must solve a system similar to (3.3a), (3.4*), (3.5a), (3.5b*), (3.5c*), (3.6a), (3.6b*) and (3.6c*) where sums of products of lower-order terms may appear added to the right-hand side of the equations. When the mode is not resonant so that $\varDelta$ given by (3.13) is not zero, one may solve the linear homogeneous system. Our task is to indicate how one should be able to select series so that order by order there is no net resonant term arising in the combination of sums of products of lower orders which form the inhomogeneous terms of the generalization of the system.
Before we address the central issue in this work we must characterize the structure of the resonance more fully. In resonance the quantities $\bar {v}$ and $\bar {f}$ are identically zero. Thus, from (3.5b*), (3.5c*) and (3.6a*) we find that
Finally, the system is closed with the two conditions
The state given by (4.1), (4.2) is the resonant mode in any order, provided
We finally turn to the inhomogeneous system where we expand the solution to some order in the amplitude, say order $P$. We obtain a system of the form (3.3a), (3.4*), (3.5) and (3.6) where the index one is replaced by $P$ and inhomogeneous terms are added to the right-hand sides of the equations. These terms are sums of products of terms of order lower than $P$. We can solve the system provided that there is no net contribution at a resonant term, for if there were such a term we could not guarantee the construction of a periodic solution of the system. We must show that the series can be arranged so that no net inhomogeneous term is present.
We modify the structure of the solution in orders $(P-2)$ and $(P-1)$ so that no net singular terms arise in order $P$. We assume that before modification there are terms at order $P$ with $m$ and $n$ satisfying (4.3). We add a resonant term at order $(P-2)$ with angular dependence $\exp {\textrm {i}(\mu y +\nu z)}$ and undetermined $x$ dependence. There are two such independent additions $B_y^{(P-2)}$ and $v^{(P-2)}$, say. We assume that there are non-resonant terms at the first order of the form (3.2). When these terms beat against each other there will be terms of order $P-1$ with angular dependence $\exp {\textrm {i}((\mu \pm m) y +(\nu \pm n) z)}$. Finally, at order $P$, terms from order $(P-1)$ beat again against the first-order terms with angular structure $\exp {\textrm {i}(m y +n z)}$. We choose the, as yet undetermined, $x$ dependence so as to satisfy the two constraints. The quantities $u^{(P)},B_x^{(P)},\varPhi ^{(P)}$, and $\varPi ^{(P)}$ which were all zero in the solution of the homogeneous data problem are all non-zero and are determined.
5. Analysis of steady flows that admit an additional symmetry vector
We have so far outlined how one can perturbatively construct non-symmetric steady plasma flows in ideal MHD. We considered flows that are nearly parallel to avoid resonances. Compared with non-symmetric MHD equilibrium without flows, the steady flow system (2.1) is far more complicated due to extra equations and nonlinearities from the flow variables. In particular, it is not clear if the resonances can be avoided. Fortunately, analytical progress can be made if the steady flows possess a symmetry vector, the Hameiri vector (Hameiri Reference Hameiri1983), which is closely related to the quasisymmetry vector (Rodriguez, Helander & Bhattacharjee Reference Rodriguez, Helander and Bhattacharjee2020; Burby et al. Reference Burby, Kallinikos and MacKay2020a). In the following, we employ a more formal approach to investigate non-symmetric 3-D flows that come equipped with the Hameiri vector. We shall now state our goals and the main results obtained in this section.
Our primary objective is to investigate the properties of a class of non-symmetric 3-D flow system that has flow components both parallel and perpendicular to the magnetic field and admits the Hameiri's symmetry vector. Under an additional assumption that the density does not change along this vector, we show that non-symmetric generalizations of Bernoulli's law, angular momentum conservation, a GGS equation and associated generalized Hamada conditions can be obtained. Therefore, our main results show that the description of these special non-symmetric flows parallels symmetric flows.
The organization of this section is as follows. We begin with a discussion of the Hameiri vector and point out its close connections with the quasisymmetry vector in § 5.1. We then discuss the key assumptions that we make in order to carry out the analysis in § 5.2. We derive the generalizations of Bernoulli's law in § 5.3. We show that a GGS equation can be constructed following the approach of Burby et al. (Reference Burby, Kallinikos and MacKay2020a) in § 5.4. We discuss the generalization of the Hamada conditions in § 5.5 and summarize the results in § 5.6.
5.1. Hameiri's vector $\textbf {C}$ and weak quasisymmetry
We briefly summarize a key result due to Hameiri, which is essential to our work. Assuming a steady flow with nested toroidal flux surfaces labelled by $\psi$, Hameiri (Reference Hameiri1983) showed that ideal MHD allows an additional cross-helicity of the form $\int \textrm {d}^3\,\boldsymbol {r}\boldsymbol {u}\boldsymbol {\cdot }\boldsymbol {C}$, where $\textrm {d}^3\boldsymbol {r}$ denotes volume integral over 3-D space, if there exists a vector $\boldsymbol {C}$ such that
Hameiri has shown that, in axisymmetry, $\boldsymbol {C}$ is $R^2\boldsymbol {\nabla } \varphi$ in the usual $(R,Z,\varphi )$ cylindrical coordinates used in the tokamak literature. Furthermore, from (5.1a–d) we can see that the vector $\boldsymbol {C}$ is analogous to a magnetic field since it is frozen in with the flow and lies on flux surfaces.
We shall make an additional assumption that the density is constant in the direction of the Hameiri vector i.e. $\boldsymbol {C}\boldsymbol {\cdot }\boldsymbol {\nabla } \rho =0$. Writing $\boldsymbol {C}=\rho \boldsymbol {Q}$, we get an equivalent set of conditions on $\boldsymbol {Q}$,
The conditions on $\boldsymbol {Q}$ from (5.2a–d) read
while the conditions on $\boldsymbol {u}$ from (5.2a–d) and (2.1) are
We want to point out that the conditions given in (5.3) are closely related to the ‘weak quasisymmetry’ conditions obtained in Rodriguez et al. (Reference Rodriguez, Helander and Bhattacharjee2020). If instead of $\boldsymbol {Q}\boldsymbol {\cdot } \boldsymbol {\nabla } \rho =0$, we assume $\boldsymbol {Q}\boldsymbol {\cdot } \boldsymbol {\nabla } B=0$ and $\rho =\rho (\psi,B)$ we get exactly the weak quasisymmetry (QS) conditions. This close connection of ideal MHD flows with QS was pointed out earlier by Helander (Reference Helander2007, Reference Helander2014).
The analogy with ‘weak QS’ becomes more evident when we look at the internal consistency of (5.3). For given $\boldsymbol {B}$ and $\rho$, the system (5.3) is an overdetermined system for $\boldsymbol {Q}$ as it represents four constraints for the three components of $\boldsymbol {Q}$. From (5.3b) we get
The component $\left ( \boldsymbol {B}\boldsymbol {\cdot } \boldsymbol {Q}\right )$ can be chosen such that (5.3a) is satisfied. The divergence-free condition (5.3c) then imposes the following condition on $\boldsymbol {B}$
which must be satisfied in order for (5.3) to have a consistent solution. Such a consistency condition indeed appears in QS systems (Burby et al. Reference Burby, Kallinikos and MacKay2020a; Rodriguez et al. Reference Rodriguez, Helander and Bhattacharjee2020) without the last term since $\boldsymbol {Q}\boldsymbol {\cdot } \boldsymbol {\nabla } B=0$.
If we assume that $\boldsymbol {Q}$ is given and we can construct a $\boldsymbol {B}$ consistent with $\boldsymbol {Q}$. Using (5.3b), we can write $\boldsymbol {B}$ in terms of $\boldsymbol {Q}$ as
which, is analogous to the symmetry flux coordinates in a tokamak (D'haeseleer et al. Reference D'haeseleer, Hitchon, Callen and Shohet1991) and the form used by Burby et al. (Reference Burby, Kallinikos and MacKay2020a) to study the strong form of QS. To ensure that $\boldsymbol {\nabla }\boldsymbol {\cdot }\boldsymbol {B}=0$ is satisfied by (5.7), we must have
In the following, we shall assume that there exists a Hameiri vector of the form $\boldsymbol {C}=\rho \boldsymbol {Q}$, where $\boldsymbol {C}\boldsymbol {\cdot } \boldsymbol {\nabla } \rho =0$ and $\boldsymbol {Q}$ satisfies (5.3) and (5.8) and is given.
5.2. Underlying assumptions
Before proceeding further, we shall discuss the various assumptions that we make in our analysis. We begin with our fundamental assumption that there exist non-symmetric plasma flows that possess nested flux surfaces. In our formal exploration, we are not forced to make any specific assumptions about the rotational transform of the magnetic field or the magnetic shear. The nestedness of flux surfaces is not guaranteed but is often made in the literature. In the following, we shall derive the necessary conditions (the Hamada conditions) required to suppress current singularities on the rational surfaces and preserve the nested surfaces. It will be apparent that the Hamada conditions cannot possibly be satisfied in a toroidal volume unless magnetic shear is weak. With the analysis of previous sections in mind, we now make a conscious choice of weak magnetic shear in the following.
Our second important assumption is the existence of Hameiri's symmetry vector $\boldsymbol {C}$. Theoretically, $\boldsymbol {C}$ can be obtained by applying Noether's theorem to the one-fluid ideal MHD Lagrangian (Ilgisonis & Pastukhov Reference Ilgisonis and Pastukhov2000). In particular, $\boldsymbol {C}$ is related to the relabelling symmetry of ideal MHD (Vladimirov & Moffatt Reference Vladimirov and Moffatt1995; Hameiri Reference Hameiri1998; Ilgisonis & Pastukhov Reference Ilgisonis and Pastukhov2000). However, the system of equations determining $\boldsymbol {C}$, (5.1a–d), being a nonlinear overdetermined system, makes it difficult to prove the existence of such a vector. Moreover, unlike QS, a systematic study of $\boldsymbol {C}$ has not been carried out to the best of our knowledge.
As discussed in detail in Hameiri (Reference Hameiri1998), an important necessary condition that needs to be satisfied in order for $\boldsymbol {C}$ to exist in a toroidal domain (Hameiri Reference Hameiri1983) is that, on each closed magnetic field line,
Here, $m(\psi )$ is assumed to be a smooth function and $\textrm {d}\ell$ denotes the differential element along the magnetic field. If condition (5.9) is not satisfied then the flow must be nearly parallel (Hameiri Reference Hameiri1998).
We then assumed that $\boldsymbol {C}\boldsymbol {\cdot }\boldsymbol {\nabla }\rho =0$, which implies that the density (and hence pressure from (2.2)) possesses a continuous symmetry along $\boldsymbol {C}$. This choice has been made solely for simplicity and analytical tractability. Let us briefly discuss some of the consequences of the assumption of $\boldsymbol {Q}\boldsymbol {\cdot } \boldsymbol {\nabla } \rho =0=\boldsymbol {Q}\boldsymbol {\cdot } \boldsymbol {\nabla } p$. Firstly, we note that (5.9) is satisfied identically. Secondly, if $\boldsymbol {Q}$ lines do not close on themselves then density and pressure must not vary on a flux surface i.e. $p=p(\psi ),\rho =\rho (\psi )$. Therefore, such systems must be subsonic (Tasso & Throumoulopoulos Reference Tasso and Throumoulopoulos1998) since large centrifugal forces lead to variations of density and pressure on flux surfaces. When $\boldsymbol {Q}$ lines are closed, flows do not have to be subsonic, and pressure variations on flux surfaces are allowed. Although we shall not make any assumption on the closure of $\boldsymbol {Q}$ lines, it is interesting to note that, in QS, the symmetry lines close helically or toroidally, depending on the type of QS.
The final assumption that $\boldsymbol {Q}$ satisfies (5.8), is due to the fact that (5.7) is divergence free only in axisymmetry (Burby et al. Reference Burby, Kallinikos and MacKay2020a). The condition (5.8) is not really an additional constraint on $\boldsymbol {Q}$ (or $\boldsymbol {C}=\rho \boldsymbol {Q}$), since it follows directly from (5.3b). Any $\boldsymbol {Q}$ satisfying (5.3b) with a divergence-free physical magnetic field necessarily satisfies (5.8). We have included the condition (5.8) only because of its usage in subsequent calculations.
5.3. Generalized Bernoulli law and angular momentum
It is well known (Hameiri Reference Hameiri1983; Throumoulopoulos et al. Reference Throumoulopoulos, Weitzner and Tasso2006) that five scalar functions are needed to describe the steady axisymmetric ideal MHD flows. Out of these, two functions are given by the electrostatic potential $\varPhi '(\psi )$ and the entropy (assumed constant). The remaining three functions, denoted by $(\varLambda (\psi ), H(\psi ), I(\psi ))$, are associated with the parallel component of the flow, the Bernoulli law and angular momentum, respectively. In the following, our goal is to derive these quantities systematically for non-symmetric flows with the help of the Hameiri vector $\boldsymbol {C}=\rho \boldsymbol {Q}$.
For a given magnetic field $\boldsymbol {B}$ and $\boldsymbol {Q}$ that satisfy (5.3) and (5.6), it can be shown (Hameiri Reference Hameiri1998) that a steady state flow satisfying (5.4) is given by
which clearly resembles the form of the axisymmetric steady flow of a tokamak (Hameiri Reference Hameiri1983). We also get the flux function $\varLambda (\psi )$ from (5.10).
From (5.10), we note that any two vectors from $\boldsymbol {u}$, $\boldsymbol {B}$ and $\boldsymbol {Q}$ can be used as basis vectors on a flux surface. We shall choose $\boldsymbol {B}$ and $\boldsymbol {Q}$ since the commutator of the two derivatives $\boldsymbol {B}\boldsymbol {\cdot }\boldsymbol {\nabla },\boldsymbol {Q}\boldsymbol {\cdot } \boldsymbol {\nabla }$ vanish (as shown in Appendix B). Choosing $\boldsymbol {B},\boldsymbol {Q},\boldsymbol {\nabla }\psi$ as basis vectors, we shall now obtain the components of the ideal MHD force balance equation (2.1b).
We find it convenient to rewrite (2.1b) in the form of a vorticity equation
Here, we have defined $h$ such that $\boldsymbol {\nabla } h= (1/\rho )\boldsymbol {\nabla } p(\rho )$, and we have used the standard vector identity $\boldsymbol {u}\boldsymbol {\cdot }\boldsymbol {\nabla } u = \boldsymbol {\nabla } (u^2/2) -\boldsymbol {u}\times (\boldsymbol {\nabla }\times \boldsymbol {u})$. Dotting with $\boldsymbol {B}$ and $\boldsymbol {Q}$ respectively, and using (5.4b), (5.10) we get
To evaluate the divergence terms that appear in (5.12) we use (5.4b) and (B2). We note that the assumption $\boldsymbol {Q}\boldsymbol {\cdot } \boldsymbol {\nabla } \rho =0$ helps us to simplify (5.12b).
Simplifying (5.12) using (5.3a), (5.6) and (B2) we get
In axisymmetry, the $\boldsymbol {Q}\boldsymbol {\cdot } \boldsymbol {\nabla }$ terms in (5.13) vanish identically, and we get two homogeneous magnetic differential equations. Solving the homogeneous magnetic differential equations, we get two flux functions which denote Bernoulli's law and momentum conservation in the symmetry ($\boldsymbol {Q}$) direction (Hameiri Reference Hameiri1983; Tasso & Throumoulopoulos Reference Tasso and Throumoulopoulos1998). The situation is similar here except that the magnetic differential equations are non-homogeneous; therefore, we need to check for consistency conditions. We postpone the discussion on the consistency conditions until § 5.5.
Assuming that the solvability conditions are satisfied, we solve (5.13) and obtain the generalized Bernoulli and the generalized angular momentum equations
where, $\mathcal {H},\mathcal {I}$ defined by
are single-valued and doubly periodic functions if Hamada conditions (see (5.24) below) are satisfied. They vanish identically in the axisymmetric case. As is well known (Hameiri Reference Hameiri1983, Reference Hameiri1998; Beskin Reference Beskin2009) availability of conserved quantities is very helpful in characterizing flow patterns.
5.4. GGS equation for flows with the Hameiri vector
We are now in a position to obtain the GGS from the $\boldsymbol {\nabla } \psi \boldsymbol {\cdot }$ component of the vorticity equation (5.11). The step-by-step derivation is provided in Appendix D. The steps to derive the GGS are identical to the ones in the derivation of the flow-modified classical GS equation in axisymmetry or helical symmetry. However, there are important differences which we shall now discuss. First, in the absence of axial/helical symmetry, the vector
which denotes the deviation of $\boldsymbol {Q}$ from the axial/helical symmetry vector (Burby et al. Reference Burby, Kallinikos and MacKay2020a) appears in the GGS. We can check that (5.8), the condition required to ensure $\boldsymbol {\nabla }\boldsymbol {\cdot }\boldsymbol {B}=0$ when $\boldsymbol {B}$ is expressed in terms of $\boldsymbol {Q}$ (5.7), is identical to
Obviously, if $\boldsymbol {Q}$ is an axial/helical symmetry both terms individually vanish and $\boldsymbol {B}$ is automatically divergence free. Secondly, the geometric quantity
which does not appear in the axisymmetric GS but appears in the helically symmetric GS, also appears in the GGS. Finally, terms appear with the force-like quantity
which vanishes identically in static MHD (Burby et al. Reference Burby, Kallinikos and MacKay2020a).
With the definitions (5.16), (5.18) and (5.19) in mind, we can write the flow-modified GGS in the following form
where,
In axisymmetry (Hameiri Reference Hameiri1983), $\mathcal {N}=0$. For helical symmetry, an additional term of the form $(I\omega _Q)/Q^4$ appears. To study a force balance other than ideal MHD force balance, we can replace $\boldsymbol {F}_{\boldsymbol {Q}}$ by other forces (Rodriguez, Sengupta & Bhattacharjee Reference Rodriguez, Sengupta and Bhattacharjee2021). For the case of flows only in the symmetry direction $\varLambda (\psi )=0$ and for strictly parallel flows $\varPhi '(\psi )=0$. One can easily extend the GGS to include other non-relativistic forces.
Unlike symmetric geometries, the coefficients and derivatives in the flow-modified GGS depend on all three spatial variables (Constantin et al. Reference Constantin, Drivas and Ginsberg2021; Burby et al. Reference Burby, Kallinikos and MacKay2020a). Hence, we need to impose the additional condition
Whether such solutions to the GGS can be constructed is a challenging open problem.
5.5. The generalized Hamada conditions
We now return to the consistency conditions that must be satisfied in order for (5.13) to have non-singular smooth solutions. When the rotational transform is irrational, we can flux-surface average (5.13). Since both $\boldsymbol {B}$ and $\boldsymbol {C}$ are tangential to the flux surfaces, there is no inconsistency. On the other hand, for rational rotational transform, the vanishing of the two equations under closed field line $\oint \textrm {d}\ell /B$ integral is not automatic and leads to non-trivial constraints (Newcomb Reference Newcomb1959) – the so-called Hamada conditions (Hamada Reference Hamada1962; Helander Reference Helander2014).
We recall that in ideal MHD equilibrium with scalar pressure $p(\psi )$, the integral constraints that must be satisfied on rational surfaces are (Grad Reference Grad1967; Solov'ev & Shafranov Reference Solov'ev and Shafranov1970; Grad Reference Grad1971)
where, $c_1(\psi ),c_2(\psi )$ are single-valued continuous functions of $\psi$ that reduce to constants in the vacuum limit. We give a straightforward proof of these conditions in Appendix C.
Hamada conditions (5.23a–c) are identically satisfied in symmetric geometries but are not satisfied in general for non-symmetric geometry with arbitrary pressure and rotational transform profiles (Boozer Reference Boozer1981; Weitzner Reference Weitzner2016), leading to singular currents on rational surfaces (Loizu et al. Reference Loizu, Hudson, Bhattacharjee and Helander2015a,Reference Loizu, Hudson, Bhattacharjee, Lazerson and Helanderb). In the following, we shall obtain the generalized Hamada conditions for steady flows with Hameiri's vector $\boldsymbol {C}=\rho \boldsymbol {Q}$.
We note that (5.9) is already in the form of a Hamada condition. To find the other conditions, we carry out the $\oint \textrm {d}\ell /B$ integrals of (5.13) (denoted by $\langle \ \rangle$ brackets). Using the commutation of the two derivatives to get two homogeneous equations of the form $\boldsymbol {Q}\boldsymbol {\cdot } \boldsymbol {\nabla } \langle A \rangle =0$ implies that $A$ must be a flux function. Thus, we obtain two consistency conditions
which show that the circulations of $\boldsymbol {u}$ and $\boldsymbol {B}$ along the closed field line are required to be flux functions in order that (5.13) be solvable.
The final Hamada condition is obtained by closed field line averaging of (5.20) on a rational flux surface, which leads to the following equation for $C_2'(\psi )$ (see Appendix E)
Equations (5.24) and (5.25), together with (5.9) describe integral constraints that need to be satisfied on each rational flux surface. These are the generalized Hamada conditions that include the effect of steady flows. In Appendix C we show that when these consistency conditions are satisfied magnetic resonances on rational surfaces are avoided. Given the complicated nature of the generalized Hamada conditions, it is doubtful that they will be satisfied by a generic steady flow. Furthermore, since the Hamada conditions need to be satisfied on each rational surface, the magnetic shear needs to be sufficiently weak to avoid passage through multiple low-order rational surfaces.
5.6. Summary of § 5
In § 5, we discussed a special class of non-symmetric steady flows endowed with Hameiri's symmetry vector $\boldsymbol {C}=\rho \boldsymbol {Q}$, where $\boldsymbol {C}\boldsymbol {\cdot } \boldsymbol {\nabla } \rho =0$ and $\boldsymbol {Q}$ satisfy the overdetermined system (5.3). If the integral constraints (5.9), (5.24) and (5.25), are satisfied then the description of such flows involves a GGS (5.20) and four flux functions similar to axisymmetric flows (Hameiri Reference Hameiri1983; Tasso & Throumoulopoulos Reference Tasso and Throumoulopoulos1998) namely, $\varLambda (\psi ),\varPhi '(\psi ),H(\psi ),I(\psi )$. Note that only the first four appear in the GGS since we have assumed a barotropic equation of state with constant entropy. The flux function $m(\psi )$ is needed to ensure the existence of Hameiri's symmetry vector. The integral constraints (5.24) are necessary to avoid current singularities (see Appendix C). They furnish two more flux functions $C_1(\psi ), C_2(\psi )$ which measure the circulations of $\boldsymbol {u}$ and $\boldsymbol {B}$ in a closed field line system. Physically, they amount to the constraint that the flows should not lead to accumulation of charges inside the closed magnetic flux tubes on rational surfaces. Finally, non-symmetric geometry introduces several extra terms grouped together as $\mathcal {N}$ such as $\mathcal {H},\mathcal {I}, \boldsymbol {w}\boldsymbol {\cdot } \boldsymbol {\nabla } \psi$ etc.
6. Discussion
We studied several classes of non-symmetric ideal MHD equilibrium with flows larger than the diamagnetic flows. Such flows could occur in stellarators with quasisymmetry. We showed that the techniques developed in carrying out perturbation methods to all orders for static MHD equilibrium (Weitzner Reference Weitzner2014) can also be extended to MHD equilibrium with flows. The basic idea is that if the field lines are closed, i.e. the rotational transform is rational and the magnetic shear is weak, one can systematically eliminate magnetic resonances at each order by utilizing the ‘free function’ that one gets by solving the magnetic differential equation from lower orders. Plasma flows introduce extra resonances, which are absent in static MHD. It is, in general, very complicated to eliminate all such resonances. However, for the class of equilibrium with nearly parallel flows, one can do so. It is to be noted that exact solutions with large parallel flows can be constructed (Kamchatnov Reference Kamchatnov1982) when the flow is considered incompressible. We have argued that compressible, nearly parallel flows are also possible in the MHD model of plasmas.
Although the present work on the analysis of magnetic resonances has been carried out in rectangular coordinates with periodic boundary conditions, extensions to polar and toroidal coordinates should be fairly direct. The principal issue is behaviour near the magnetic axis. To deal with this complication, one must require that components of the fields and flows vanish sufficiently rapidly near the axis. The interactions through the nonlinear terms typically do not destroy these properties. We leave the detailed near-axis analysis for the future.
We have studied a particular class of perpendicular flows by taking a more formal approach. Our approach is based on utilizing a symmetry vector found by Hameiri (Hameiri Reference Hameiri1981, Reference Hameiri1998). The generalizations of Bernoulli's law and conservation of angular momentum were also obtained thanks to Hameiri's symmetry vector. We obtain a flow-modified GGS equation for such flows and the constraints resulting from the closed magnetic field lines. It is clear from the nature of the constraints that if such non-parallel flows exist, they must be special. There exists a close connection between Hameiri's vector and the QS vector that will be further discussed elsewhere (Rodriguez et al. Reference Rodriguez, Sengupta and Bhattacharjee2021).
We have not attempted to estimate the damping of the flows due to neoclassical effects, which is an essential but difficult question (Simakov & Helander Reference Simakov and Helander2009). We shall address this in detail in a forthcoming paper (Rodriguez et al. Reference Rodriguez, Sengupta and Bhattacharjee2021). We only make a few observations here. First, the neoclassical effects are tied to the parallel electric field, which does not appear in ideal MHD. If the parallel electric field can be shown to be self-consistently small or higher order in the Larmor radius, then the MHD description prevails, and flow damping would be minimized. Therefore, our results are valid within the framework of ideal MHD. Second, the techniques developed here can be easily extended to more sophisticated MHD models that are relevant to astrophysical plasmas.
Our results here do not necessarily contradict earlier results (Simakov & Helander Reference Simakov and Helander2009; Sugama et al. Reference Sugama, Watanabe, Nunami and Nishimura2011) that argue that large flows cannot be supported even if the stellarator is quasisymmetric but not completely axisymmetric. We show that it is challenging to avoid magnetic resonances unless the flow is nearly parallel. Except for very special flows, the generalized Hamada conditions on rational surfaces cannot be satisfied, severely limiting the possibility of steady non-symmetric flows in generic stellarators. We expect that kinetic constraints will further restrict the class of available solutions.
Acknowledgements
The authors would like to sincerely thank the anonymous reviewer for their excellent, meticulous and constructive suggestions that have significantly improved the organization and presentation of the results in the paper. In addition, W.S. would like to thank E. Rodriguez, A. Bhattacharjee, A.B. Hassam, E.J. Paul and J. Juno for helpful discussions.
Editor William Dorland thanks the referees for their advice in evaluating this article.
Funding
This research was partly funded by the US DOE grant no. DEFG02-86ER53223 and Simons Foundation/SFARI (560651, AB).
Declaration of interests
The authors report no conflict of interest.
Appendix A. Analysis of the dispersion relation (2.7)
Our goal here is to show that real solutions of the dispersion relation (2.7) exist. We shall rewrite the dispersion relation in terms of dimensionless physical quantities and analyse the parametric space of these parameters, which supports real solutions of (2.7). We shall also look at some simple, special cases which will allow us to gain insight into the parametric space.
Let us consider a generic flow of the form
Associated with this flow, we can define the parallel and perpendicular Mach numbers as follows
$c$ being the sound speed.
As noted in the text, the dispersion relation is independent of the magnitude of $\boldsymbol {k}$. Therefore, we can divide (2.7) by $(\rho c^2 k^2)^2$ to get
where,
We will not consider any components of ${\hat {\boldsymbol {k}}}$ orthogonal to $\boldsymbol {u}$ and $\boldsymbol {B}$ since they trivially satisfy (A3). Hence, we will assume that the unit vector $\hat {k}$ is coplanar with $\boldsymbol {u},\boldsymbol {B}$ and therefore, admits the following orthogonal decomposition
Furthermore, (A4a–c) together with (A2a,b) implies that
To determine $\theta$ we substitute (A6a,b) into the dispersion relation (A3). which leads to the following quartic equation for $\xi =\tan \theta$
We note that there are only three physical parameters that come into the problem: plasma beta, $\beta =\rho c^2/B^2$, parallel and perpendicular Mach numbers $(M_\parallel, M_\perp )$ given by (A2a,b). In the following, we shall assume that all three physical quantities are real and positive.
The theory of quartic equations is well developed, and the necessary conditions for a quartic equation to have real roots can be determined straightforwardly (Rees Reference Rees1922). When the quartic discriminant is negative, the equation has two real and two imaginary roots. We have plotted the associated parametric region in figures 1(a) and 1(b). The region where four real roots can exist is shown in figure 1(c). These two regions are complimentary and hence the equation always has a real root. The details are not particularly illuminating and will be omitted here. Fortunately, a lot of insight can be gained by studying some simple cases, which we shall do next.
We begin with the special case of purely parallel flows i.e. $M_\perp =0$. In this case, the equation simplifies considerably and we get
As $M_\parallel ^2\to 1/(\beta +1), \xi \to \infty$ and the angle $\theta =\arctan {\xi }\to {\rm \pi}/2$, which implies that $\boldsymbol {k}$ becomes orthogonal to $\hat {B}$. At the other extreme, when $M_\parallel ^2=1$ or $M_\parallel ^2=1/\beta$, $\xi =0$, i.e. $\boldsymbol {k}$ is parallel to $\hat {B}$. The case $M_\parallel ^2=1/\beta$ is discussed in Kamchatnov (Reference Kamchatnov1982). These special solutions appear as lines of discontinuity in figures 1(b) and 1(c). The quartic determinant is zero along these curves leading to multiple roots. The solutions corresponding to $M_\parallel ^2=1/(\beta +1),1/\beta$ for flows of the type (3.12) were described as pathological solutions in (3.13) and lead to $\varDelta =0$. Physically, the system can develop shocks along these curves.
Omitting the three special cases discussed above, we find that real solutions with $\xi >0$ exist for $M_\parallel$ less (greater) than one, provided $M_\parallel ^2$ lies inside (outside) of the interval $(1/(\beta +1),1/\beta )$. We note that for large values of $\beta$ subsonic parallel ($M_\parallel ^2<1$) flows are rare compared with the supersonic parallel flows, in accordance with figure 1(a).
Next, we consider the case of purely perpendicular flows, i.e. $M_\parallel =0$. If
Since $\xi ^2$ is positive definite, $\xi$ is real and the flow is perpendicular and supersonic.
If $M_\parallel =0$ but $M_\perp$ does not satisfy (A9a–c), (A7) reduces to the following bi-quadratic equation
For the roots to be real, the discriminant of (A10) must be positive, which implies
We find that $\beta =1$ represents a discontinuity in the solution, as depicted in figures 1(a) and 1(b). For values of $\beta$ close to zero we find that perpendicular flows can exist only if $M_\perp \approx 1$. For large $\beta$ we find that the first condition is easily satisfied but not the second one.
Furthermore, for $\xi$ to have four real roots we need
On the other hand, when $M_\perp ^2<1+1/\beta$, $\xi$ has at least two real roots provided (A11) is satisfied. Therefore, both supersonic and subsonic perpendicular flows can exist.
Finally, we note that ${\hat {\boldsymbol {k}}},\hat {\boldsymbol {B}}, {\hat {\boldsymbol {u}}_\perp }$ can be plotted on a unit sphere, as shown in figure 2. We choose $\hat {\boldsymbol {B}}$ to be in the $z$ direction and, therefore, $\boldsymbol {\hat {u}_\perp }$ must lie in the $x\text {--}y$ plane. The coplanar vector ${\hat {\boldsymbol {k}}}$ makes an angle $\theta$ with $\boldsymbol {B}/B$, where $\theta =\arctan {\xi }$. Since the equation for $\xi$, (A7), depends on $\boldsymbol {u}$ only through the combinations $M_\parallel \propto \boldsymbol {u}\boldsymbol {\cdot } \boldsymbol {B}$ and $M_\perp \propto |\boldsymbol {u}_\perp |$, it is independent of the direction of $\hat {\boldsymbol {u}}_\perp$. As a result, the curve that $\boldsymbol {k}$ traces out is symmetric around $\hat {B}$ and hence a circle on the unit sphere.
Appendix B. Useful identities and expressions
We collect some useful identities and expressions in this appendix that have been used in deriving some of the key results in the text.
A useful expression is that of the commutator $[\boldsymbol {B}\boldsymbol {\cdot }\boldsymbol {\nabla },\boldsymbol {Q}\boldsymbol {\cdot } \boldsymbol {\nabla } ]$. Employing standard Einstein summation convention, $\boldsymbol {B}\boldsymbol {\cdot }\boldsymbol {\nabla }=B_i\partial _i$, $\boldsymbol {Q}\boldsymbol {\cdot } \boldsymbol {\nabla } =Q_j\partial _j$, together with (5.4b) implies that $\epsilon _{klm}B_l Q_m = \partial _k \psi$. It is straightforward to show that
The following identity is useful to simplify the divergence of $\boldsymbol {u}\times \boldsymbol {\nabla }\psi$
Appendix C. Closed field line consistency conditions and current singularities
The close relation between current singularities that can develop on rational surfaces in a non-symmetric toroidal domain in ideal MHD, and the necessary conditions in the form of integrals over closed magnetic field lines such as $\oint \textrm {d}\ell /B = f(\psi )$ are well known (Hamada Reference Hamada1962; Grad Reference Grad1967; Solov'ev & Shafranov Reference Solov'ev and Shafranov1970; Boozer Reference Boozer1981; Loizu et al. Reference Loizu, Hudson, Bhattacharjee and Helander2015a). We show here that, to avoid current singularities on a rational surface, the integral constraints of the form (5.23a–c) for ideal MHD equilibrium, and (5.24), (5.25) for steady flows, are sufficient conditions.
We shall use a generic $(\psi,\alpha,\varphi )$ coordinate system, where $\psi$ is the flux label, $\alpha$ is the field line label and $\varphi$ denotes another angle like coordinate. We will assume that the magnetic field lines are closed and shear is weak such that both $\boldsymbol {\nabla } \alpha$ and $\boldsymbol {\nabla } \varphi$ are single valued. A more general derivation for finite shear can be carried out using generalized Boozer coordinates.
In the $(\psi,\alpha,\varphi )$ coordinates, we can represent $\boldsymbol {B}$ in the following covariant and contravariant forms
The functions $B_\psi,B_\alpha, B_\varphi$ are all single-valued functions since $\boldsymbol {B}$ is single valued and so are the gradients of $(\psi,\alpha,\varphi )$. The Jacobian $\mathcal {J}^{-1}= \boldsymbol {\nabla }\psi \times \boldsymbol {\nabla } \alpha \boldsymbol {\cdot } \boldsymbol {\nabla } \varphi = \boldsymbol {B}\boldsymbol {\cdot }\boldsymbol {\nabla } \varphi$ connects $B_\varphi$ to $B^2$ through
The current obtained by taking the curl of the covariant form of $\boldsymbol {B}$ satisfies
We now consider a general force-balance equation of the form
Note that $\boldsymbol {F}$ does not have a $\boldsymbol {\nabla } \varphi$ component because of parallel force balance condition $\boldsymbol {F}\cdot \boldsymbol {B}=0$.
If the current $\boldsymbol {J}$ simultaneously satisfies $\boldsymbol {J}=\boldsymbol {\nabla }\times \boldsymbol {B}$ and $\boldsymbol {J}\times \boldsymbol {B}=\boldsymbol {F}$, we can equate $\boldsymbol {J}\times \boldsymbol {B}$ from (C5a,b) to (C4) to get the following magnetic differential equations for $B_\psi$ and $B_\alpha$ that enforce force balance
Eliminating $B_\varphi$ we obtain another magnetic differential equation
Since $B_\psi$ and $B_\alpha$ are single valued, the well-known necessary and sufficient conditions for the solvability of magnetic differential equations (Newcomb Reference Newcomb1959) lead to the following constraints
Using (C2) and $\textrm {d}\ell /B = \mathcal {J} \,\textrm {d}\varphi$ we see that these constraints can be written as
We note that the last consistency condition is not independent of the first two since it is obtained by eliminating $B_\phi$ from the first two conditions. However, once the integral constraints (C9) are satisfied, we can solve the Magnetic differential equation (MDE) for $B_\psi$ and $B_\alpha$. Therefore, the constraints (C9) are necessary and sufficient for force balance.
For ideal MHD equilibrium $\boldsymbol {F}=p'(\psi )\boldsymbol {\nabla }\psi$ and we recover (5.23a–c) from (C9). For steady flows both $F_\psi$ and $F_\alpha$ are non-zero and we recover (5.24b) and (5.25).
To show the sufficiency of the integral conditions, we start with the solution of the force-balance condition and impose the divergence-free condition on the current. Writing $\boldsymbol {J}=\boldsymbol {J}_\perp + \boldsymbol {B} (j_{||}/B)$, and taking the divergence, we obtain a magnetic differential equation for the parallel component of current,
Therefore, if the Newcomb condition
is not satisfied, then current singularity follows (Loizu et al. Reference Loizu, Hudson, Bhattacharjee and Helander2015a).
Using (C5a,b) and (C1a,b) we can obtain an expression for the perpendicular component of the current
where $\boldsymbol {e}_\psi = \mathcal {J} \boldsymbol {\nabla }\alpha \times \boldsymbol {\nabla } \varphi$ etc. are the basis vectors. The divergence of $\boldsymbol {J}_\perp$ is
If the integral constraints (C9) are satisfied, (C7) is solvable and we can rewrite (C13) as
which satisfies (C11). The MDE (C10) can then be solved for the Pfirsch–Schlüter currents.
Therefore, the integral constraints (C9) are sufficient for the current $\boldsymbol {J}$ to satisfy force balance and be divergence free.
Appendix D. Derivation of the flow-modified GGS
Since the steps in deriving the flow-modified GGS are a bit tricky and one can easily get lost in the plethora of terms, we give a step-by-step derivation in this Appendix. Starting with the flow $\boldsymbol {u}$ given by (5.10) we obtain the following expression for the vorticity $\boldsymbol {\omega }$ and $\boldsymbol {u}\times \boldsymbol {\omega }$
With the help of (D1) we obtain the right-hand side of $\boldsymbol {\nabla }\psi$ dotted with the vorticity equation (5.11)
Now, using (5.14a), the left side of $\boldsymbol {\nabla }\psi$ dotted with the vorticity equation yields
We can see that the $\varPhi ''(\psi )$ term cancels when we equate the left (D3) to the right side (D2).
Let us first evaluate the first term on the right-hand side of (D2) by replacing $\boldsymbol {\nabla } \psi$ by $\boldsymbol {B}\times \boldsymbol {Q}$ such that
The components of $\boldsymbol {J}$ can be calculated as follows
Substituting (D5) in (D4) and using
we get
where,
In simplifying the last term of (D2) we use instead of $\boldsymbol {u}\times \boldsymbol {\nabla } \psi$
We now rearrange all the terms in the $\boldsymbol {\nabla }\psi \cdot$ component of the vorticity equation in the form
where $Y_3$ is given by
with $C=\left ( \boldsymbol {B}\boldsymbol {\cdot } \boldsymbol {Q}\right )$. We simplify $Y_3$ by using $\boldsymbol {\nabla } (\varLambda ^2/\rho )=\varLambda \boldsymbol {\nabla } (\varLambda /\rho )+(\varLambda /\rho )\boldsymbol {\nabla } \varLambda$ and (D6a,b)
Substituting (D12) into $Y_3$ and simplifying we get
The final step is to use the generalized momentum equation (5.14b) to replace the term, $C\boldsymbol {\nabla } C$ in (D13) by $C\boldsymbol {\nabla } \left ( I(\psi )+\mathcal {I} +\varLambda \left ( \boldsymbol {u}\boldsymbol {\cdot } \boldsymbol {Q} \right )\right )$. After some straightforward algebra we get
Appendix E. Derivation of the final Hamada condition
To derive the final Hamada condition (5.25), it is convenient to set up flux coordinates $(\psi,\alpha,\beta )$ (discussed in Hameiri Reference Hameiri1998) such that
For closed field line systems with zero magnetic shear, $\boldsymbol {\nabla }\alpha$ is a single-valued quantity. Since, $\boldsymbol {B}$ and $\boldsymbol {Q}$ must be single valued, $\boldsymbol {\nabla }\beta$ in general is either single valued or at least the multi-valued part of $\boldsymbol {\nabla }\beta$ is in the $\boldsymbol {\nabla }\psi$ direction (Hameiri Reference Hameiri1998). For closed $\boldsymbol {Q}$ lines, $\beta$ is single valued.
Although we use the $(\psi,\alpha,\beta )$ coordinates to derive (5.25), the result is independent of the choice and details of the flux coordinates. In particular, any flux coordinates that allow radial currents, e.g. generalized Boozer coordinates (Rodriguez et al. Reference Rodriguez, Helander and Bhattacharjee2020), can be used.
The condition $\boldsymbol {B}\times \boldsymbol {Q}=\boldsymbol {\nabla } \psi$ implies that the Jacobian $\mathcal {J}^{-1}=\boldsymbol {\nabla }\psi \times \boldsymbol {\nabla }\alpha \boldsymbol {\cdot } \boldsymbol {\nabla } \beta$ is unity. Therefore,
In these coordinates, the divergence of a vector of the form $f\boldsymbol {\nabla } \psi$ is given by
Averaging (E3) along the closed magnetic field line we get
where, $[\beta ]$ denotes the jump in $\beta$ along the closed field line. If $\boldsymbol {Q}$ lines are closed $[\beta ]=0$.
We obtain (5.25) for $f=(1-\varLambda ^2/\rho )/Q^2$ upon using (5.14), (5.24) and the identity