text:
# Grad-Shafranov equation in noncircular stationary axisymmetric spacetimes

## I Introduction

While most neutron stars have magnetic fields of $`10^{12}`$-$`10^{13}`$ G, studies of soft gamma-ray repeaters and anomalous X-ray pulsars suggest that a significant fraction $`(10\%)`$ of neutron stars is born with larger magnetic fields $`10^{14}`$-$`10^{15}`$ G k98 ; dt92 ; m02 ; t01 . The internal magnetic field of a new born neutron star may be even larger $`10^{16}`$ G if it is generated by the helical dynamo dt92 ; td93 . A magnetic field of nearly maximum strength allowed by the virial theorem $`10^{17}`$ G may be achieved if the central engine of gamma-ray bursts are neutron stars u92 ; n98 ; kr98 ; wyhw00 . In such super-magnetized neutron stars, so-called magnetars dt92 ; t01 , the magnetic fields have substantial effects on their internal stellar structure. Especially the deformation due to the magnetic stress becomes non-negligible b95 ; bg96 ; i01 ; kok99 . Since the deformation affects the precession, oscillations and the gravitational wave emission of neutron stars i01 ; bg96 ; m99 ; c02 ; it00 , it is important to investigate equilibrium configurations of magnetars.

General relativistic effects are sizable in the interior of a neutron star, so that any quantitative investigation of the magnetars has to be based on general relativistic magnetohydrodynamics (MHD) l67 ; nt73 ; bo78 ; bo79 . Therefore, we have to solve the matter and electromagnetic field configurations in a curved spacetime, and have to take account of the electromagnetic energy-momentum as a source of the gravitational field. So far several works have been devoted to equilibrium configurations of a magnetized star in a stationary axisymmetric spacetime b93 ; b95 ; kok99 ; cpl01 . However these works consider only poloidal magnetic fields for simplicity, since the existence of only a poloidal field is compatible with the circularity of the spacetime gb93 ; o02 . In a circular spacetime, there exists a family of two-surfaces everywhere orthogonal to the plane defined by the two Killing vectors associated with stationarity $`\eta ^\mu =(/t)^\mu `$ and axisymmetry $`\xi ^\mu =(/\phi )^\mu `$ p66 ; c69 ; c72 . Thus one may choose the coordinates $`\left(x^\mu \right)=(t,x^1,x^2,\phi )`$ such that the metric components $`g_{01}`$, $`g_{02}`$, $`g_{31}`$ and $`g_{32}`$ are identically zero. As a consequence, the problem is simplified dramatically.

However non-negligible toroidal magnetic fields are likely to exist in nature. Differential rotation generated during the gravitational collapse s03 ; s00  or in the binary coalescence su00  may wind up the frozen-in magnetic field to amplify the toroidal component kr98 ; wyhw00 . A toroidal magnetic field may be generated by the $`\alpha `$-$`\mathrm{\Omega }`$ dynamo during the first few seconds after the formation of a millisecond pulsar dt92 ; td93 . In addition, convective motion may also exist in the interior of a neutron star meas76 , which also violates the circularity of the spacetime gb93 ; o02 . Thus, we have to consider noncircular spacetimes. The degree of noncircularity of the spacetime in a neutron star with mass $`M_{}`$ and radius $`R_{}`$, will be about $`(M_{}/R_{})v_{\mathrm{mf}}`$ and $`0.1`$-$`0.01(M_{}/R_{})_M`$, where $`v_{\mathrm{mf}}`$ is the velocity of the meridional flow and $`_M`$ is the ratio of the magnetic energy to the gravitational energy Ioka2 .

The problem to obtain an equilibrium configuration of a magnetar can be separated into two parts. The first part is the Einstein equations which determine the spacetime geometry under a given configuration of matter and electromagnetic fields. The second part is the matter and electromagnetic field equations in a given spacetime geometry. A $`(2+1)+1`$ formalism to solve the Einstein equations under the presence of a spatial Killing vector was developed by Maeda et al. m80 ; nok87  and by Sasaki s84 . This formalism is similar to the well-known $`3+1`$ formalism sy78 . Later Gourgoulhon and Bonazzola gb93  developed a similar but different $`(2+1)+1`$ formalism which is more suited for stationary axisymmetric spacetimes. So here we focus our attention on the second problem, i.e., to formulate the equations of motion of matter and electromagnetic fields in a curved spacetime.

It is well known that the basic equations for a stationary axisymmetric ideal MHD system can be reduced to a single second-order, nonlinear partial differential equation, the so-called Grad-Shafranov (GS) equation, for a quantity called the flux function, $`\mathrm{\Psi }`$. The GS equation was derived in the Newtonian case l86 , the Schwarzschild spacetime case ml86 , and the Kerr spacetime case n91 ; b97 . The flux function $`\mathrm{\Psi }`$ is such that it is constant over each surface generated by rotating the magnetic field lines (or equivalently the flow lines) about the axis of symmetry and the GS equation determines the transfield equilibrium. Any physical quantities can be calculated from the solution $`\mathrm{\Psi }`$ of the GS equation. However, the GS equation in noncircular spacetimes has never been derived explicitly.

In this paper we derive the GS equation explicitly in noncircular (i.e., the most general) stationary axisymmetric spacetimes. This is a first step toward the study of equilibrium configurations of magnetars. An attempt to solve the GS equation will be discussed in a subsequent paper. This paper is organized as follows. In Sec. II, we briefly review the conservation laws in stationary axisymmetric general relativistic ideal MHD systems bo78 ; bo79  that are used to characterize the matter and electromagnetic field configurations. We neglect dissipative effects, which is a reasonable assumption because of the high conductivity and the low viscosity in neutron stars. In Sec. III, we derive the GS equation in an un-elucidated form. At this stage it is not clear if the GS equation is a second-order differential equation for the flux function $`\mathrm{\Psi }`$. In Sec. IV, we briefly review the $`(2+1)+1`$ formalism by Gourgoulhon and Bonazzola gb93  to describe the geometry of noncircular stationary axisymmetric spacetimes in a transparent way. We do not, however, discuss the Einstein equations but assume the geometry to be given. In Sec. V, we explicitly demonstrate that all physical quantities except for the metric can be evaluated from the flux function $`\mathrm{\Psi }`$. In Sec. VI, we write down the GS equation in the covariant form projected onto the 2-surface $`\mathrm{\Sigma }_{t\phi }`$ spanned by $`t=`$const. and $`\phi =`$const.. We also discuss various limits of the GS equation. Finally, we summarize our result in Sec. VII. The energy-momentum tensor decomposed in the $`(2+1)+1`$ form is given in Appendix A, and notation and symbols are summarized in Appendix B.

We use the units $`c=G=k_B=1`$. Greek indices ($`\mu ,\nu ,\alpha ,\beta ,\mathrm{}`$) run from $`0`$ to $`3`$, small Latin indices ($`i,j,k,\mathrm{}`$) from $`1`$ to $`3`$, and capital Latin indices ($`A,B,C,\mathrm{}`$) from $`1`$ to $`2`$, where $`x^0=t`$ and $`x^3=\phi `$. The signature of the 4-metric is $`(,+,+,+)`$.

## II Basic equations and conservation laws

### II.1 Basic equations for general relativistic magnetohydrodynamics

The basic equations governing a general relativistic ideal MHD system are as follows l67 ; nt73 . Baryons are conserved,

$`(\rho u^\mu )_{;\mu }=u^\mu \rho _{,\mu }+\rho u_{}^{\mu }{}_{;\mu }{}^{}=0,`$ (1)
where $`\rho `$ is the rest mass density (i.e., the baryon mass times the baryon number density) and $`u^\mu `$ is the fluid 4-velocity with

$`u_\mu u^\mu =1.`$ (2)
The electromagnetic field is governed by the Maxwell equations,

$`F_{[\mu \nu ;\alpha ]}`$ $`=`$ $`0,`$ (3)
$`F_{}^{\mu \nu }{}_{;\nu }{}^{}`$ $`=`$ $`4\pi J^\mu ,`$ (4)
where $`F_{\mu \nu }`$ and $`J^\mu `$ are the field strength tensor and the electric current 4-vector, respectively. Equation (3) implies the existence of the vector potential $`A_\mu `$,

$`F_{\mu \nu }=A_{\nu ,\mu }A_{\mu ,\nu }.`$ (5)
The electric and magnetic fields in the fluid rest frame are defined as

$`E_\mu `$ $`=`$ $`F_{\mu \nu }u^\nu ,`$ (6)
$`B^\mu `$ $`=`$ $`{\displaystyle \frac{1}{2}}ϵ^{\mu \nu \alpha \beta }u_\nu F_{\alpha \beta },`$ (7)
where $`ϵ_{\mu \nu \alpha \beta }`$ is the Levi-Civita antisymmetric tensor with $`ϵ_{0123}=\sqrt{g}`$. Equations (6) and (7) are inverted to give

$`F_{\mu \nu }=u_\mu E_\nu u_\nu E_\mu +ϵ_{\mu \nu \alpha \beta }u^\alpha B^\beta ,`$ (8)
with $`E_\mu u^\mu =B_\mu u^\mu =0`$. In the ideal MHD, we assume the perfect conductivity, so that

$`E_\mu `$ $`=`$ $`F_{\mu \nu }u^\nu =0.`$ (9)
The equations of motion for the fluid are given by $`T_{}^{\mu \nu }{}_{;\nu }{}^{}=0`$, where $`T^{\mu \nu }`$ is the total energy-momentum tensor of the fluid and electromagnetic fields,

$`T^{\mu \nu }=(\rho +\rho ϵ+p)u^\mu u^\nu +pg^{\mu \nu }+{\displaystyle \frac{1}{4\pi }}\left(F^{\mu \alpha }F_{}^{\nu }{}_{\alpha }{}^{}{\displaystyle \frac{1}{4}}g^{\mu \nu }F^{\alpha \beta }F_{\alpha \beta }\right).`$ (10)
Then we obtain the conservation of the fluid energy,

$`u^\mu (\rho +\rho ϵ)_{,\mu }+(\rho +\rho ϵ+p)u_{}^{\mu }{}_{;\mu }{}^{}=0,`$ (11)
and the Euler equations,

$`(\rho +\rho ϵ+p)u_{}^{\mu }{}_{;\nu }{}^{}u^\nu +(g^{\mu \nu }+u^\mu u^\nu )p_{,\nu }F^{\mu \nu }J_\nu =0,`$ (12)
where $`ϵ`$ and $`p`$ are the internal energy per unit mass and pressure, respectively. Eliminating $`u_{}^{\mu }{}_{;\mu }{}^{}`$ from Eqs. (1) and (11) gives

$`u^\mu (\rho +\rho ϵ)_{,\mu }=\mu u^\mu \rho _{,\mu },`$ (13)
where

$`\mu =1+ϵ+{\displaystyle \frac{p}{\rho }},`$ (14)
is the enthalpy per unit mass. Assuming local thermodynamic equilibrium, the first law of thermodynamics is given by

$`dϵ=pd\left({\displaystyle \frac{1}{\rho }}\right)+TdS,`$ (15)
where $`S`$ and $`T`$ are the entropy per unit mass and the temperature. Then Eqs. (14) and (15) imply

$`d\mu ={\displaystyle \frac{dp}{\rho }}+TdS.`$ (16)
Finally we supply the equation of state,

$`p=p(\rho ,S).`$ (17)
### II.2 Conservation laws in a stationary axisymmetric spacetime

