Computation of Equilibria of a Flowing Two-fluid in Two Dimensions

Abstract

Two-fluid equilibria encompass more physics than the Grad–Shafranov (GS) equation for static equilibria, or even the flowing MHD model. However the two-fluid system is more complicated, and, worse yet, is a singular perturbation problem. The latter difficulty is overcome using the “nearby-fluids” ordering. A “1.5D” solution method has been used to interpret results from the TCS experiment. These results, summarized here, exhibit trends indicative of the improved stability and transport observed experimentally. An algorithm for solving 2D equilibria has been developed based on a relaxation method. The magnetic flux function, governed by the extended GS equation, is updated by successive-overrelaxation, while the toroidal field and flow components and the density are updated using a Newton–Raphson-like method.

Appendix: Equations for nearby-fluids equilibria

The equations that follow use dimensionless variables based on the reference quantities for the magnetic field B _R (scale of the poloidal field), the density n _R (peak density), and the system length scale L (e.g. the separatrix radius). From these follow reference scales for other variables: the Alfven speed $V_{\rm AR} \equiv B_{\rm R}/(4\pi m_{i}n_{\rm R})^{1/2}$; temperature $T_{\rm R}\equiv B_{\rm R}^2 /4\pi n_{\rm R} k$; and the electrostatic potential ϕ_R≡ kT _R/e. Cylindrical coordinates (r,z) are used.

In the nearby-fluids ordering, the “arbitrary” surface functions for the ion and electron stream functions $\bar{\psi}_i$ and $\bar{\psi}_e$ are replaced by F, G:

$$ \bar{\psi}_e (x) = \int\limits_0^x F({x}^{\prime})d{x}^{\prime} \quad \bar{\psi}_i (x) = \int\limits_0^x F({x}^{\prime})d{x}^{\prime}+\varepsilon G(x) $$

(A1,2)

The other “arbitrary” functions are unchanged: total enthalpies H _i (...), H _e (...) and entropies S _i (...), S _e (...) for each species. The extended GS equation for the magnetic flux variable ψ (r,z) is

$$ \Delta^\ast \psi = \frac{1}{\varepsilon}\left[rB_\theta F(\psi)-rnu_\theta \right]-nr^2\left[H^{\prime}_e (\psi)-T_e S^{\prime}_e (\psi)\right] $$

(A3)

where $\Delta^\ast \equiv r^2\nabla \cdot [(1/r)\nabla]$ is the Grad–Shafranov operator, and H _e , S _e are surface functions representing the electron total enthalpy and entropy. The ion surface variable is

$$ Y = \psi + \varepsilon ru_{\theta} $$

(A4)

The toroidal field and flow components are

$$ B_\theta = \frac{1}{r}\frac{G(\hbox{Y})-\varepsilon \left\{\bar{F}r^2\left[H^{\prime}_i (Y)-T_i S^{\prime}_i (Y)\right]-M_{A2}^2 \Delta_n^\ast\psi_{i} \right\}}{1-M_{A2}^2} $$

(A5)

$$ u_\theta = \frac{1}{rn}\left[\bar{\psi}^{\prime} rB_\theta+ \varepsilon \left\{nr^2\left[H^{\prime}_i (Y)-T_i S^{\prime}_i (Y)\right]- \bar{\psi}^{\prime}_i \Delta_n^\ast \psi_i \right\}\right] $$

(A6)

where and H _i , S _i are surface functions representing the ion total enthalpy and entropy. The “averaged” $\bar{F}$ and the two-fluid Alfven number (“Mach number” using the Alfven speed) are

$$ \bar{F} \equiv (1/2)\left[F(\psi +0.789\varepsilon ru_\theta)+ F(\psi +0.211\varepsilon ru_\theta)\right]+O(\varepsilon^4) \quad M_{A2} \equiv \sqrt{\left|\bar{F}\bar{\psi}^{\prime}_{i}/n\right|}, $$

(A7,8)

respectively, were ${\Delta_n^\ast \equiv r^2\nabla \cdot [(1/nr)\nabla]}$ is the density-weighted form of the GS operator. The density is governed by the combined Bernoulli equation (sum of ion and electron Bernoulli equations): it is updated by a Newton–Raphson iterative method as follows

$$ n_{\rm new} = n\left[1+\frac{H_i (Y)+H_e (\psi)-f_B} {\gamma (T_i +T_e)(1-M_{s2}^2)}\right] $$

(A9)

where the two-fluid Mach number and the Bernoulli function are

$$ M_{s2} \equiv \sqrt{\left|\left\{u_p^2 +u_\theta^2 + \varepsilon ru_\theta \left[H^{\prime}_e (\psi)-T_e S^{\prime}_e (\psi)\right]\right\}/\gamma (T_i +T_e)\right|} $$

(A10)

$$ f_B \equiv \gamma (T_i +T_e)/(\gamma-1)-u^2/2, $$

(A11)

respectively. Equations of state give the temperatures

$$ T_i = n^{\gamma -1}e^{(\gamma -1)S_i (Y)} \quad T_e = n^{\gamma -1}e^{(\gamma -1)S_e (\psi)} $$

(A12,13)

The electrostatic potential is given by the electron Bernoulli equation

$$ \phi = \gamma T_e/(\gamma -1)-H_e (\psi) $$

(A14)

The poloidal field and flow are

$$ B_p = \left|\nabla \psi\right|/r \quad u_p = \left|\nabla \psi_i (Y) \right|/ nr $$

(A15,16)

The solution method is as follows. The arbitrary functions F, G,H _i , H _e , S _i , S _e are pre-specified. The first guess for ψ (r,z), n(r,z) is made. In the iteration cycle B _θ, u _θ, and n are updated first (A5,6,9) after which ψ is updated by successive-overrelaxation with a finite-difference expression for Δ^* (A3). The iteration is repeated until the solution converges. If the intent is to model experiment, then the arbitrary functions are adjusted and the procedure repeated.