mean molecular weight. We fix $\mu_e=2$ in all our computations here, which corresponds to a fully ionized pure hydrogen gas. This choice of parameters is same as that in [@Hachisu_1986].

The form of the current density and the boundary condition {#sec:bc}
----------------------------------------------------------

From the assumptions of axisymmetry and stationarity, we introduce magnetic flux function $\Psi$ as follows: $$\begin{aligned}
   H_R \equiv -\frac{1}{R}\P{\Psi}{z} , \hspace{10pt}
   H_z \equiv \frac{1}{R}\P{\Psi}{R} \ ,\end{aligned}$$ where $H_R$ and $H_z$ are magnetic field components in the $R$-direction and $z$-direction, respectively. We assume this flux function is positive in the entire space. By introducing this magnetic flux function, equation (\[Eq:divH\]) can be automatically satisfied. It should be noted that the magnetic flux function $\Psi$ can be expressed as: $$\begin{aligned}
    \Psi = r \sin \theta A_{\varphi} \ , \end{aligned}$$ where $A_{\varphi}$ is the $\varphi$-component of the vector potential $\Vec{A} = (A_R, A_{\varphi}, A_z)$.

As shown in [@Tomimura_Eriguchi_2005], for axisymmetric and stationary barotropes with rigid rotation we can constrain the form of the electric current density by using an integrability condition of the equations of motion, equation (\[Eq:eular\]): $$\begin{aligned}
  \frac{\Vec{j}}{c} = \frac{1}{4\pi} \D{\kappa(\Psi)}{\Psi} \Vec{H} + r \sin \theta \rho 
   \mu(\Psi) \Vec{e}_\varphi \ ,
   \label{Eq:current}\end{aligned}$$ where $\kappa(\Psi)$ and $\mu(\Psi)$ are arbitrary functions 