Here we recapitulate the conservation laws derived from the basic equations in the previous subsection in a stationary axisymmetric spacetime. There exists two Killing vectors associated with stationarity and axisymmetry, which we denote by $`\eta ^\mu `$ and $`\xi ^\mu `$, respectively. The Lie derivatives of all physical quantities along the Killing vectors must vanish, e.g., $`_\xi u^\mu =\xi ^\nu u_{}^{\mu }{}_{;\nu }{}^{}u^\nu \xi _{}^{\mu }{}_{;\nu }{}^{}=0`$. We take $`\eta ^\mu =(/t)^\mu `$ and $`\xi ^\mu =(/\phi )^\mu `$ so that $`x^0=t`$ and $`x^3=\phi `$ are the time and azimuthal coordinates associated with the Killing vectors $`\eta ^\mu `$ and $`\xi ^\mu `$, respectively. Thus all physical quantities are independent of $`t`$ and $`\phi `$.

Bekenstein and Oron bo78 ; bo79  showed that a stationary axisymmetric system has several conserved quantities along each flow line. This is a general relativistic generalization of Ferraro’s integrability condition fe37 ; fe54 ; wo59 . By exploiting the gauge freedom to make $`A_{\mu ,\nu }\eta ^\nu =A_{\mu ,0}=0`$ and $`A_{\mu ,\nu }\xi ^\nu =A_{\mu ,3}=0`$, we can show that the magnetic potential $`\mathrm{\Psi }:=A_\mu \xi ^\mu =A_3`$ as well as the electric potential $`\mathrm{\Phi }:=A_\mu \eta ^\mu =A_0`$ are constant along each flow line, i.e., $`u^\mu (\xi ^\nu A_\nu )_{,\mu }=u^\mu (\eta ^\nu A_\nu )_{,\mu }=0`$ be85 . Henceforth we label the flow line by $`\mathrm{\Psi }`$, which we refer to as the flux function as in the non-relativistic case l86 . The $`\mathrm{\Psi }=\mathrm{const}.`$ surfaces are called the flux surfaces, which are generated by rotating the magnetic field lines (or the flow lines) about the axis of symmetry.

According to Bekenstein and Oron bo78 ; bo79 , one can show that

$`F_{03}`$ $`=`$ $`0,`$ (18)
$`F_{0A}`$ $`=`$ $`\mathrm{\Omega }F_{A3},`$ (19)
$`F_{31}`$ $`=`$ $`\mathrm{\Psi }_{,1}=C\sqrt{g}\rho u^2,`$ (20)
$`F_{23}`$ $`=`$ $`\mathrm{\Psi }_{,2}=C\sqrt{g}\rho u^1,`$ (21)
$`F_{12}`$ $`=`$ $`C\sqrt{g}\rho (u^3\mathrm{\Omega }u^0),`$ (22)
where $`\mathrm{\Omega }(\mathrm{\Psi })`$ and $`C(\mathrm{\Psi })`$ are conserved along each flow line and hence are functions of the flux function $`\mathrm{\Psi }`$. The above equations are effectively first integrals of the Maxwell equations (3). It may be useful to rewrite the above equations as

$`B^\mu =C\rho \left[(u_0+\mathrm{\Omega }u_3)u^\mu +\eta ^\mu +\mathrm{\Omega }\xi ^\mu \right].`$ (23)
Note that $`\mathrm{\Omega }(\mathrm{\Psi })`$ is the $`\mathrm{\Psi }`$-derivative of the electric potential, $`\mathrm{\Omega }(\mathrm{\Psi })=d\mathrm{\Phi }/d\mathrm{\Psi }`$, and coincides with the angular velocity $`d\phi /dt=u^3/u^0=\mathrm{\Omega }`$ if there is no toroidal field $`F_{12}=0`$ and $`C0`$, from Eq. (22). In addition, one can show that $`E(\mathrm{\Psi })`$, $`L(\mathrm{\Psi })`$ and $`D(\mathrm{\Psi })`$ are also conserved along each flow line bo78 , where

$`D`$ $`=`$ $`\mu (u_0+\mathrm{\Omega }u_3),`$ (24)
$`E`$ $`=`$ $`\left(\mu +{\displaystyle \frac{B^2}{4\pi \rho }}\right)u_0+C(u_0+\mathrm{\Omega }u_3){\displaystyle \frac{B_0}{4\pi }},`$ (25)
$`L`$ $`=`$ $`\left(\mu +{\displaystyle \frac{B^2}{4\pi \rho }}\right)u_3+C(u_0+\mathrm{\Omega }u_3){\displaystyle \frac{B_3}{4\pi }},`$ (26)
and $`B^2:=B_\mu B^\mu `$. Equations (24), (25) and (26) are first integrals of the equations of motion $`T^{\mu \nu }{}_{;\nu }{}^{}=0`$ along $`B_\mu `$, $`\eta _\mu `$ and $`\xi _\mu `$, respectively. In particular, Eq. (25) is a generalized Bernouilli’s equation. These conserved quantities are not mutually independent but there is a relation among them bo78 ,

$`D=E\mathrm{\Omega }L.`$ (27)
Together with Eqs. (24) – (26), this implies that

$`{\displaystyle \frac{B^2}{\rho }}+C(B_0+\mathrm{\Omega }B_3)=0.`$ (28)
With Eq. (28), we can rewrite $`E(\mathrm{\Psi })`$ and $`L(\mathrm{\Psi })`$ in Eqs. (25) and (26) as

$`E`$ $`=`$ $`\mu u_0{\displaystyle \frac{1}{4\pi }}C\mathrm{\Omega }(u_0B_3u_3B_0),`$ (29)
$`L`$ $`=`$ $`\mu u_3+{\displaystyle \frac{1}{4\pi }}C(u_0B_3u_3B_0).`$ (30)
Finally, from Eqs. (13) – (15) one finds that the entropy per unit mass $`S`$ is conserved along each flow line, $`u^\mu S_{,\mu }=0`$, as a result of the perfect fluid form of the energy-momentum tensor. For a stationary axisymmetric spacetime, this implies that $`S`$ is a function of $`\mathrm{\Psi }`$, $`S=S(\mathrm{\Psi })`$.

In summary, for a give flux function $`\mathrm{\Psi }`$, there exists five conserved quantities, $`E(\mathrm{\Psi })`$, $`L(\mathrm{\Psi })`$, $`\mathrm{\Omega }(\mathrm{\Psi })`$, $`C(\mathrm{\Psi })`$ and $`S(\mathrm{\Psi })`$. Except for $`S(\mathrm{\Psi })`$, there are no perfectly relevant physical interpretations of these quantities. Nevertheless, by considering several limiting cases, we may associate them with terms that describe their qualitative nature. We may call $`E(\mathrm{\Psi })`$ the total energy per unit mass, $`L(\mathrm{\Psi })`$ the total angular momentum per unit mass, $`\mathrm{\Omega }(\mathrm{\Psi })`$ the angular velocity, and $`C(\mathrm{\Psi })`$ the magnetic field strength relative to the magnitude of meridional flow. Since these conserved quantities are essentially the first integrals of the equations of motion, specification of these functions characterizes the configuration of the fluid flow and the electromagnetic field. As we will see in Sec. V, all physical quantities are completely determined once these conserved quantities are given as functions of the flux function $`\mathrm{\Psi }`$, provided that the spatial configuration of $`\mathrm{\Psi }`$ is known. Therefore the problem reduces to solving an equation for the flux function that determines the spatial configuration of $`\mathrm{\Psi }`$, that is, the GS equation.

## III Grad-Shafranov equation in the component expression

The GS equation is given by the transfield component of the Euler equations (12). In this section, we consider the $`x^A`$-derivative of the flux function in the Euler equations, and factorize the resulting equation to derive the GS equation. We express equations in terms of their explicit coordinate components, since it is the most straightforward way to incorporate the symmetry, e.g., $`(\mathrm{})_{,0}=(\mathrm{})_{,3}=0`$. Accordingly, the GS equation is given in the component expression. A covariant form of the GS equation based on $`(2+1)+1`$ formalism, which may be useful for numerical calculations, will be given in Sec. VI.

Using Eq. (16), the Euler equations (12) can be written as

$`\rho \mu u_{\mu ;\nu }u^\nu +\rho \mu _{,\mu }+\rho u_\mu u^\nu \mu _{,\nu }F_{\mu \nu }J^\nu \rho TS_{,\mu }=0,`$ (31)
where we have used that $`u^\mu S_{,\mu }=0`$. First, let us consider the $`x^1`$-component of the Euler equations (31). The first term can be expressed as

$`\rho \mu u_{1;\mu }u^\mu =\rho \mu (u^1u_{1;1}+u^2u_{1;2})\rho \mu \left(u^0\mathrm{\Gamma }_{10}^\mu u_\mu +u^3\mathrm{\Gamma }_{13}^\mu u_\mu \right),`$ (32)
where $`\mathrm{\Gamma }_{\mu \nu }^\alpha `$ is the Christoffel symbol. With Eq. (2), the third term in Eq. (31) can be transformed as

$`\rho u_1u^\mu \mu _{,\mu }`$ $`=`$ $`\rho u_1u^1\mu _{,1}+\rho u_1u^2\mu _{,2}`$ (33)
$`=`$ $`\rho \mu _{,1}(1+u_0u^0+u_2u^2+u_3u^3)+\rho u_1u^2\mu _{,2}`$
$`=`$ $`\rho u^2(u_1\mu _{,2}u_2\mu _{,1})\rho \mu _{,1}\rho u^0(\mu u_0)_{,1}\rho u^3(\mu u_3)_{,1}`$
$`+\rho \mu u^0(u_{0;1}+\mathrm{\Gamma }_{01}^\mu u_\mu )+\rho \mu u^3(u_{3;1}+\mathrm{\Gamma }_{31}^\mu u_\mu )`$
$`=`$ $`\rho u^2[(\mu u_1)_{;2}(\mu u_2)_{;1})]\rho \mu _{,1}\rho u^0(\mu u_0)_{,1}\rho u^3(\mu u_3)_{,1}`$
$`\rho \mu (u^1u_{1;1}+u^2u_{1;2})+\rho \mu (u^0\mathrm{\Gamma }_{01}^\mu u_\mu +u^3\mathrm{\Gamma }_{13}^\mu u_\mu ),`$
where we have used the fact $`u^\mu u_{\mu ;\nu }=0`$. Thus putting Eqs. (32) and (33) together, and using Eq. (20), we find

$`\rho (\mu u_{1;\mu }u^\mu +\mu _{,1}+u_1u^\mu \mu _{,\mu })`$ $`=`$ $`\rho u^2\left[(\mu u_1)_{;2}(\mu u_2)_{;1}\right]\rho u^0(\mu u_0)_{,1}\rho u^3(\mu u_3)_{,1},`$ (34)
$`=`$ $`{\displaystyle \frac{1}{C\sqrt{g}}}\left[(\mu u_1)_{;2}(\mu u_2)_{;1}\right]\mathrm{\Psi }_{,1}\rho u^0(\mu u_0)_{,1}\rho u^3(\mu u_3)_{,1}.`$
Next, the fourth term in Eq. (31) can be transformed as

