Grad–Shafranov equation

The Grad–Shafranov equation (H. Grad and H. Rubin (1958); Vitalii Dmitrievich Shafranov (1966)) is the equilibrium equation in ideal magnetohydrodynamics (MHD) for a two dimensional plasma, for example the axisymmetric toroidal plasma in a tokamak. This equation is a two-dimensional, nonlinear, elliptic partial differential equation obtained from the reduction of the ideal MHD equations to two dimensions, often for the case of toroidal axisymmetry (the case relevant in a tokamak). The flux function \psi is both a dependent and an independent variable in this equation:

\Delta^{*}\psi = -\mu_{0}R^{2}\frac{dp}{d\psi}-\frac{1}{2}\frac{dF^2}{d\psi}

where \mu_0 is the magnetic permeability, p(\psi) is the pressure, F(\psi)=RB_{\phi}

and the magnetic field and current are given by

\vec{B}=\frac{1}{R}\nabla\psi\times \hat{e}_{\phi}+\frac{F}{R}\hat{e}_{\phi}
\mu_0\vec{J}=\frac{1}{R}\frac{dF}{d\psi}\nabla\psi\times \hat{e}_{\phi}-\frac{1}{R}\Delta^{*}\psi \hat{e}_{\phi}

The elliptic operator

\Delta^{*} is given by

\Delta^{*}\psi \equiv R^{2} \vec{\nabla} \cdot \left( \frac{1}{R^{2}} \vec{\nabla} \psi \right) = R\frac{\partial}{\partial R}\left(\frac{1}{R}\frac{\partial \psi}{\partial R}\right)+\frac{\partial^2 \psi}{\partial Z^2}.

The nature of the equilibrium, whether it be a tokamak, reversed field pinch, etc. is largely determined by the choices of the two functions F(\psi) and p(\psi) as well as the boundary conditions.

Derivation (in slab coordinates):
To begin we assume that the system is 2-dimensional with z as the invariant axis, i.e. \partial /\partial z = 0 for all quantities. Then the magnetic field can be written in cartesian coordinates as

 \bold{B} = (\partial A/\partial y,-\partial A /\partial x,B_z(x,y))

or more compactly,

 \bold{B} =\nabla A \times \hat{\bold{z}} + B_z \hat{\bold{z}},

where A(x,y)\hat{\bold{z}} is the vector potential for the in-plane (x and y components) magnetic field. Note that based on this form for B we can see that A is constant along any given magnetic field line, since \nabla A is everywhere perpendicular to B. (Also note that -A is the flux function \psi mentioned above.)

Two dimensional, stationary, magnetic structures are described by the balance of pressure forces and magnetic forces, i.e.:

\nabla p = \bold{j} \times \bold{B},

where p is the plasma pressure and j is the electric current. Note from the form of this equation that we also know p is a constant along any field line, (again since \nabla p is everywhere perpendicular to B. Additionally, the two-dimensional assumption (\partial / \partial z ) means that the z- component of the left hand side must be zero, so the z-component of the magnetic force on the right hand side must also be zero. This means that \bold{j}_\perp \times \bold{B}_\perp = 0, i.e. \bold{j}_\perp is parallel to \bold{B}_\perp.

We can break the right hand side of the previous equation into two parts:

\bold{j} \times \bold{B} = j_z (\hat{\bold{z}} \times \bold{B_\perp}) +\bold{j_\perp} \times \hat{\bold{z}}B_z ,

where the \perp subscript denotes the component in the plane perpendicular to the z-axis. The z component of the current in the above equation can be written in terms of the one-dimensional vector potential as j_z = -\nabla^2 A/\mu_0. . The in plane field is

\bold{B}_\perp = \nabla A \times \hat{\bold{z}} ,

and using Maxwell–Ampère's equation, the in plane current is given by

\bold{j}_\perp = (1/\mu_0)\nabla B_z \times \hat{\bold{z}}.

In order for this vector to be parallel to \bold{B}_\perp as required, the vector \nabla B_z must be perpendicular to \bold{B}_\perp, and B_z must therefore, like p be a field-line invariant.

Rearranging the cross products above, we see that

\hat{\bold{z}} \times \bold{B}_\perp = \nabla A -(\bold{\hat z} \cdot\nabla A) \bold{\hat z} = \nabla A,

and

\bold{j}_\perp \times B_z\bold{\hat{z}} = -(1/\mu_0)B_z\nabla B_z +(B_z/\mu_0)(\bold{\hat z}\cdot\nabla B_z)\bold{\hat z}=-(1/\mu_0) B_z\nabla B_z.

These results can be substituted into the expression for \nabla p to yield:

\nabla p = -\left[(1/\mu_0) \nabla^2 A\right]\nabla A-(1/\mu_0)B_z\nabla B_z.

Now, since p and B_z are constants along a field line, and functions only of A, we note that \nabla p = (d p /dA)\nabla A and  \nabla B_z = (d B_z/dA)\nabla A. Thus, factoring out \nabla A and rearraging terms we arrive at the Grad–Shafranov equation:

\nabla^2 A = -\mu_0 \frac{d}{dA}\left(p + \frac{B_z^2}{2\mu_0}\right)

References