$`F_{1\mu }J^\mu `$ $`=`$ $`F_{10}J^0F_{13}J^3{\displaystyle \frac{1}{4\pi \sqrt{g}}}F_{1A}(\sqrt{g}F^{AB})_{,B}`$ (35)
$`=`$ $`\left(\mathrm{\Omega }J^0J^3\right)\mathrm{\Psi }_{,1}+{\displaystyle \frac{1}{4\pi \sqrt{g}}}F_{12}(\sqrt{g}F^{12})_{,1},`$
where the second line follows from Eqs. (19) and (20). From Eqs. (8) and (9), we have $`F^{12}=(u_0B_3u_3B_0)/\sqrt{g}`$. Therefore, together with Eq. (22) that gives $`F_{12}`$, the last term in Eq. (35) is expressed as

$`{\displaystyle \frac{1}{4\pi \sqrt{g}}}F_{12}(\sqrt{g}F^{12})_{,1}`$ $`=`$ $`{\displaystyle \frac{1}{4\pi }}C\rho (u^3\mathrm{\Omega }u^0)(u_0B_3u_3B_0)_{,1}`$ (36)
$`=`$ $`{\displaystyle \frac{1}{4\pi }}\rho u^0(u_0B_3u_3B_0)(C\mathrm{\Omega })^{}\mathrm{\Psi }_{,1}+{\displaystyle \frac{1}{4\pi }}\rho u^3(u_0B_3u_3B_0)C^{}\mathrm{\Psi }_{,1},`$
$``$ $`\rho u^0\left[{\displaystyle \frac{1}{4\pi }}C\mathrm{\Omega }(u_0B_3u_3B_0)\right]_{,1}\rho u^3\left[{\displaystyle \frac{1}{4\pi }}C(u_0B_3u_3B_0)\right]_{,1},`$
where primes denote differentiation with respect to $`\mathrm{\Psi }`$. Finally, the last term in Eq. (31) gives

$`\rho TS_{,1}=\rho TS^{}\mathrm{\Psi }_{,1}.`$ (37)
Combining Eqs. (34) – (37), we have

$`\{{\displaystyle \frac{1}{C\sqrt{g}}}[(\mu u_1)_{,2}(\mu u_2)_{,1}](J^3\mathrm{\Omega }J^0)`$
$`{\displaystyle \frac{1}{4\pi }}\rho u^0(u_0B_3u_3B_0)\left(C\mathrm{\Omega }\right)^{}+{\displaystyle \frac{1}{4\pi }}\rho u^3(u_0B_3u_3B_0)C^{}\rho TS^{}\}\mathrm{\Psi }_{,1}`$
$`\rho u^0\left[\mu u_0{\displaystyle \frac{1}{4\pi }}C\mathrm{\Omega }\left(u_0B_3u_3B_0\right)\right]_{,1}\rho u^3\left[\mu u_3+{\displaystyle \frac{1}{4\pi }}C\left(u_0B_3u_3B_0\right)\right]_{,1}=0.`$ (38)
Recalling the expressions for $`E`$ and $`L`$ given by Eqs. (29) and (30), respectively, we see that the last two terms are just their derivatives. Therefore, we can factor out the $`x^1`$-derivative of the flux function in Eq. (38) to obtain

$`\{{\displaystyle \frac{1}{C\sqrt{g}}}[(\mu u_1)_{,2}(\mu u_2)_{,1}](J^3\mathrm{\Omega }J^0)`$
$`+\rho u^0[E^{}{\displaystyle \frac{1}{4\pi }}(u_0B_3u_3B_0)\left(C\mathrm{\Omega }\right)^{}]\rho u^3[L^{}{\displaystyle \frac{1}{4\pi }}(u_0B_3u_3B_0)C^{}]\rho TS^{}\}\mathrm{\Psi }_{,1}=0.`$ (39)
The same analysis applies to the $`x^2`$-component of Eq. (31), and one finds the above equation (39) with the replacement of $`\mathrm{\Psi }_{,1}`$ by $`\mathrm{\Psi }_{,2}`$. Therefore, by assuming $`\mathrm{\Psi }_{,A}0`$ ($`A=1,2`$), the GS equation is given by

$`J^3\mathrm{\Omega }J^0+{\displaystyle \frac{1}{C\sqrt{g}}}\left[(\mu u_1)_{,2}(\mu u_2)_{,1}\right]\rho u^0\left[E^{}\mathrm{\Lambda }\left(C\mathrm{\Omega }\right)^{}\right]+\rho u^3\left[L^{}\mathrm{\Lambda }C^{}\right]+\rho TS^{}=0,`$ (40)
where, for convenience, we have introduced an auxiliary quantity $`\mathrm{\Lambda }`$ defined by

$`\mathrm{\Lambda }={\displaystyle \frac{1}{4\pi }}(u_0B_3u_3B_0).`$ (41)
At this stage, however, it is not clear if Eq. (40) gives a second-order, nonlinear partial differential equation for the flux function $`\mathrm{\Psi }`$, since the dependence on the flux function is unknown. In Sec. V, we explicitly demonstrate that all the physical quantities appearing in the above equation can be expressed in terms of $`\mathrm{\Psi }`$ and its $`x^A`$-derivatives, and in Sec. VI we derive the GS equation in the covariant form and we make it explicit that it is indeed a second-order, non-linear differential equation for $`\mathrm{\Psi }`$.

## IV $`(2+1)+1`$ DECOMPOSITION

In this section, we briefly review the $`(2+1)+1`$ formalism of the Einstein equations for stationary axisymmetric spacetimes developed by Gourgoulhon and Bonazzola gb93 , in order to describe our metric in a covariant fashion. Note that this formalism is different from the $`(2+1)+1`$ formalism by Maeda, Sasaki, Nakamura, and Miyama m80 ; nok87  and Sasaki s84 , which is suitable to the axisymmetric gravitational collapse. Here we adopt the formalism by Gourgoulhon and Bonazzola because it is more convenient for a spacetime which is not only axisymmetric but also stationary.

Let $`n^\mu `$ be the unit timelike 4-vector orthogonal to the $`t=\mathrm{const}.`$ hypersurface $`\mathrm{\Sigma }_t`$ and oriented in the direction of increasing $`t`$,

$`n_\mu =Nt_{,\mu }.`$ (42)
The lapse function $`N`$ is determined by the requirement

$`n_\mu n^\mu =1.`$ (43)
The 3-metric induced by $`g_{\mu \nu }`$ on $`\mathrm{\Sigma }_t`$ is given by

$`h_{\mu \nu }=g_{\mu \nu }+n_\mu n_\nu .`$ (44)
Similarly, let $`m_\mu `$ be the unit spacelike 4-vector orthogonal to the $`t=\mathrm{const}.`$ and $`\phi =\mathrm{const}.`$ 2-surface $`\mathrm{\Sigma }_{t\phi }`$ and oriented in the direction of increasing $`\phi `$,

$`m_\mu =Mh_{\mu }^{}{}_{}{}^{\nu }\phi _{,\nu }=M\phi _{|\mu },`$ (45)
where the vertical stroke $`|`$ denotes the covariant derivative associated with the 3-metric $`h_{\mu \nu }`$. The coefficient $`M`$ is determined by

$`m_\mu m^\mu =1.`$ (46)
The induced 2-metric on $`\mathrm{\Sigma }_{t\phi }`$ is given by

$`H_{\mu \nu }=h_{\mu \nu }m_\mu m_\nu =g_{\mu \nu }+n_\mu n_\nu m_\mu m_\nu .`$ (47)
The covariant derivative associated with the 2-metric $`H_{\mu \nu }`$ is denoted by a double vertical stroke $``$. There is a relation between the determinants as

$`\sqrt{g}=N\sqrt{h}=NM\sqrt{H}.`$ (48)
Any 4-vector can be decomposed into its projection onto $`\mathrm{\Sigma }_{t\phi }`$, the component parallel to $`n_\mu `$ and that to $`m_\mu `$. The Killing vectors are decomposed as

$`\eta ^\mu `$ $`=`$ $`Nn^\mu N^\mu =Nn^\mu MN^\phi m^\mu N_{\mathrm{\Sigma }}^{}{}_{}{}^{\mu },`$ (49)
$`\xi ^\mu `$ $`=`$ $`Mm^\mu M_{\mathrm{\Sigma }}^{}{}_{}{}^{\mu },`$ (50)
where the shift vector $`N^\mu `$ is (minus) the projection of $`\eta ^\mu `$ onto $`\mathrm{\Sigma }_t`$, $`M_\mathrm{\Sigma }^\mu `$ is (minus) the projection of $`\xi ^\mu `$ onto $`\mathrm{\Sigma }_{t\phi }`$, and $`N_\mathrm{\Sigma }^\mu `$ is the projection of $`N^\mu `$ onto $`\mathrm{\Sigma }_{t\phi }`$. For our choice of coordinates, i.e., for $`x^0=t`$ and $`x^3=\phi `$, the component expressions for $`n^\mu `$, $`n_\mu `$, $`m^\mu `$ and $`m_\mu `$ are

$`n_\mu `$ $`=`$ $`(N,0,0,0),`$ (51)
$`n^\mu `$ $`=`$ $`({\displaystyle \frac{1}{N}},{\displaystyle \frac{N^1}{N}},{\displaystyle \frac{N^2}{N}},{\displaystyle \frac{N^\phi }{N}}),`$ (52)
$`m_\mu `$ $`=`$ $`(MN^\phi ,0,0,M),`$ (53)
$`m^\mu `$ $`=`$ $`(0,{\displaystyle \frac{M_{\mathrm{\Sigma }}^{}{}_{}{}^{1}}{M}},{\displaystyle \frac{M_{\mathrm{\Sigma }}^{}{}_{}{}^{2}}{M}},{\displaystyle \frac{1}{M}}).`$ (54)
Note that $`N_\mathrm{\Sigma }{}_{}{}^{\mu }=(0,N_\mathrm{\Sigma }{}_{}{}^{1},N_\mathrm{\Sigma }{}_{}{}^{2},0)`$ and $`N^A=N_\mathrm{\Sigma }{}_{}{}^{A}+N^\phi M_\mathrm{\Sigma }^A`$.

The explicit component expressions of $`g^{\mu \nu }`$, $`g_{\mu \nu }`$, $`h^{\mu \nu }`$, and $`h_{\mu \nu }`$ are given by

$`\left(\begin{array}{cc}g_{00}& g_{0j}\\ g_{i0}& g_{ij}\end{array}\right)`$ $`=`$ $`\left(\begin{array}{cc}N_kN^kN^2& N_j\\ & \\ N_i& h_{ij}\end{array}\right),`$ (60)
$`\left(\begin{array}{cc}g^{00}& g^{0j}\\ g^{i0}& g^{ij}\end{array}\right)`$ $`=`$ $`\left(\begin{array}{cc}{\displaystyle \frac{1}{N^2}}& {\displaystyle \frac{N^j}{N^2}}\\ & \\ {\displaystyle \frac{N^i}{N^2}}& h^{ij}{\displaystyle \frac{N^iN^j}{N^2}}\end{array}\right),`$ (66)
$`\left(\begin{array}{cc}h_{AB}& h_{A3}\\ h_{3B}& h_{33}\end{array}\right)`$ $`=`$ $`\left(\begin{array}{cc}H_{AB}& M_{\mathrm{\Sigma }}^{}{}_{A}{}^{}\\ & \\ M_{\mathrm{\Sigma }}^{}{}_{B}{}^{}& M^2+M_{\mathrm{\Sigma }}^{}{}_{A}{}^{}M_{\mathrm{\Sigma }}^{}{}_{}{}^{A}\end{array}\right),`$ (72)
$`\left(\begin{array}{cc}h^{AB}& h^{A3}\\ h^{3B}& h^{33}\end{array}\right)`$ $`=`$ $`\left(\begin{array}{cc}H^{AB}+{\displaystyle \frac{M_{\mathrm{\Sigma }}^{}{}_{}{}^{A}M_{\mathrm{\Sigma }}^{}{}_{}{}^{B}}{M^2}}& {\displaystyle \frac{M_{\mathrm{\Sigma }}^{}{}_{}{}^{A}}{M^2}}\\ & \\ {\displaystyle \frac{M_{\mathrm{\Sigma }}^{}{}_{}{}^{B}}{M^2}}& {\displaystyle \frac{1}{M^2}}\end{array}\right).`$ (78)
where $`i,j,k=1,2,3`$ and $`A,B=1,2`$. We can express the 4-metric $`g_{\mu \nu }`$ in terms of $`N`$, $`N^\phi `$, $`N_\mathrm{\Sigma }^A`$, $`M`$, $`M_\mathrm{\Sigma }^A`$ and $`H_{AB}`$ as

$`g_{\mu \nu }dx^\mu dx^\nu `$ $`=`$ $`[N^2M^2\left(N^\phi \right)^2N_\mathrm{\Sigma }{}_{A}{}^{}N_{\mathrm{\Sigma }}^{}{}_{}{}^{A}]dt^22(M^2N^\phi N_\mathrm{\Sigma }{}_{}{}^{A}M_{\mathrm{\Sigma }}^{}{}_{A}{}^{})dtd\phi `$ (79)
$``$ $`2N_\mathrm{\Sigma }{}_{A}{}^{}dtdx^A2M_\mathrm{\Sigma }{}_{A}{}^{}d\phi dx^A+H_{AB}dx^Adx^B+(M^2+M_\mathrm{\Sigma }{}_{A}{}^{}M_{\mathrm{\Sigma }}^{}{}_{}{}^{A})d\phi ^2,`$
where the functions $`N`$, $`N^\phi `$, $`N_\mathrm{\Sigma }^A`$, $`M`$, $`M_\mathrm{\Sigma }^A`$ and $`H_{AB}`$ depend only on the coordinate $`(x^1,x^2)`$. Since we only assume that physical quantities are independent of $`x^0=t`$ and $`x^3=\phi `$, the metric $`g_{\mu \nu }`$ in Eq. (79) has some freedom in the choice of coordinates. We will leave the coordinate freedom unspecified. In Sec. VI the covariant GS equation will be given as an equation projected onto $`\mathrm{\Sigma }_{t\phi }`$.

## V Physical quantities from flux function $`\mathrm{\Psi }`$

Provided that the metric $`g_{\mu \nu }`$ is given and the conserved quantities $`E(\mathrm{\Psi })`$, $`L(\mathrm{\Psi })`$, $`\mathrm{\Omega }(\mathrm{\Psi })`$, $`C(\mathrm{\Psi })`$ and $`S(\mathrm{\Psi })`$ are given as functions of $`\mathrm{\Psi }`$, all the physical quantities can be evaluated once the (effectively 2-dimensional) configuration of the flux function $`\mathrm{\Psi }`$ is known (see c86b  for the circular case). In this section, we explicitly demonstrate this fact for the most general case of noncircular spacetimes.

### V.1 Fluid 4-velocity $`u^\mu `$

First let us consider the fluid 4-velocity $`u^\mu `$. It is useful to prepare two vectors $`\eta ^\mu +\mathrm{\Omega }\xi ^\mu `$ and $`\xi ^\mu +\mathrm{\Theta }\eta ^\mu `$ constructed from two Killing vectors $`\eta ^\mu `$ and $`\xi ^\mu `$, and to make them orthogonal to each other $`(\eta ^\mu +\mathrm{\Omega }\xi ^\mu )(\xi _\mu +\mathrm{\Theta }\eta _\mu )=0`$ by taking

$`\mathrm{\Theta }={\displaystyle \frac{\xi _\nu (\eta ^\nu +\mathrm{\Omega }\xi ^\nu )}{\eta _\mu (\eta ^\mu +\mathrm{\Omega }\xi ^\mu )}}={\displaystyle \frac{g_{03}+\mathrm{\Omega }g_{33}}{g_{00}+\mathrm{\Omega }g_{03}}}.`$ (80)
Then we can decompose the fluid 4-velocity in the coordinate bases as

$`u^\mu =u_\eta (\eta ^\mu +\mathrm{\Omega }\xi ^\mu )+u_\xi (\xi ^\mu +\mathrm{\Theta }\eta ^\mu )+\stackrel{~}{u}_\mathrm{\Sigma }{}_{}{}^{\mu },`$ (81)
where $`\eta ^\mu +\mathrm{\Omega }\xi ^\mu =(1,0,0,\mathrm{\Omega })`$, $`\xi ^\mu +\mathrm{\Theta }\eta ^\mu =(\mathrm{\Theta },0,0,1)`$ and $`\stackrel{~}{u}_\mathrm{\Sigma }{}_{}{}^{\mu }=(0,u^1,u^2,0)`$ in the component expressions, and hence $`\stackrel{~}{u}_\mathrm{\Sigma }{}_{}{}^{\mu }n_{\mu }^{}=\stackrel{~}{u}_\mathrm{\Sigma }{}_{}{}^{\mu }m_{\mu }^{}=0`$ from Eqs. (51) and (53), and

$`u^0`$ $`=`$ $`u_\eta +\mathrm{\Theta }u_\xi ,`$ (82)
$`u^3`$ $`=`$ $`u_\xi +\mathrm{\Omega }u_\eta .`$ (83)
The decomposition in Eq. (81) is not conforming to the spirit of the $`(2+1)+1`$ formalism but makes it easy to obtain the coefficients $`u_\eta `$ and $`u_\xi `$ as shown below.

From Eqs. (20) and (21), the term $`\stackrel{~}{u}_\mathrm{\Sigma }^\mu `$ is given by

$`\stackrel{~}{u}_\mathrm{\Sigma }{}_{}{}^{\mu }={\displaystyle \frac{1}{NMC\rho }}ϵ^{\mu \nu }\mathrm{\Psi }_{,\nu },`$ (84)
where the antisymmetric tensor $`ϵ^{\mu \nu }`$ is defined by

$`ϵ^{\mu \nu }=ϵ^{\mu \nu \alpha \beta }n_\alpha m_\beta .`$ (85)
With Eqs. (23) – (26), (49) and (50), the coefficients $`u_\eta `$ and $`u_\xi `$ are expressed as

$`u_\eta `$ $`=`$ $`{\displaystyle \frac{E\mathrm{\Omega }L}{G_\eta \mu }}{\displaystyle \frac{\stackrel{~}{N}_\mathrm{\Sigma }}{G_\eta }}={\displaystyle \frac{D}{G_\eta \mu }}{\displaystyle \frac{\stackrel{~}{N}_\mathrm{\Sigma }}{G_\eta }},`$ (86)
$`u_\xi `$ $`=`$ $`{\displaystyle \frac{(L\mathrm{\Theta }E)}{G_\xi \mu }}\left({\displaystyle \frac{4\pi \mu }{G_\eta C^2\rho }}\right)\left(1{\displaystyle \frac{4\pi \mu }{G_\eta C^2\rho }}\right)^1+{\displaystyle \frac{\stackrel{~}{M}_\mathrm{\Sigma }}{G_\xi }},`$ (87)
where $`G_\eta `$ and $`G_\xi `$ are defined by

$`G_\eta `$ $`=`$ $`(\eta _\mu +\mathrm{\Omega }\xi _\mu )(\eta ^\mu +\mathrm{\Omega }\xi ^\mu )=(g_{00}+2\mathrm{\Omega }g_{03}+\mathrm{\Omega }^2g_{33}),`$ (88)
$`G_\xi `$ $`=`$ $`(\xi _\mu +\mathrm{\Theta }\eta _\mu )(\xi ^\mu +\mathrm{\Theta }\eta ^\mu )=g_{33}+2\mathrm{\Theta }g_{03}+\mathrm{\Theta }^2g_{00},`$ (89)
and $`\stackrel{~}{N}_\mathrm{\Sigma }`$ and $`\stackrel{~}{M}_\mathrm{\Sigma }`$ are defined by

$`\stackrel{~}{N}_\mathrm{\Sigma }`$ $`=`$ $`\stackrel{~}{u}_\mathrm{\Sigma }{}_{}{}^{\mu }(N_\mathrm{\Sigma }{}_{\mu }{}^{}+\mathrm{\Omega }M_\mathrm{\Sigma }{}_{\mu }{}^{}),`$ (90)
$`\stackrel{~}{M}_\mathrm{\Sigma }`$ $`=`$ $`\stackrel{~}{u}_\mathrm{\Sigma }{}_{}{}^{\mu }(M_\mathrm{\Sigma }{}_{\mu }{}^{}+\mathrm{\Theta }N_\mathrm{\Sigma }{}_{\mu }{}^{}).`$ (91)
Note that, using Eq. (84), $`\stackrel{~}{N}_\mathrm{\Sigma }`$ and $`\stackrel{~}{M}_\mathrm{\Sigma }`$ are expressed as

$`\stackrel{~}{N}_\mathrm{\Sigma }={\displaystyle \frac{1}{NMC\rho }}ϵ^{AB}\mathrm{\Psi }_B(N_\mathrm{\Sigma }{}_{A}{}^{}+\mathrm{\Omega }M_\mathrm{\Sigma }{}_{A}{}^{}),`$ (92)
$`\stackrel{~}{M}_\mathrm{\Sigma }={\displaystyle \frac{1}{NMC\rho }}ϵ^{AB}\mathrm{\Psi }_B(M_\mathrm{\Sigma }{}_{A}{}^{}+\mathrm{\Theta }N_\mathrm{\Sigma }{}_{A}{}^{}).`$ (93)
Note also that

$`M_{\mathrm{Alf}}^{}{}_{}{}^{2}`$ $`:=`$ $`{\displaystyle \frac{4\pi \mu }{G_\eta C^2\rho }}`$ (94)
$`=`$ $`{\displaystyle \frac{4\pi \mu \rho }{B^2}}(G_\xi u_\xi {}_{}{}^{2}+\stackrel{~}{u}_\mathrm{\Sigma }{}_{}{}^{A}\stackrel{~}{u}_{\mathrm{\Sigma }}^{}{}_{A}{}^{}2u_\xi \stackrel{~}{M}_\mathrm{\Sigma }+{\displaystyle \frac{\stackrel{~}{N}_\mathrm{\Sigma }^2}{G_\eta }}),`$
is the square of the effective Alfv$`\stackrel{´}{\mathrm{e}}`$n Mach number $`M_{\mathrm{Alf}}`$ (where the second equality follows from Eq. (104) below). At the Alfv$`\stackrel{´}{\mathrm{e}}`$n point $`M_{\mathrm{Alf}}=1`$, the numerator $`L\mathrm{\Theta }E`$ should vanish to keep the velocity $`u_\xi `$ finite c86a .

Thus, from Eqs. (81), (86), (87) and (84), given the metric $`g_{\mu \nu }`$ and the conserved functions $`E(\mathrm{\Psi })`$, $`L(\mathrm{\Psi })`$, $`\mathrm{\Omega }(\mathrm{\Psi })`$, $`C(\mathrm{\Psi })`$ and $`S(\mathrm{\Psi })`$, the fluid 4-velocity $`u^\mu `$ can be obtained from the flux function $`\mathrm{\Psi }`$ and its first derivatives $`\mathrm{\Psi }_{,A}`$ if the density $`\rho `$ and the enthalpy $`\mu `$ are additionally known. The expression for $`\rho `$ will be given in the next subsection. The enthalpy $`\mu `$ is then determined as a function of $`\rho `$ and $`\mathrm{\Psi }`$ as will be discussed also in the next subsection.

Once the components of $`u^\mu `$ are known, the $`(2+1)+1`$ decomposition of the fluid 4-velocity is easily performed. With the help of Eqs. (81), (49) and (50), we have

$`u^\mu =u_nn^\mu +u_mm^\mu +u_{\mathrm{\Sigma }}^{}{}_{}{}^{\mu },`$ (95)
where

$`u_n`$ $`=`$ $`N(u_\eta +\mathrm{\Theta }u_\xi ),`$ (96)
$`u_m`$ $`=`$ $`M\left[\left(\mathrm{\Omega }N^\phi \right)u_\eta +\left(1N^\phi \mathrm{\Theta }\right)u_\xi \right],`$ (97)
$`u_{\mathrm{\Sigma }}^{}{}_{}{}^{\mu }`$ $`=`$ $`\stackrel{~}{u}_\mathrm{\Sigma }{}_{}{}^{\mu }(u_\eta +\mathrm{\Theta }u_\xi )N_{\mathrm{\Sigma }}^{}{}_{}{}^{\mu }(u_\xi +\mathrm{\Omega }u_\eta )M_{\mathrm{\Sigma }}^{}{}_{}{}^{\mu }.`$ (98)
### V.2 Density $`\rho `$ and other thermodynamical quantities $`p`$, $`ϵ`$, $`\mu `$ and $`T`$

The pressure $`p`$, the internal energy $`ϵ`$, the enthalpy $`\mu `$ and the temperature $`T`$ are functions of the density $`\rho `$ and the entropy $`S`$ from Eqs. (14), (15) and (17). Hence, given $`S`$ as a function of $`\mathrm{\Psi }`$, the only remaining quantity to be known is the density $`\rho `$.

The density $`\rho `$ is determined by the normalization of the 4-velocity $`u_\mu u^\mu =G_\eta u_{\eta }^{}{}_{}{}^{2}+G_\xi u_{\xi }^{}{}_{}{}^{2}+\stackrel{~}{u}_\mathrm{\Sigma }{}_{A}{}^{}\stackrel{~}{u}_{\mathrm{\Sigma }}^{}{}_{}{}^{A}+2u_\eta \stackrel{~}{u}_\mathrm{\Sigma }{}_{}{}^{\mu }(\eta _\mu +\mathrm{\Omega }\xi _\mu )+2u_\xi \stackrel{~}{u}_\mathrm{\Sigma }{}_{}{}^{\mu }(\xi _\mu +\mathrm{\Theta }\eta _\mu )=1`$, that is,

$`{\displaystyle \frac{(E\mathrm{\Omega }L)^2}{G_\eta \mu ^2}}+{\displaystyle \frac{(4\pi )^2(L\mathrm{\Theta }E)^2}{G_\xi G_{\eta }^{}{}_{}{}^{2}C^4\rho ^2}}\left(1{\displaystyle \frac{4\pi \mu }{G_\eta C^2\rho }}\right)^2+{\displaystyle \frac{H^{AB}\mathrm{\Psi }_{,A}\mathrm{\Psi }_{,B}}{N^2M^2C^2\rho ^2}}+{\displaystyle \frac{\stackrel{~}{N}_\mathrm{\Sigma }^2}{G_\eta }}{\displaystyle \frac{\stackrel{~}{M}_\mathrm{\Sigma }^2}{G_\xi }}=1.`$ (99)
This equation is what is called the wind equation (see c86b  for the circular case). Note that $`\rho `$ contains the first-order derivatives $`\mathrm{\Psi }_{,A}`$ through this equation.

### V.3 Magnetic field $`B^\mu `$

The magnetic field is also calculated from the flux function. With Eqs. (23), (95), (49) and (50), the $`(2+1)+1`$ decomposition of the magnetic field is given by

$`B^\mu =B_nn^\mu +B_mm^\mu +B_\mathrm{\Sigma }{}_{}{}^{\mu },`$ (100)
where

$`B_n`$ $`=`$ $`C\rho N\left[\left(G_\eta u_\eta +\stackrel{~}{N}_\mathrm{\Sigma }\right)\left(u_\eta +\mathrm{\Theta }u_\xi \right)1\right],`$ (101)
$`B_m`$ $`=`$ $`C\rho M\left[\left(G_\eta u_\eta +\stackrel{~}{N}_\mathrm{\Sigma }\right)\left\{\left(\mathrm{\Omega }N^\phi \right)u_\eta +\left(1N^\phi \mathrm{\Theta }\right)u_\xi \right\}+N^\phi \mathrm{\Omega }\right],`$ (102)
$`B_\mathrm{\Sigma }^\mu `$ $`=`$ $`C\rho [(G_\eta u_\eta +\stackrel{~}{N}_\mathrm{\Sigma })u_\mathrm{\Sigma }{}_{}{}^{\mu }+N_\mathrm{\Sigma }{}_{}{}^{\mu }+\mathrm{\Omega }M_\mathrm{\Sigma }{}_{}{}^{\mu }].`$ (103)
The magnetic strength is given by

$`B^2=B^\mu B_\mu =C^2\rho ^2[(G_\eta u_\eta +\stackrel{~}{N}_\mathrm{\Sigma })^2G_\eta ]=C^2\rho ^2[G_\eta (G_\xi u_{\xi }^{}{}_{}{}^{2}+\stackrel{~}{u}_\mathrm{\Sigma }{}_{A}{}^{}\stackrel{~}{u}_{\mathrm{\Sigma }}^{}{}_{}{}^{A}2u_\xi \stackrel{~}{M}_\mathrm{\Sigma })+\stackrel{~}{N}_\mathrm{\Sigma }{}_{}{}^{2}].`$ (104)
The $`(2+1)+1`$ decomposition of the energy-momentum tensor is given in Appendix A.

### V.4 Electric current $`J^\mu `$

Let us consider the following components of the electric current,

$`J^0={\displaystyle \frac{1}{4\pi NM\sqrt{H}}}\left(NM\sqrt{H}F^{0A}\right)_{,A}={\displaystyle \frac{1}{4\pi NM}}\left(NMF^{0A}\right)_{||A},`$ (105)
$`J^3={\displaystyle \frac{1}{4\pi NM\sqrt{H}}}\left(NM\sqrt{H}F^{3A}\right)_{,A}={\displaystyle \frac{1}{4\pi NM}}\left(NMF^{3A}\right)_{||A}.`$ (106)
The field strength tensor components appearing in the above equations are also expressed in terms of the flux function $`\mathrm{\Psi }`$ as

$`F^{0A}`$ $`=`$ $`\left(g^{00}g^{AB}g^{0B}g^{A0}\right)\mathrm{\Omega }\mathrm{\Psi }_{,B}+\left(g^{0B}g^{A3}g^{03}g^{AB}\right)\mathrm{\Psi }_{,B}+\left(g^{01}g^{A2}g^{02}g^{A1}\right)F_{12}`$ (107)
$`=`$ $`{\displaystyle \frac{1}{N^2}}\left[H^{AB}+{\displaystyle \frac{M_\mathrm{\Sigma }^AM_\mathrm{\Sigma }^B}{M^2}}\right]\mathrm{\Omega }\mathrm{\Psi }_{,B}+{\displaystyle \frac{1}{N^2}}\left[N^\phi H^{AB}{\displaystyle \frac{M_\mathrm{\Sigma }^AN_{\mathrm{\Sigma }}^{}{}_{}{}^{B}}{M^2}}\right]\mathrm{\Psi }_{,B}`$
$`\left[\left({\displaystyle \frac{N_{\mathrm{\Sigma }}^{}{}_{}{}^{B}+N^\phi M_\mathrm{\Sigma }^B}{N^2}}\right)ϵ_{B}^{}{}_{}{}^{A}+{\displaystyle \frac{M_\mathrm{\Sigma }^AN_{\mathrm{\Sigma }}^{}{}_{}{}^{B}M_\mathrm{\Sigma }^Cϵ_{BC}}{N^2M^2}}\right]{\displaystyle \frac{F_{12}}{\sqrt{H}}},`$
$`F^{3A}`$ $`=`$ $`\left(g^{30}g^{AB}g^{3B}g^{A0}\right)\mathrm{\Omega }\mathrm{\Psi }_{,B}+\left(g^{3B}g^{A3}g^{33}g^{AB}\right)\mathrm{\Psi }_{,B}+\left(g^{31}g^{A2}g^{32}g^{A1}\right)F_{12}`$ (108)
$`=`$ $`{\displaystyle \frac{1}{N^2}}\left[N^\phi H^{AB}{\displaystyle \frac{N_{\mathrm{\Sigma }}^{}{}_{}{}^{A}M_\mathrm{\Sigma }^B}{M^2}}\right]\mathrm{\Omega }\mathrm{\Psi }_{,B}\left[\left({\displaystyle \frac{1}{M^2}}\left({\displaystyle \frac{N^\phi }{N}}\right)^2\right)H^{AB}{\displaystyle \frac{N_{\mathrm{\Sigma }}^{}{}_{}{}^{A}N_{\mathrm{\Sigma }}^{}{}_{}{}^{B}}{N^2M^2}}\right]\mathrm{\Psi }_{,B}`$
$`+\left[\left\{\left({\displaystyle \frac{1}{M^2}}\left({\displaystyle \frac{N^\phi }{N}}\right)^2\right)M_\mathrm{\Sigma }^B{\displaystyle \frac{N^\phi N_{\mathrm{\Sigma }}^{}{}_{}{}^{B}}{N^2}}\right\}ϵ_{B}^{}{}_{}{}^{A}+{\displaystyle \frac{N_{\mathrm{\Sigma }}^{}{}_{}{}^{A}N_{\mathrm{\Sigma }}^{}{}_{}{}^{B}M_\mathrm{\Sigma }^Cϵ_{BC}}{N^2M^2}}\right]{\displaystyle \frac{F_{12}}{\sqrt{H}}},`$
where $`F_{12}`$ is expressed as

$`{\displaystyle \frac{F_{12}}{\sqrt{H}}}=CNM\rho (u^3\mathrm{\Omega }u^0)=CNM\rho u_\xi \left(1\mathrm{\Omega }\mathrm{\Theta }\right).`$ (109)
In the above, the first equality follows from Eq. (22), and the second from Eqs. (82) and (83). Thus, $`J^0`$ and $`J^3`$ are expressed in terms of $`\mathrm{\Psi }`$ and its first and second derivatives.

### V.5 Auxiliary quantity $`\mathrm{\Lambda }`$

We also need to evaluate the auxiliary quantity $`\mathrm{\Lambda }`$ defined by Eq. (41), that is,

$`\mathrm{\Lambda }={\displaystyle \frac{1}{4\pi }}(u_0B_3u_3B_0).`$ (110)
From the expression of $`B^\mu `$ given by Eq. (23), we have

$`\mathrm{\Lambda }={\displaystyle \frac{1}{4\pi }}C\rho \left[u_0\left(g_{03}+\mathrm{\Omega }g_{33}\right)u_3\left(g_{00}+\mathrm{\Omega }g_{03}\right)\right].`$ (111)
Using the component expressions of $`u^\mu `$ given by Eq. (81), this is rewritten as

$`\mathrm{\Lambda }`$ $`=`$ $`{\displaystyle \frac{1}{4\pi }}C\rho [u_\xi \left(g_{03}{}_{}{}^{2}g_{00}g_{33}\right)(1\mathrm{\Omega }\mathrm{\Theta })+\stackrel{~}{u}_\mathrm{\Sigma }{}_{}{}^{\mu }\eta _{\mu }^{}(g_{03}+\mathrm{\Omega }g_{33})\stackrel{~}{u}_\mathrm{\Sigma }{}_{}{}^{\mu }\xi _{\mu }^{}(g_{00}+\mathrm{\Omega }g_{03})]`$ (112)
$`=`$ $`{\displaystyle \frac{1}{4\pi }}C\rho [u_\xi \left(g_{03}{}_{}{}^{2}g_{00}g_{33}\right)(1\mathrm{\Omega }\mathrm{\Theta })+\stackrel{~}{M}_\mathrm{\Sigma }(g_{00}+\mathrm{\Omega }g_{03})],`$
where the second line follows from Eqs. (49), (50), (80) and (91). The above form is sufficient for $`\mathrm{\Lambda }`$ to be obtained from the flux function, but it can be further simplified if we use Eqs. (80) and (89). From these equations, we find

$$(1\mathrm{\Omega }\mathrm{\Theta })\left(g_{00}g_{33}g_{03}^{}{}_{}{}^{2}\right)=G_\xi \left(g_{00}+\mathrm{\Omega }g_{03}\right).$$

Therefore, we obtain

$`\mathrm{\Lambda }={\displaystyle \frac{1}{4\pi }}(u_0B_3u_3B_0)={\displaystyle \frac{1}{4\pi }}C\rho \left(G_\xi u_\xi \stackrel{~}{M}_\mathrm{\Sigma }\right)\left(g_{00}+\mathrm{\Omega }g_{03}\right).`$ (113)
## VI Grad-Shafranov equation in the covariant form

Now we are ready to show that the GS equation (40) is indeed a second-order differential equation for the flux function $`\mathrm{\Psi }`$. At the same time, following the spirit of the $`(2+1)+1`$ formalism, we express the GS equation in the covariant form with respect to the geometry of $`\mathrm{\Sigma }_{t\phi }`$.

The covariant expression for the GS equation is readily obtained as

$`J^3\mathrm{\Omega }J^0+{\displaystyle \frac{1}{NMC}}ϵ^{AB}(\mu u_{\mathrm{\Sigma }}^{}{}_{A}{}^{})_{||B}\rho (u_\eta +\mathrm{\Theta }u_\xi )\left[E^{}\mathrm{\Lambda }(C\mathrm{\Omega })^{}\right]+\rho (u_\xi +\mathrm{\Omega }u_\eta )\left[L^{}\mathrm{\Lambda }C^{}\right]+\rho TS^{}=0,`$ (114)
where we have replaced $`\sqrt{g}`$, $`u^0`$ and $`u^3`$ in the original GS equation (40) by their $`(2+1)+1`$ type expressions (48), (82) and (83), respectively, and, as before, a double vertical stroke $``$ denotes the covariant differentiation with respect to the 2-metric $`H_{AB}`$.

In the previous section, we have seen that $`u_\eta `$, $`u_\xi `$, $`u_\mathrm{\Sigma }^A`$ and $`\mathrm{\Theta }`$ (section V A), $`\rho `$, $`\mu `$ and $`T`$ (section V B), $`J^0`$ and $`J^3`$ (section V D), and $`\mathrm{\Lambda }`$ (section V E) are all expressed in terms of $`\mathrm{\Psi }`$ and its derivatives, given the conserved functions $`E(\mathrm{\Psi })`$, $`L(\mathrm{\Psi })`$, $`\mathrm{\Omega }(\mathrm{\Psi })`$, $`C(\mathrm{\Psi })`$ and $`S(\mathrm{\Psi })`$, and the metric $`g_{\mu \nu }`$. In particular, we have seen that $`J^0`$ and $`J^3`$ contain the second-order derivatives of $`\mathrm{\Psi }`$, while $`\rho `$ (hence $`\mu `$) as well as $`u_\mathrm{\Sigma }_A`$ contain the first-order derivatives of $`\mathrm{\Psi }`$. Thus, the GS equation (114) is a second-order, non-linear differential equation for $`\mathrm{\Psi }`$, where the first three terms contain the second-order derivatives.

### VI.1 No toroidal field limit

From Eqs. (20) – (22), we find that the toroidal field and the meridional flow vanish if $`|C|\mathrm{}`$. Here note that $`u^3\mathrm{\Omega }u^0=u_\xi (1\mathrm{\Omega }\mathrm{\Theta })C^20`$ in Eq. (22) from Eqs. (82), (83) and (87) (and hence $`\mathrm{\Omega }`$ coincides with the angular velocity $`d\phi /dt=u^3/u^0=\mathrm{\Omega }`$). In the absence of the toroidal field and the meridional flow, a spacetime is circular. The circular limit is expressed as gb93

$`N_\mathrm{\Sigma }{}_{}{}^{A}0,M_\mathrm{\Sigma }{}_{}{}^{A}0.`$ (115)
Therefore, in the $`|C|\mathrm{}`$ limit, the GS equation (114) reduces to

$`J^3\mathrm{\Omega }J^0\rho u_\eta \left[E^{}\mathrm{\Omega }L^{}C\mathrm{\Lambda }\mathrm{\Omega }^{}\right]+\rho TS^{}=0,`$ (116)
where the density is determined by $`\left(E\mathrm{\Omega }L\right)^2/G_\eta \mu ^2=G_\eta u_{\eta }^{}{}_{}{}^{2}=1`$ from Eq. (99), and

$`C\mathrm{\Lambda }={\displaystyle \frac{L\mathrm{\Theta }E}{G_\eta }}\left(g_{00}+\mathrm{\Omega }g_{03}\right)\left(1{\displaystyle \frac{4\pi \mu }{G_\eta C^2\rho }}\right)^1.`$ (117)
Here we regard $`|C|\mathrm{}`$ as the limit of a sequence of models with $`|C|<\mathrm{}`$. The last term $`\left[1(4\pi \mu )/(G_\eta C^2\rho )\right]^1`$ in Eq. (117) can be neglected if the density $`\rho `$ is finite. However, in the case when there is a surface with $`\rho =0`$ like a star and the flux function $`\mathrm{\Psi }`$ is not constant on that surface, the last term $`\left[1(4\pi \mu )/(G_\eta C^2\rho )\right]^1`$ diverges near the surface. Unless one can fine-tune the term $`L\mathrm{\Theta }E`$ so that its zero point cancels this divergence, which seems unlikely to be possible, we should demand the rigid rotation $`\mathrm{\Omega }^{}=0`$ in Eq. (116). This is consistent with Bonazzola et al. b93 ; b95 . Note that if the flux function $`\mathrm{\Psi }`$ is constant on the $`\rho =0`$ surface, we may find $`C`$ that satisfies $`C^2\rho \mathrm{}`$ on the surface.

### VI.2 No poloidal field limit

The poloidal field vanishes if we let $`\mathrm{\Psi }\delta \widehat{\mathrm{\Psi }}`$ and take a limit $`\delta 0`$, as we can see from Eqs. (20) and (21). In this process we relabel the flow lines by $`\widehat{\mathrm{\Psi }}`$ and replace the conserved quantities as $`EE(\widehat{\mathrm{\Psi }})`$, $`LL(\widehat{\mathrm{\Psi }})`$, $`\mathrm{\Omega }\mathrm{\Omega }(\widehat{\mathrm{\Psi }})`$, $`CC(\widehat{\mathrm{\Psi }})`$ and $`SS(\widehat{\mathrm{\Psi }})`$. In the limit $`\delta 0`$, the meridional flow vanishes $`\stackrel{~}{u}_\mathrm{\Sigma }{}_{}{}^{A}0`$ from Eq. (84). Then we can show that the spacetime is circular as expressed in Eq. (115) o02 . Therefore, in the $`\delta 0`$ limit, the GS equation (114) reduces to

$`\left(u_\eta +\mathrm{\Theta }u_\xi \right)\left[E^{}\mathrm{\Lambda }\left(C\mathrm{\Omega }\right)^{}\right]\left(u_\xi +\mathrm{\Omega }u_\eta \right)\left[L^{}\mathrm{\Lambda }C^{}\right]TS^{}=0,`$ (118)
where primes now denote differentiation with respect to $`\widehat{\mathrm{\Psi }}`$. This is an algebraic equation. Here we regard $`\mathrm{\Psi }0`$ as the limit of a sequence of models with $`\mathrm{\Psi }0`$. If $`\mathrm{\Psi }`$ is exactly zero, the transfield components of the Euler equations (12) are satisfied regardless of the GS equation (see Sec. III). Therefore there may exist ’isolated’ solutions which cannot be obtained by the limit discussed here.

### VI.3 No magnetic field limit

There are two limits for configurations with no magnetic field. The first way to obtain such configurations is to let $`\mathrm{\Psi }\delta _1\widehat{\mathrm{\Psi }}`$ and $`C\widehat{C}/\delta _2`$ and take the limit $`\delta _10`$ and $`\delta _20`$, as we can see from Eqs. (20) – (22). Here note that $`u^3\mathrm{\Omega }u^0=u_\xi (1\mathrm{\Omega }\mathrm{\Theta })C^20`$ in Eq. (22) from Eqs. (82), (83) and (87) (and hence $`\mathrm{\Omega }`$ coincides with the angular velocity $`d\phi /dt=u^3/u^0=\mathrm{\Omega }`$). In this process we relabel the flow lines by $`\widehat{\mathrm{\Psi }}`$ and replace the conserved quantities as $`EE(\widehat{\mathrm{\Psi }})`$, $`LL(\widehat{\mathrm{\Psi }})`$, $`\mathrm{\Omega }\mathrm{\Omega }(\widehat{\mathrm{\Psi }})`$, $`CC(\widehat{\mathrm{\Psi }})`$ and $`SS(\widehat{\mathrm{\Psi }})`$. In the limit $`\delta _10`$ and $`\delta _20`$, the meridional flow vanishes $`\stackrel{~}{u}_\mathrm{\Sigma }{}_{}{}^{A}0`$ from Eqs. (20) and (21), so that the spacetime becomes circular as in Eq. (115). Therefore, the GS equation (114) reduces to

$`u_\eta \left[\left(E\mathrm{\Omega }L\right)^{}\mathrm{\Omega }^{}\left(C\mathrm{\Lambda }L\right)\right]TS^{}=0,`$ (119)
where primes denote differentiation with respect to $`\widehat{\mathrm{\Psi }}`$, $`C\mathrm{\Lambda }`$ is given by Eq. (117), and the density is determined from Eq. (99) as

$`{\displaystyle \frac{(E\mathrm{\Omega }L)^2}{G_\eta \mu ^2}}=G_\eta u_\eta {}_{}{}^{2}=1.`$ (120)
Let us see the relation between this limit ($`\mathrm{\Psi }0`$ and $`|C|\mathrm{}`$) and the case of a rotating star b70 ; b93 . From Eq. (120) we have the Bernouilli’s equation for a rotating fluid as

$`\mathrm{ln}\mu \mathrm{ln}u_\eta \mathrm{ln}\left(E\mathrm{\Omega }L\right)=0.`$ (121)
From the definition of $`L`$, Eq. (30), we have $`C\mathrm{\Lambda }L=\mu u_3=u_3u_\eta (E\mathrm{\Omega }L)`$. Then for an isentropic star $`S^{}=0`$, the GS equation in Eq. (119) is written as

$`u_3u_\eta ={\displaystyle \frac{d}{d\mathrm{\Omega }}}\left[\mathrm{ln}(E\mathrm{\Omega }L)\right],`$ (122)
where we have used the fact that $`\mathrm{ln}(E\mathrm{\Omega }L)`$ can be regarded as a function of $`\mathrm{\Omega }`$, since $`E`$, $`L`$ and $`\mathrm{\Omega }`$ are functions of $`\widehat{\mathrm{\Psi }}`$ only. Equation (122) is the well-known integrability condition for a rotating fluid b70 ; b93 .

Note that if we regard $`|C|\mathrm{}`$ as the limit of a sequence of models with $`|C|<\mathrm{}`$, we should demand the rigid rotation $`\mathrm{\Omega }^{}=0`$, as discussed in Sec. VI.1.

The second way to obtain configurations with no magnetic field is to let $`\mathrm{\Psi }\delta \widehat{\mathrm{\Psi }}`$ and $`C\delta \widehat{C}`$ and take the limit $`\delta 0`$. In this process we relabel the flow lines by $`\widehat{\mathrm{\Psi }}`$ and replace the conserved quantities as $`EE(\widehat{\mathrm{\Psi }})`$, $`LL(\widehat{\mathrm{\Psi }})`$, $`\mathrm{\Omega }\mathrm{\Omega }(\widehat{\mathrm{\Psi }})`$, $`\widehat{C}\widehat{C}(\widehat{\mathrm{\Psi }})`$ and $`SS(\widehat{\mathrm{\Psi }})`$. As we can see from Eqs. (20) and (21), there exists a meridional flow in this case. The toroidal field vanishes since in Eq. (22) $`u^3\mathrm{\Omega }u^0=u_\xi (1\mathrm{\Omega }\mathrm{\Theta })`$ and $`u_\xi (L\mathrm{\Theta }E)/G_\xi \mu +\stackrel{~}{M}_\mathrm{\Sigma }/G_\xi `$ in the limit $`\delta 0`$, from Eqs. (82), (83) and (87). The GS equation in Eq. (114) reduces to

$`\rho u_\eta (E^{}\mathrm{\Omega }L^{})+\rho u_\xi (L^{}\mathrm{\Theta }E^{})+{\displaystyle \frac{1}{NM\widehat{C}}}ϵ^{AB}(\mu u_\mathrm{\Sigma }{}_{A}{}^{})_{||B}+\rho TS^{}=0,`$ (123)
where primes denote differentiation with respect to $`\widehat{\mathrm{\Psi }}`$. One can introduce a yet new flux function $`\stackrel{~}{\mathrm{\Psi }}`$ by $`d\stackrel{~}{\mathrm{\Psi }}=d\widehat{\mathrm{\Psi }}/\widehat{C}`$ to absorb the function $`\widehat{C}`$ into the definition of the new flux function $`\stackrel{~}{\mathrm{\Psi }}`$. The resulting equation may be directly obtained from the Euler equations for a perfect fluid.

### VI.4 Newtonian limit

In the Newtonian limit, all physical quantities are expanded in power series of the typical fluid velocity w72 . The metric reduces to

$`g_{\mu \nu }dx^\mu dx^\nu =(1+2\varphi )dt^2+(12\varphi )\left(dZ^2+dR^2+R^2d\phi ^2\right),`$ (124)
where the Newtonian potential $`\varphi `$ is of order $`O(v^2)`$, and we adopt the cylindrical coordinate $`(t,Z,R,\phi )`$. We denote the 3-dimensional velocity by

$`v^i:={\displaystyle \frac{u^i}{u^0}}={\displaystyle \frac{dx^i}{dt}},`$ (125)
where

$`u^0=1\varphi +{\displaystyle \frac{v^2}{2}},`$ (126)
and $`v^2=v^iv_i`$. We regard the internal energy $`ϵ`$ and the pressure $`p`$ to be $`O(v^2)`$. To make the energy density of the electromagnetic field $`O(v^2)`$, we demand

$`B^iO(v),B^0O(v^2),\mathrm{\Psi }O(v),\mathrm{\Omega }O(v),`$ (127)
from Eqs. (8) and (18) – (22).

From Eqs. (20) – (22), we find

$`B^A`$ $`=`$ $`C\rho v^A,`$ (128)
$`B^{\widehat{\phi }}`$ $`=`$ $`C\rho \left(v^{\widehat{\phi }}R\mathrm{\Omega }\right),`$ (129)
where $`B^{\widehat{\phi }}:=RB^3`$ and $`v^{\widehat{\phi }}:=Rv^3`$. From Eqs. (30), (113) and (24), we also have

$`L`$ $`=`$ $`Rv^{\widehat{\phi }}+C\mathrm{\Lambda }=R\left(v^{\widehat{\phi }}{\displaystyle \frac{CB^{\widehat{\phi }}}{4\pi }}\right)O(v),`$ (130)
$`D1`$ $`=`$ $`ϵ+{\displaystyle \frac{p}{\rho }}+{\displaystyle \frac{v^2}{2}}+\varphi R\mathrm{\Omega }v^{\widehat{\phi }}O(v^2).`$ (131)
These results are to be compared with the Newtonian results (note the correspondences between our notation and that of l86  as $`\mathrm{\Omega }G`$, $`C4\pi /F`$, $`LH/F`$, and $`D1J`$).

Let us obtain the Newtonian GS equation, which is of order $`O(v)`$. First consider the $`J^3`$ and $`\mathrm{\Omega }J^0`$ terms in Eq. (114). In the Newtonian order, these terms may be evaluated on the flat background with $`N=1`$ and $`M=R^2`$. Then $`\mathrm{\Omega }J^0`$ is found to be $`O(v^3)`$, and the term $`J^3`$ is given by

$`J^3={\displaystyle \frac{1}{4\pi NM}}\left({\displaystyle \frac{N}{M}}H^{AB}\mathrm{\Psi }_{,B}\right)_{||A}={\displaystyle \frac{1}{4\pi R^2}}\mathrm{\Delta }^{}\mathrm{\Psi },`$ (132)
where

$`\mathrm{\Delta }^{}=R{\displaystyle \frac{}{R}}{\displaystyle \frac{1}{R}}{\displaystyle \frac{}{R}}+{\displaystyle \frac{^2}{Z^2}}.`$ (133)
Next consider the $`ϵ^{AB}(\mu u_\mathrm{\Sigma }{}_{A}{}^{})_B`$ term in Eq. (114). By using Eqs. (98) and (84), we have

$`{\displaystyle \frac{1}{NMC}}ϵ^{AB}(\mu u_\mathrm{\Sigma }{}_{A}{}^{})_{||B}={\displaystyle \frac{1}{NMC}}\left({\displaystyle \frac{\mu }{NMC\rho }}\mathrm{\Psi }_{,A}\right)_{||B}H^{AB}={\displaystyle \frac{1}{4\pi R^2}}[{\displaystyle \frac{4\pi }{C^2\rho }}\mathrm{\Delta }^{}\mathrm{\Psi }+{\displaystyle \frac{4\pi }{C}}\left({\displaystyle \frac{1}{C\rho }}\right)\mathrm{\Psi }].`$ (134)
Finally consider the terms proportional to $`\rho `$ in Eq. (114). To the lowest order, from Eqs. (80), (86), (87) and (113), we have

$`\mathrm{\Theta }`$ $`=`$ $`R^2\mathrm{\Omega }O(v),`$ (135)
$`u_\eta `$ $`=`$ $`1O(1),`$ (136)
$`u_\xi `$ $`=`$ $`{\displaystyle \frac{(LR^2\mathrm{\Omega })}{R^2}}\left({\displaystyle \frac{4\pi }{C^2\rho }}\right)\left(1{\displaystyle \frac{4\pi }{C^2\rho }}\right)^1O(v),`$ (137)
$`E^{}`$ $`=`$ $`(D1+\mathrm{\Omega }L)^{}O(v),`$ (138)
$`L^{}`$ $``$ $`O(1),`$ (139)
$`\mathrm{\Lambda }`$ $`=`$ $`{\displaystyle \frac{C\rho u_\xi R^2}{4\pi }}O(v).`$ (140)
Then we can show

$`\rho u_\eta \left[E^{}\mathrm{\Lambda }C\mathrm{\Omega }^{}\mathrm{\Omega }L^{}\right]=\rho \left[(D1)^{}+Rv^{\widehat{\phi }}\mathrm{\Omega }^{}\right].`$ (141)
From Eqs. (140) and (130), we may express $`\rho u_\xi `$ as

$`\rho u_\xi ={\displaystyle \frac{4\pi \mathrm{\Lambda }}{CR^2}}={\displaystyle \frac{4\pi }{CR^2}}{\displaystyle \frac{LRv^{\widehat{\phi }}}{C}}.`$ (142)
Using this expression, we can show

$`\rho u_\xi \left[L^{}\mathrm{\Lambda }C^{}\right]`$ $`=`$ $`{\displaystyle \frac{1}{CR^2}}\left({\displaystyle \frac{4\pi L}{C}}Rv^{\widehat{\phi }}{\displaystyle \frac{4\pi }{C}}\right){\displaystyle \frac{C}{4\pi }}\left[\left({\displaystyle \frac{4\pi L}{C}}\right)^{}L\left({\displaystyle \frac{4\pi }{C}}\right)^{}+C\mathrm{\Lambda }\left({\displaystyle \frac{4\pi }{C}}\right)^{}\right]`$ (143)
$`=`$ $`{\displaystyle \frac{1}{4\pi R^2}}\left({\displaystyle \frac{4\pi L}{C}}Rv^{\widehat{\phi }}{\displaystyle \frac{4\pi }{C}}\right)\left[\left({\displaystyle \frac{4\pi L}{C}}\right)^{}Rv^{\widehat{\phi }}\left({\displaystyle \frac{4\pi }{C}}\right)^{}\right],`$
where the second equality follows from Eq. (130). It is also easy to show that $`\rho \mathrm{\Theta }u_\xi \left[E^{}\mathrm{\Lambda }(C\mathrm{\Omega })^{}\right]O(v^3)`$.

Therefore from Eqs. (132), (134), (141) and (143), the GS equation in the Newtonian limit is given by

$`\left(1{\displaystyle \frac{4\pi }{C^2\rho }}\right)\mathrm{\Delta }^{}\mathrm{\Psi }{\displaystyle \frac{4\pi }{C}}\left({\displaystyle \frac{1}{C\rho }}\right)\mathrm{\Psi }`$
$`=4\pi \rho R^2\left[(D1)^{}+Rv^{\widehat{\phi }}\mathrm{\Omega }^{}\right]\left({\displaystyle \frac{4\pi L}{C}}Rv^{\widehat{\phi }}{\displaystyle \frac{4\pi }{C}}\right)\left[\left({\displaystyle \frac{4\pi L}{C}}\right)^{}Rv^{\widehat{\phi }}\left({\displaystyle \frac{4\pi }{C}}\right)^{}\right]+4\pi R^2\rho TS^{}.`$ (144)
This is equivalent to the Newtonian GS equation in l86 .

## VII Summary

We have derived the GS equation (114) in noncircular (the most general) stationary axisymmetric spacetimes. The GS equation has been given in the covariant form projected onto the $`t=`$const. and $`\phi =`$const. 2-surface $`\mathrm{\Sigma }_{t\phi }`$. We have also derived the wind equation (99) in noncircular spacetimes. We have discussed various limits of the GS equation (no toroidal field limit, no poloidal field limit, no magnetic field limit and Newtonian limit).

To obtain equilibrium configurations of magnetars, we have to solve the GS equation (114). As first glance, it looks formidable to solve it. One possibility is to take a perturbative approach to solve the GS equation. Unless the magnetic field is as strong as the maximum magnetic field allowed by the virial theorem $`10^{18}`$ G b95 ; bg96 , we may assume weak magnetic fields compared with gravity. Then the magnetic field may be treated as a small perturbation on an already-known non-magnetized configuration. This approach is similar to that developed for slowly rotating stars c33 ; h67 , in which the perturbation parameter is the angular velocity. Work in this direction is in progress. The preliminary study indicates that the degree of noncircularity of the spacetime in a neutron star with mass $`M_{}`$ and radius $`R_{}`$, is about $`(N_{\mathrm{\Sigma }}^{}{}_{}{}^{\mu }N_{\mathrm{\Sigma }}^{}{}_{\mu }{}^{})^{1/2}(M_{}/R_{})v_{\mathrm{mf}}`$ and $`(M_{\mathrm{\Sigma }}^{}{}_{}{}^{\mu }M_{\mathrm{\Sigma }}^{}{}_{\mu }{}^{})^{1/2}0.1`$-$`0.01(M_{}/R_{})_M`$, where $`v_{\mathrm{mf}}`$ is the velocity of the meridional flow, $`_M`$ is the ratio of the magnetic energy to the gravitational energy, and the length scale is normalized by the mass Ioka2 .

###### Acknowledgements.

KI is grateful to F. Takahara for useful discussions. This work was supported in part by the Monbukagaku-sho Grant-in-Aid for Scientific Research, Nos. 00660 (KI) and 14047214 (MS).

## Appendix A $`(2+1)+1`$ decomposition of the energy-momentum tensor

The $`(2+1)+1`$ decomposition of the energy-momentum tensor is

$`T^{\mu \nu }`$ $`=`$ $`en^\mu n^\nu +j(n^\mu m^\nu +m^\mu n^\nu )+j^A\left(n^\mu H_A{}_{}{}^{\nu }+H_A{}_{}{}^{\mu }n_{}^{\nu }\right)`$ (145)
$`+sm^\mu m^\nu +s^A\left(m^\mu H_A{}_{}{}^{\nu }+H_A{}_{}{}^{\mu }m_{}^{\nu }\right)+s^{AB}H_A{}_{}{}^{\mu }H_{B}^{}{}_{}{}^{\nu },`$
where $`n_\mu `$ and $`m_\mu `$ are the unit timelike and spacelike normals to the 2-surface $`\mathrm{\Sigma }_{t\phi }`$, respectively, and $`A,B=1,2`$. For the electromagnetic field in an ideal MHD system, we have from Eqs. (8) and (9), $`F^{\mu \alpha }F^\nu {}_{\alpha }{}^{}=(u^\mu u^\nu +g^{\mu \nu })B^2B^\mu B^\nu `$ and $`F^{\mu \nu }F_{\mu \nu }=2B^2`$, and hence

$`{\displaystyle \frac{1}{4\pi }}\left(F^{\mu \alpha }F_{}^{\nu }{}_{\alpha }{}^{}{\displaystyle \frac{1}{4}}g^{\mu \nu }F^{\alpha \beta }F_{\alpha \beta }\right)={\displaystyle \frac{1}{4\pi }}\left[\left(u^\mu u^\nu +{\displaystyle \frac{1}{2}}g^{\mu \nu }\right)B^2B^\mu B^\nu \right].`$ (146)
Then, using Eqs (95) and (100) the components of the energy-momentum tensor in Eq. (145) are obtained as

$`e`$ $`=`$ $`T^{\mu \nu }n_\mu n_\nu =\left(\rho +\rho ϵ+p\right)\left(u_n\right)^2p+{\displaystyle \frac{1}{4\pi }}\left[\left\{\left(u_n\right)^2{\displaystyle \frac{1}{2}}\right\}B^2\left(B_n\right)^2\right],`$ (147)
$`j`$ $`=`$ $`T^{\mu \nu }n_\mu m_\nu =\left(\rho +\rho ϵ+p\right)u_nu_m+{\displaystyle \frac{1}{4\pi }}\left(u_nu_mB^2B_nB_m\right),`$ (148)
$`j^A`$ $`=`$ $`H^A{}_{\mu }{}^{}T_{}^{\mu \nu }n_\nu =(\rho +\rho ϵ+p)u_nu_\mathrm{\Sigma }{}_{}{}^{A}+{\displaystyle \frac{1}{4\pi }}(B^2u_nu_\mathrm{\Sigma }{}_{}{}^{A}B_nB_\mathrm{\Sigma }{}_{}{}^{A}),`$ (149)
$`s`$ $`=`$ $`T^{\mu \nu }m_\mu m_\nu =\left(\rho +\rho ϵ+p\right)\left(u_m\right)^2+p+{\displaystyle \frac{1}{4\pi }}\left[\left\{\left(u_m\right)^2+{\displaystyle \frac{1}{2}}\right\}B^2\left(B_m\right)^2\right],`$ (150)
$`s^A`$ $`=`$ $`H^A{}_{\mu }{}^{}T_{}^{\mu \nu }m_\nu =(\rho +\rho ϵ+p)u_mu_\mathrm{\Sigma }{}_{}{}^{A}+{\displaystyle \frac{1}{4\pi }}(B^2u_mu_\mathrm{\Sigma }{}_{}{}^{A}B_mB_\mathrm{\Sigma }{}_{}{}^{A}),`$ (151)
$`s^{AB}`$ $`=`$ $`H^A{}_{\mu }{}^{}H_{}^{B}{}_{\nu }{}^{}T_{}^{\mu \nu }=(\rho +\rho ϵ+p)u_\mathrm{\Sigma }{}_{}{}^{A}u_{\mathrm{\Sigma }}^{}{}_{}{}^{B}+pH^{AB}+{\displaystyle \frac{1}{4\pi }}[\left(u_\mathrm{\Sigma }{}_{}{}^{A}u_{\mathrm{\Sigma }}^{}{}_{}{}^{B}+{\displaystyle \frac{1}{2}}H^{AB}\right)B^2B_\mathrm{\Sigma }{}_{}{}^{A}B_{\mathrm{\Sigma }}^{}{}_{}{}^{B}],`$ (152)
where $`u_n`$, $`u_m`$, $`u_\mathrm{\Sigma }^\mu `$, $`B_n`$, $`B_m`$, $`B_\mathrm{\Sigma }^\mu `$ and $`B^2`$ are given by Eqs. (96) – (98) and (101) – (104), respectively.

## Appendix B Symbols

Here, we summarize definitions of some of the symbols we use, which may not be commonly used, with the equation numbers where they are defined or introduced.

$``$ Quantities conserved along each flow line:

* Flux function, $`\mathrm{\Psi }:=A_\mu \xi ^\mu =A_3`$.
* Electric potential, $`\mathrm{\Phi }:=A_\mu \eta ^\mu =A_0`$.
* ‘Magnetic field strength’ relative to the magnitude of meridional flow; Eqs. (20) – (22).
* ‘Fluid energy’ per unit mass; Eqs. (24) and (27).
* ‘Energy’ per unit mass; Eqs. (25) and (29).
* ‘Angular momentum’ per unit mass; Eqs. (26) and (30).
* ‘Angular velocity’ of the magnetic field line, $`\mathrm{\Omega }=d\mathrm{\Phi }/d\mathrm{\Psi }`$; Eq. (19).
* Entropy per unit mass; Eqs. (15).

$``$ Quantities associated with the metric:

* Killing vector associated with stationarity, $`\eta ^\mu =(/t)^\mu `$; Eq. (49).
* Killing vector associated with axisymmetry, $`\xi ^\mu =(/\phi )^\mu `$; Eq. (50).
* Unit timelike 4-vector orthogonal to $`t=`$const. hypersurface $`\mathrm{\Sigma }_t`$; Eqs. (42).
* 3-metric on $`\mathrm{\Sigma }_t`$; Eqs. (44).
* Unit spacelike 4-vector orthogonal to $`t=`$const. and $`\phi =`$const. hypersurface $`\mathrm{\Sigma }_{t\phi }`$; Eqs. (45).
* 2-metric on $`\mathrm{\Sigma }_{t\phi }`$; Eq. (47).
* Norm of $`\eta ^\mu +\mathrm{\Omega }\xi ^\mu `$; Eq. (88).
* Quantity such that $`\xi ^\mu +\mathrm{\Theta }\eta ^\mu `$ is orthogonal to $`\eta ^\mu +\mathrm{\Omega }\xi ^\mu `$; Eq. (80).
* Norm of $`\xi ^\mu +\mathrm{\Theta }\eta ^\mu `$; Eq. (89).
* Lapse function, $`N=\eta ^\mu n_\mu `$; Eqs. (42).
* Shift vector, $`N^\mu =Nn^\mu \eta ^\mu `$; Eq. (49).
* Projection of $`N^\mu `$ onto $`\mathrm{\Sigma }_{t\phi }`$; Eq. (49).
* $`=m_\mu \xi ^\mu `$; Eqs. (45).
* $`=Mm^\mu \xi ^\mu `$; Eq. (50).

$``$ Quantities associated with the fluid:

* Fluid 4-velocity; Eqs. (10).
* Rest mass density; Eqs. (10).
* Enthalpy per unit mass; Eq. (14).
* Internal energy per unit mass; Eq. (10).
* Pressure; Eq. (10).
* Temperature; Eq. (15).
* $`=(u^0\mathrm{\Theta }u^3)/(1\mathrm{\Omega }\mathrm{\Theta })`$; Eq. (81), or Eqs. (82) and (83).
* $`=(u^3\mathrm{\Omega }u^0)/(1\mathrm{\Omega }\mathrm{\Theta })`$; Eq. (81), or Eqs. (82) and (83).
* $`=n_\mu u^\mu `$; Eqs. (95) and (96).
* $`=m_\mu u^\mu `$; Eqs. (95) and (97).
* Projection of $`u^\mu `$ onto $`\mathrm{\Sigma }_{t\phi }`$; Eqs. (95) and (98).
* $`u^A`$ ($`A=1,2`$) components of $`u^\mu `$; Eq. (81).
* A component of $`\stackrel{~}{u}_\mathrm{\Sigma }^\mu `$ defined by Eq. (90).
* A component of $`\stackrel{~}{u}_\mathrm{\Sigma }^\mu `$ defined by Eq. (91).

$``$ Quantities associated with the electromagnetic field:

* Electric field in the fluid rest frame; Eq. (6).
* Magnetic field in the fluid rest frame; Eq. (7).
* $`=n_\mu B^\mu =n_\mu u_\nu ϵ^{\mu \nu \alpha \beta }F_{\alpha \beta }/2`$; Eq. (100).
* $`=m_\mu B^\mu =u_\mu m_\nu ϵ^{\mu \nu \alpha \beta }F_{\alpha \beta }/2`$; Eq. (100).
* Projection of $`B^\mu `$ onto $`\mathrm{\Sigma }_{t\phi }`$, Eq. (100).
* Electromagnetic current 4-vector; Eq. (4).

$``$ Others:

* Alfv$`\stackrel{´}{\mathrm{e}}`$n Mach number; Eq. (94).
* An auxiliary quantity defined by Eq. (41).