Lecture 13 The Shafranov Shift

← Previous Lecture Contents Next Lecture →

Lecture 13
The Shafranov Shift

Overview

Why this is a classic problem. The Shafranov-shift calculation is the first place where the Grad–Shafranov equation of Lecture 12 becomes a genuinely tokamak problem rather than a dressed-up cylinder. In the large-aspect-ratio limit,
1.
the zeroth-order equilibrium is the screw-pinch force balance already encountered in Lecture 11,
2.
the first toroidal correction displaces the magnetic axis outward by an amount of order \(a^2/R_0\),
3.
the same calculation shows why an external vertical field is required to hold the plasma in toroidal force balance, and
4.
the key dimensionless combinations \(\beta _p\) and \(\ell _i\) emerge naturally.

The Shafranov shift is one of the best examples in MHD of a calculation that is at once elegant, practical, and historically formative. It is elegant because the algebra can be carried explicitly in a controlled expansion in inverse aspect ratio. It is practical because the result connects directly with magnetic measurements, diamagnetic loops, and the vertical-field coils of a tokamak. And it is historically formative because it taught the field, very early, that toroidal curvature produces an outward “hoop” force that cannot be balanced by plasma currents alone.

Historical Perspective

Shafranov’s analysis of toroidal plasma equilibrium in the late 1950s and early 1960s made clear that a torus is not just a bent cylinder: the toroidal curvature of the current channel creates an outward force that must be balanced by externally applied fields Shafranov [1960, 1966]. The modern language of “Shafranov shift,” “vertical field,” and the combinations \(\beta _p+\ell _i/2\) all trace back to this line of work. In present-day tokamaks the same physics is embedded in equilibrium solvers and reconstruction codes, but the large-aspect-ratio calculation remains the cleanest place to see where it comes from.

We begin from the isotropic MHD equilibrium implied by the static limit of the momentum equation, Eq. (1.8), and specialize to the axisymmetric Grad–Shafranov equation already derived in Lecture 12,

\Delta ^\star \psi = -\muo R^2 \pp {p}{\psi } - F(\psi )\pp {F}{\psi }. \tag{13.1}

The purpose of this lecture is to solve Eq. (13.1) asymptotically for a large-aspect-ratio tokamak with circular flux surfaces.

Caution

Fixed boundary versus free boundary. Two related but distinct questions are mixed together in many first presentations of this problem.
1.
In a fixed-boundary calculation one imagines the plasma edge held at a prescribed circular boundary and solves for the internal shift of the flux surfaces.
2.
In a free-boundary calculation one asks what external vertical field is required so that the plasma does not simply drift radially outward.

The algebra below treats both, and keeping the distinction straight is important.

13.1 Large-Aspect-Ratio Expansion of the Grad–Shafranov Equation

Local coordinates about the magnetic axis. We work in cylindrical coordinates \((R,\phi ,Z)\) and introduce local polar coordinates \((r,\theta )\) centered on the magnetic axis,

R = R_0 + r\cos \theta , \qquad Z = r\sin \theta , \qquad \epsilon \equiv \frac {r}{R_0}\ll 1. \tag{13.2}

The coordinate derivatives needed below are

\begin{aligned}\pp {R}{r} &= \cos \theta , & \pp {Z}{r} &= \sin \theta , \\ \pp {R}{\theta } &= -r\sin \theta , & \pp {Z}{\theta } &= r\cos \theta ,\end{aligned}

and therefore

\begin{aligned}\pp {r}{R} &= \cos \theta , & \pp {r}{Z} &= \sin \theta , \\ \pp {\theta }{R} &= -\frac {\sin \theta }{r}, & \pp {\theta }{Z} &= \frac {\cos \theta }{r}.\end{aligned}

The Grad–Shafranov operator in local form. The operator in Eq. (13.1) is

\Delta ^\star \psi = \pp {^2 \psi }{R^2} - \frac {1}{R}\pp {\psi }{R} + \pp {^2 \psi }{Z^2}. \tag{13.7}

The first and third terms form the planar Laplacian. Transforming them to the \((r,\theta )\) coordinates gives

\pp {^2 \psi }{R^2} + \pp {^2 \psi }{Z^2} = \frac {1}{r}\pp {}{r}\!\left (r\pp {\psi }{r}\right ) + \frac {1}{r^2}\pp {^2 \psi }{\theta ^2}. \tag{13.8}

For the remaining term we use the chain rule:

\pp {\psi }{R} = \pp {\psi }{r}\pp {r}{R} + \pp {\psi }{\theta }\pp {\theta }{R} = \cos \theta \,\pp {\psi }{r} - \frac {\sin \theta }{r}\pp {\psi }{\theta }. \tag{13.9}

Substituting Eqs. (13.8) and (13.9) into Eq. (13.7) yields the exact local form

\boxed { \Delta ^\star \psi = \frac {1}{r}\pp {}{r}\!\left (r\pp {\psi }{r}\right ) + \frac {1}{r^2}\pp {^2 \psi }{\theta ^2} - \frac {1}{R_0(1+\epsilon \cos \theta )} \left ( \cos \theta \,\pp {\psi }{r} - \frac {\sin \theta }{r}\pp {\psi }{\theta } \right ). } \tag{13.10}

Ordering of the equilibrium quantities. We now assume a conventional large-aspect-ratio tokamak ordering. If the safety factor is order unity, \[ q \sim \frac {r B_\phi }{R_0 B_\theta }\sim 1, \] then

B_\theta \sim \epsilon B_\phi . \tag{13.11}

If the poloidal beta is order unity, \[ \beta _p \sim \frac {2\muo p}{B_\theta ^2} \sim 1, \] then the pressure scales like

p \sim \frac {B_\theta ^2}{2\muo } \sim \epsilon ^2\frac {B_\phi ^2}{2\muo }. \tag{13.12}

Thus the equilibrium is close to a cylindrical screw pinch at leading order, with toroidal curvature entering as an \(\mathcal O(\epsilon )\) correction.

Expansion of the flux function. We seek a regular expansion

\psi (r,\theta ) = \psi _0(r) + \psi _1(r,\theta ) + \mathcal O(\epsilon ^2), \qquad \psi _1=\mathcal O(\epsilon ), \tag{13.13}

with a circular plasma boundary to lowest order. The free functions are expanded about \(\psi _0\):

\begin{aligned}\pp {p}{\psi } &= \left .\pp {p}{\psi }\right |_{\psi _0} + \left .\pp {^2 p}{\psi ^2}\right |_{\psi _0}\psi _1 +\cdots , \\ F\pp {F}{\psi } &= \left .\frac {d}{d\psi }\left (\frac {F^2}{2}\right )\right |_{\psi _0} + \left .\frac {d^2}{d\psi ^2}\left (\frac {F^2}{2}\right )\right |_{\psi _0}\psi _1 +\cdots .\end{aligned}

13.2 Zeroth Order: The Cylindrical Equilibrium

Order \(\mathcal O(1)\) Grad–Shafranov equation. Collecting the \(\mathcal O(1)\) terms of Eq. (13.10) gives

\frac {1}{r}\dd {}{r}\!\left (r\dd {\psi _0}{r}\right ) = -\muo R_0^2 \left .\pp {p}{\psi }\right |_{\psi _0} - \left .F\pp {F}{\psi }\right |_{\psi _0}. \tag{13.16}

Define the leading-order poloidal field by

B_{\theta 0}(r) \equiv \frac {1}{R_0}\dd {\psi _0}{r}. \tag{13.17}

At the same order the toroidal field is approximately cylindrical, so

F(\psi _0)=R B_\phi \simeq R_0 B_z(r). \tag{13.18}

Using

\frac {d}{d\psi _0} = \frac {dr}{d\psi _0}\dd {}{r} = \frac {1}{R_0 B_{\theta 0}}\dd {}{r}, \tag{13.19}

Eq. (13.16) becomes

\begin{aligned}\frac {R_0}{r}\dd {}{r}\!\left (rB_{\theta 0}\right ) &= -\muo R_0^2 \dd {p}{r}\frac {dr}{d\psi _0} -\frac {1}{2}\dd {F^2}{r}\frac {dr}{d\psi _0} \nonumber \\ &= -\muo R_0 \frac {1}{B_{\theta 0}}\dd {p}{r} -\frac {R_0}{2B_{\theta 0}}\dd {B_z^2}{r}.\end{aligned}

Multiplying by \(B_{\theta 0}/R_0\) gives

\frac {B_{\theta 0}}{r}\dd {}{r}\!\left (rB_{\theta 0}\right ) = -\muo \dd {p}{r} - \dd {}{r}\!\left (\frac {B_z^2}{2}\right ). \tag{13.21}

Expand the derivative on the left:

\frac {B_{\theta 0}}{r}\dd {}{r}(rB_{\theta 0}) = \frac {B_{\theta 0}^2}{r} + B_{\theta 0}\dd {B_{\theta 0}}{r} = \frac {B_{\theta 0}^2}{r} + \dd {}{r}\!\left (\frac {B_{\theta 0}^2}{2}\right ).

Substitute this into Eq. (13.21) and collect total derivatives:

\boxed { \dd {}{r}\!\left ( p+\frac {B_{\theta 0}^2+B_z^2}{2\muo } \right ) = -\frac {B_{\theta 0}^2}{\muo r}. } \tag{13.23}

This is the cylindrical screw-pinch force balance; toroidicity has not yet appeared.

A useful integrated identity. Equation (13.23) becomes especially instructive after one integration. Define

Q(r) \equiv p(r)+\frac {B_{\theta 0}^2(r)+B_z^2(r)}{2\muo }.

Then Eq. (13.23) is \(dQ/dr=-B_{\theta 0}^2/(\muo r)\). Multiply by \(r^2\) and integrate from \(0\) to \(a\):

\int _0^a r^2\dd {Q}{r}\,dr = -\int _0^a \frac {r B_{\theta 0}^2}{\muo }\,dr. \tag{13.25}

Integrating the left-hand side by parts gives

\begin{aligned}\left [r^2 Q(r)\right ]_0^a - 2\int _0^a rQ(r)\,dr &= -\int _0^a \frac {r B_{\theta 0}^2}{\muo }\,dr, \nonumber \\ a^2\left [p(a)+\frac {B_{\theta 0}^2(a)+B_z^2(a)}{2\muo }\right ] - 2\int _0^a r \left ( p+\frac {B_{\theta 0}^2+B_z^2}{2\muo } \right )\,dr &= -\int _0^a \frac {r B_{\theta 0}^2}{\muo }\,dr.\end{aligned}

The \(B_{\theta 0}^2\) terms cancel between the last two terms, leaving

a^2\left [p(a)+\frac {B_{\theta 0}^2(a)+B_z^2(a)}{2\muo }\right ] - 2\int _0^a r \left ( p+\frac {B_z^2}{2\muo } \right )\,dr = 0. \tag{13.27}

Volume averages, poloidal beta, and diamagnetism. For a circular cylinder the natural volume average is

\langle f\rangle = \frac {2}{a^2}\int _0^a f(r)\,r\,dr. \tag{13.28}

Using Eq. (13.28), Eq. (13.27) becomes

2\muo \left [\langle p\rangle -p(a)\right ] + \langle B_z^2\rangle - B_z^2(a) - B_{\theta 0}^2(a) = 0. \tag{13.29}

If the pressure vanishes at the boundary, \(p(a)=0\), then

\langle p\rangle = \frac {B_{\theta 0}^2(a)}{2\muo } + \frac {B_z^2(a)-\langle B_z^2\rangle }{2\muo }. \tag{13.30}

This motivates the standard definition of poloidal beta,

\beta _p \equiv \frac {\langle p\rangle }{B_{\theta 0}^2(a)/(2\muo )} = \frac {4\muo }{a^2 B_{\theta 0}^2(a)}\int _0^a p(r)\,r\,dr. \tag{13.31}

If \(\langle B_z^2\rangle < B_z^2(a)\), the plasma is diamagnetic and \(\beta _p>1\). If \(\langle B_z^2\rangle > B_z^2(a)\), the plasma is paramagnetic and \(\beta _p<1\).

What does a diamagnetic loop measure? Let the vacuum toroidal field at the plasma boundary be \(B_0\equiv B_z(a)\) and write

B_z(r)=B_0+\Delta B_z(r).

The change in toroidal flux through the plasma cross-section is

\Delta \Phi = 2\pi \int _0^a \Delta B_z(r)\,r\,dr. \tag{13.33}

To leading order in \(\Delta B_z\),

\begin{aligned}B_0^2-\langle B_z^2\rangle &= B_0^2 - \frac {2}{a^2}\int _0^a \left (B_0+\Delta B_z\right )^2 r\,dr \nonumber \\ &\simeq -\frac {4 B_0}{a^2}\int _0^a \Delta B_z(r)\,r\,dr = -\frac {2B_0}{\pi a^2}\Delta \Phi .\end{aligned} \tag{13.34}

Combining Eqs. (13.30) and (13.34) shows that a diamagnetic loop measures the plasma energy content through the pressure-driven reduction in enclosed toroidal flux. This is one of the oldest and most robust equilibrium diagnostics.

13.3 First Toroidal Correction and the Shift Equation

Order \(\mathcal O(\epsilon )\) equation. Collect the first-order terms of Eq. (13.10). A convenient form is

\begin{aligned}\left ( \frac {1}{r}\pp {}{r}r\pp {}{r} + \frac {1}{r^2}\pp {^2}{\theta ^2} \right )\psi _1 - \frac {\cos \theta }{R_0}\dd {\psi _0}{r} &= \frac {d}{d\psi _0}\!\left ( -\muo R_0^2 \pp {p}{\psi } -\frac {1}{2}\dd {F^2}{\psi } \right )_{\psi _0}\psi _1 - 2\muo R_0 r \dd {p}{\psi _0}\cos \theta .\end{aligned} \tag{13.35}

The right-hand side comes from the first-order correction to the source term \(-\muo R^2 \pp {p}{\psi }-F\pp {F}{\psi }\), including the geometric factor \(R^2=R_0^2+2R_0 r\cos \theta +\mathcal O(\epsilon ^2)\).

Shift ansatz. The dominant first-order correction is a radial displacement of otherwise circular flux surfaces. We therefore write

\psi (r,\theta ) = \psi _0\!\left (r-\Delta (r)\cos \theta \right ) + \mathcal O(\epsilon ^2). \tag{13.36}

Expand to first order in \(\Delta \):

\psi (r,\theta ) = \psi _0(r) - \Delta (r)\cos \theta \,\dd {\psi _0}{r} + \mathcal O(\epsilon ^2).

Comparing with Eq. (13.13) gives

\boxed { \psi _1(r,\theta ) = -\Delta (r)\cos \theta \,\dd {\psi _0}{r} = - R_0 B_{\theta 0}(r)\,\Delta (r)\cos \theta . } \tag{13.38}

If the plasma boundary is held fixed at \(r=a\), then

\Delta (a)=0, \qquad \left .\dd {\Delta }{r}\right |_{r=0}=0, \tag{13.39}

where the second condition enforces regularity at the magnetic axis.

Angular operator identity. The left-hand side of Eq. (13.35) simplifies because

\left ( \frac {1}{r}\pp {}{r}r\pp {}{r} + \frac {1}{r^2}\pp {^2}{\theta ^2} \right )\!\left [f(r)\cos \theta \right ] = \cos \theta \left [ \frac {1}{r}\dd {}{r}\!\left (r\dd {f}{r}\right )-\frac {f}{r^2} \right ]. \tag{13.40}

Setting \(f(r)=-R_0 B_{\theta 0}(r)\Delta (r)\) and substituting Eq. (13.38) into Eq. (13.35) yields

\begin{aligned}-\frac {1}{r}\dd {}{r}\!\left [r\dd {}{r}\!\left (R_0 B_{\theta 0}\Delta \right )\right ] + \frac {R_0 B_{\theta 0}\Delta }{r^2} - B_{\theta 0} &= - R_0 B_{\theta 0}\Delta \frac {d}{d\psi _0}\!\left ( -\muo R_0^2 \pp {p}{\psi } -\frac {1}{2}\dd {F^2}{\psi } \right )_{\psi _0} - 2\muo R_0 r \dd {p}{\psi _0}.\end{aligned} \tag{13.41}

The cancellation that makes the problem solvable. The leading-order cylindrical equilibrium, Eq. (13.16), may be written as

\frac {R_0}{r}\dd {}{r}\!\left (rB_{\theta 0}\right ) = -\muo R_0^2 \pp {p}{\psi _0} - \frac {1}{2}\dd {F^2}{\psi _0}. \tag{13.42}

Differentiate Eq. (13.42) with respect to \(\psi _0\):

\frac {d}{d\psi _0}\!\left ( -\muo R_0^2 \pp {p}{\psi } -\frac {1}{2}\dd {F^2}{\psi } \right )_{\psi _0} = \frac {1}{R_0 B_{\theta 0}}\dd {}{r}\!\left [ \frac {R_0}{r}\dd {}{r}\!\left (rB_{\theta 0}\right ) \right ]. \tag{13.43}

Insert Eq. (13.43) into Eq. (13.41). The terms proportional to \(\Delta \) then cancel after the radial derivatives on the left are expanded. What survives is the combination containing only \(d\Delta /dr\):

-\frac {1}{rR_0 B_{\theta 0}} \dd {}{r}\!\left [r(R_0 B_{\theta 0})^2\dd {\Delta }{r}\right ] - B_{\theta 0} = -2\muo \frac {r}{B_{\theta 0}}\dd {p}{r}. \tag{13.44}

Multiply by \(-rB_{\theta 0}/R_0\) and use \((R_0 B_{\theta 0})^2=R_0^2 B_{\theta 0}^2\) to obtain

\boxed { \dd {}{r}\left (r B_{\theta 0}^2 \dd {\Delta }{r}\right ) = \frac {r}{R_0} \left ( 2\muo r\dd {p}{r} - B_{\theta 0}^2 \right ). } \tag{13.45}

This is the large-aspect-ratio Shafranov-shift equation.

Integrating once. Apply the regularity condition at \(r=0\) and integrate Eq. (13.45) from \(0\) to \(r\):

\begin{aligned}r B_{\theta 0}^2(r)\dd {\Delta }{r} &= \frac {1}{R_0} \int _0^r \left ( 2\muo x^2\dd {p}{x} - x B_{\theta 0}^2(x) \right )\,dx.\end{aligned}

Integrate the pressure term by parts:

\begin{aligned}\int _0^r 2\muo x^2 \dd {p}{x}\,dx &= 2\muo \left [x^2 p(x)\right ]_0^r - 4\muo \int _0^r x p(x)\,dx \nonumber \\ &= 2\muo r^2 p(r) - 4\muo \int _0^r x p(x)\,dx.\end{aligned}

Hence

\boxed { \dd {\Delta }{r} = \frac {2\muo r\,p(r)}{R_0 B_{\theta 0}^2(r)} - \frac {1}{R_0 r B_{\theta 0}^2(r)} \int _0^r \left [ 4\muo x p(x)+x B_{\theta 0}^2(x) \right ]dx. } \tag{13.48}

Equation (13.48) is usually the most useful working formula. Once \(p(r)\) and \(B_{\theta 0}(r)\) are specified, one more integration gives the shift \(\Delta (r)\) and in particular the magnetic-axis displacement \(\Delta _s=\Delta (0)\) relative to the fixed boundary.

Sign of the shift. At the plasma edge, where \(p(a)=0\),

\left .\dd {\Delta }{r}\right |_{r=a} = -\frac {1}{R_0 a B_{\theta 0}^2(a)} \int _0^a \left [ 4\muo x p(x)+x B_{\theta 0}^2(x) \right ]dx <0. \tag{13.49}

Thus a fixed boundary at \(\Delta (a)=0\) implies \(\Delta (r)>0\) in the interior: the magnetic axis is shifted outward, toward larger major radius, exactly as expected.

13.4 Edge Asymmetry, Internal Inductance, and the Vacuum Continuation

Internal inductance and edge derivative. Define the internal inductance

\ell _i \equiv \frac {2}{a^2 B_{\theta 0}^2(a)} \int _0^a B_{\theta 0}^2(r)\,r\,dr, \tag{13.50}

together with the poloidal beta of Eq. (13.31). Then Eq. (13.49) becomes

\boxed { \left .\dd {\Delta }{r}\right |_{r=a} = -\frac {a}{R_0} \left ( \beta _p+\frac {\ell _i}{2} \right ). } \tag{13.51}

The combination \(\beta _p+\ell _i/2\) therefore appears before any specific profile has been chosen.

Boundary poloidal-field asymmetry. The first-order flux is \[ \psi (r,\theta )=\psi _0(r)-R_0 B_{\theta 0}(r)\Delta (r)\cos \theta . \] At the boundary, where \(\Delta (a)=0\),

\begin{aligned}B_\theta (a,\theta ) &= \frac {1}{R}\pp {\psi }{r}\bigg |_{r=a} \nonumber \\ &= \frac {1}{R_0\left (1+\frac {a}{R_0}\cos \theta \right )} \left [ R_0 B_{\theta 0}(a) - R_0 B_{\theta 0}(a)\left .\dd {\Delta }{r}\right |_{r=a}\cos \theta \right ] + \mathcal O(\epsilon ^2) \nonumber \\ &= B_{\theta 0}(a) \left [ 1 - \left ( \frac {a}{R_0} + \left .\dd {\Delta }{r}\right |_{r=a} \right )\cos \theta \right ].\end{aligned}

Using Eq. (13.51),

\boxed { B_\theta (a,\theta ) = B_{\theta 0}(a) \left [ 1 + \frac {a}{R_0} \left ( \beta _p+\frac {\ell _i}{2}-1 \right )\cos \theta \right ]. } \tag{13.53}

This boundary asymmetry can be measured by arrays of magnetic probes. It is one of the simplest experimental signatures of the same large-aspect-ratio equilibrium physics.

Continuation into the vacuum region. Outside the plasma, \(p=0\) and the enclosed current is constant, so

B_{\theta 0}(r) = B_{\theta 0}(a)\frac {a}{r}, \qquad r>a. \tag{13.54}

Equation (13.45) becomes

\dd {}{r}\left [ \frac {a^2 B_{\theta 0}^2(a)}{r}\dd {\Delta }{r} \right ] = -\frac {a^2 B_{\theta 0}^2(a)}{R_0 r}. \tag{13.55}

Integrate from \(a\) to \(r\):

\begin{aligned}\frac {a^2 B_{\theta 0}^2(a)}{r}\dd {\Delta }{r} - a B_{\theta 0}^2(a)\left .\dd {\Delta }{r}\right |_{r=a} &= -\frac {a^2 B_{\theta 0}^2(a)}{R_0}\ln \frac {r}{a}.\end{aligned}

Divide through by \(a^2 B_{\theta 0}^2(a)\) and insert Eq. (13.51):

\dd {\Delta }{r} = -\frac {r}{R_0} \left ( \beta _p+\frac {\ell _i}{2} + \ln \frac {r}{a} \right ). \tag{13.57}

Impose \(\Delta (a)=0\) and integrate once more:

\Delta (r) = -\frac {1}{R_0} \left [ \left ( \beta _p+\frac {\ell _i}{2}-\frac 12 \right )\frac {r^2-a^2}{2} + \frac {r^2}{2}\ln \frac {r}{a} \right ]. \tag{13.58}

Because \(\Delta (r)\) is negative in the vacuum continuation, the first-order flux perturbation \[ \psi _1=-R_0 B_{\theta 0}\Delta \cos \theta \] is positive there.

Vacuum form of the first-order flux. Using \(B_{\theta 0}(r)=\muo I_p/(2\pi r)\) in the vacuum region, Eq. (13.58) gives

\begin{aligned}\psi _1(r,\theta ) &= - R_0 B_{\theta 0}(r)\Delta (r)\cos \theta \nonumber \\ &= \frac {\muo I_p}{4\pi } \left [ \left ( \beta _p+\frac {\ell _i}{2}-\frac 12 \right ) \left ( 1-\frac {a^2}{r^2} \right ) + \ln \frac {r}{a} \right ]r\cos \theta .\end{aligned} \tag{13.59}

This is the far-vacuum field generated by the plasma itself, still written in the local large-aspect-ratio variables.

13.5 Required Vertical Field

Why a plasma current is not enough. The outward curvature force of the toroidal current channel cannot be balanced by plasma currents alone. The plasma-generated vacuum field must be supplemented by an externally applied vertical field, traditionally supplied by poloidal-field coils. In the large-aspect-ratio ordering, the corresponding flux function is

\psi _{\rm ext}(r,\theta )=R_0 B_v r\cos \theta . \tag{13.60}

The full first-order flux in the vacuum is therefore

\psi _1^{\rm total}=\psi _1^{\rm plasma}+\psi _{\rm ext}.

Toroidal current filament. To determine the required \(B_v\) it is convenient to compare Eq. (13.59) with the exact vacuum field of a toroidal current filament at major radius \(R_0\). The poloidal flux of a ring current \(I_p\) is

\psi _{\rm ring}(R,Z) = \frac {\muo I_p}{\pi }\sqrt {R R_0}\, \frac {(2-k^2)K(k)-2E(k)}{2k}, \tag{13.62}

where \(K\) and \(E\) are complete elliptic integrals and

k^2 = \frac {4R R_0}{(R+R_0)^2+Z^2}. \tag{13.63}

Set \[ R=R_0+r\cos \theta , \qquad Z=r\sin \theta , \qquad \epsilon =\frac {r}{R_0}\ll 1. \] Then

1-k^2 = \frac {r^2}{4R_0^2+4R_0 r\cos \theta +r^2} = \frac {r^2}{4R_0^2} \left [ 1-\frac {r}{R_0}\cos \theta +\mathcal O(\epsilon ^2) \right ].

Define the complementary modulus \(k'=\sqrt {1-k^2}\). To first order,

k' = \frac {r}{2R_0} \left ( 1-\frac {r}{2R_0}\cos \theta \right ) +\mathcal O(\epsilon ^3). \tag{13.65}

We also need

\sqrt {R R_0} = R_0\sqrt {1+\epsilon \cos \theta } = R_0\left (1+\frac {\epsilon }{2}\cos \theta \right )+\mathcal O(\epsilon ^2),

and the standard asymptotic forms

K(k) = \ln \frac {4}{k'} + \mathcal O(k'^2\ln k'), \qquad E(k)=1+\mathcal O(k'^2\ln k'). \tag{13.67}

Since \(k'^2=\mathcal O(\epsilon ^2)\), Eq. (13.67) gives, to first order,

\begin{aligned}K(k) &= \ln \frac {8R_0}{r} + \frac {r}{2R_0}\cos \theta + \mathcal O(\epsilon ^2), \\ E(k) &= 1+\mathcal O(\epsilon ^2).\end{aligned}

Substituting these expansions into Eq. (13.62) yields

\boxed { \psi _{\rm ring}(r,\theta ) = -\frac {\muo I_p R_0}{2\pi } \left [ \ln \frac {8R_0}{r}-2 + \frac {r}{2R_0} \left ( \ln \frac {8R_0}{r}-1 \right )\cos \theta \right ]. } \tag{13.70}

The leading term immediately gives the expected cylindrical poloidal field,

\frac {1}{R_0}\dd {\psi _{\rm ring}^{(0)}}{r} = \frac {\muo I_p}{2\pi r} = B_{\theta 0}(r).

Matching the first-order terms. The \(\cos \theta \) part of Eq. (13.70) is

\psi _{\rm ring}^{(1)} = -\frac {\muo I_p}{4\pi } \left ( \ln \frac {8R_0}{r}-1 \right ) r\cos \theta . \tag{13.72}

The total first-order vacuum field must satisfy

\psi _{\rm ring}^{(1)} + R_0 B_v r\cos \theta = \psi _1(r,\theta ),

with \(\psi _1(r,\theta )\) given by Eq. (13.59). Divide by the common factor \(r\cos \theta \):

\begin{aligned}-\frac {\muo I_p}{4\pi } \left ( \ln \frac {8R_0}{r}-1 \right ) + R_0 B_v &= \frac {\muo I_p}{4\pi } \left [ \beta _p+\frac {\ell _i}{2}-\frac 12+\ln \frac {r}{a} \right ].\end{aligned}

The \(r\) dependence cancels, leaving

\boxed { B_v = \frac {\muo I_p}{4\pi R_0} \left [ \beta _p + \frac {\ell _i-3}{2} + \ln \frac {8R_0}{a} \right ]. } \tag{13.75}

This is Shafranov’s famous vertical-field formula. It states that toroidal equilibrium requires an externally applied field whose magnitude depends on both the internal pressure/current distribution and the global toroidal geometry through \(\ln (8R_0/a)\).

Interactive Shafranov Profile Explorer

Open a browser companion to Problem 13.2. The app uses the lecture’s circular pressure and current profiles, solves for the peaking exponent from \(q(0)\) and \(q(a)\), and reports the resulting \(q(r)\), \(\ell_i\), \(\beta_p\), \(\Delta(r)\), required vertical field, and boundary asymmetry factor.

Open the profile explorer

Experimental perspective. Equation (13.75) is one reason the Shafranov-shift problem has remained so central. The combination \(\beta _p+\ell _i/2\) can be inferred from magnetic measurements and diamagnetic flux, and the required vertical field is set by the same force balance that governs plasma position control. In modern tokamaks this physics is buried inside free-boundary equilibrium solvers, feedback controllers, and reconstruction codes, but the large-aspect-ratio derivation remains the cleanest route to the intuition. It also highlights the limits of external magnetic information: global quantities such as \(\beta _p\), \(\ell _i\), and overall radial force balance are accessible magnetically, whereas the detailed internal \(q\) profile generally requires additional internal constraints, as discussed in the later equilibrium-reconstruction lecture.

Takeaways

Takeaway. The Shafranov shift is not an optional geometric correction. It is the leading signature of toroidicity in a tokamak equilibrium. Starting from the Grad–Shafranov equation, one finds
1.
a zeroth-order cylindrical force balance,
2.
a first-order shift equation, Eq. (13.45),
3.
the appearance of the profile measures \(\beta _p\) and \(\ell _i\), and
4.
the necessity of an external vertical field, Eq. (13.75).

That is why this calculation remains a classic problem: it turns the geometry of a torus directly into measurable equilibrium physics.

Bibliography

Hannes Alfvén. Cosmical Electrodynamics. Clarendon Press, Oxford, 1950.

Hannes Alfvén. Existence of electromagnetic-hydrodynamic waves. Nature, 150:405–406, 1942. doi:10.1038/150405d0.

W. A. Newcomb. Convective instability induced by gravity in a plasma with a frozen-in magnetic field. Physics of Fluids, 4:391–396, 1961. doi:10.1063/1.1706342.

W. A. Newcomb. Lagrangian and hamiltonian methods in magnetohydrodynamics. Nuclear Fusion Supplement, Part 2, pages 451–463, 1962.

H. K. Moffatt. Magnetic Field Generation in Electrically Conducting Fluids. Cambridge University Press, Cambridge, England; New York, NY, USA, 1978. ISBN 9780521216401. Cambridge Monographs on Mechanics and Applied Mathematics.

Problems

Problem 13.1. Vertical field and hoop force

Starting from Eq. (13.75), show that Shafranov’s formula for the vertical field may be written in the intuitive form

2\pi R_0 I_p B_v = \frac {I_p^2}{2}\dd {}{R_0}(L_e+L_i) + 4\pi ^2\int _0^a \left ( p-\frac {B_0\,\delta B_\phi }{\muo } \right )r\,dr,

where \(\delta B_\phi (r)=B_\phi (r)-B_0\), \(L_i=2\pi R_0 \ell _i/2\), and the external inductance of a thin circular current channel is \[ L_e=\muo R_0\left [\ln \frac {8R_0}{a}-2\right ]. \] Then evaluate the net outward hoop force by integrating the radial pressure force around a torus,

F_R = -2\pi ^2 \int _0^a \int _0^{2\pi } \dd {p}{r}\cos ^2\theta \, r^2\,dr\,d\theta ,

and explain physically how each term in the “intuitive” form above contributes to toroidal force balance.

Problem 13.2. Shift of the magnetic axis

Use Eq. (13.48) to compute the magnetic-axis displacement for explicit profiles. Let

\begin{aligned}p(r) &= \hat p\left (1-\frac {r^2}{a^2}\right ), \\ j_\phi (r) &= \hat j\left (1-\frac {r^2}{a^2}\right )^\nu .\end{aligned}

Assume \[ \beta _p=\frac 12, \qquad q(0)\approx 1, \qquad q(a)=3. \]

(a): Determine \(\hat p\), \(\hat j\), and \(\nu \).
(b): Use Ampère’s law to obtain \(B_{\theta 0}(r)\).
(c): Integrate Eq. (13.48) inward from the boundary condition \(\Delta (a)=0\) to determine the Shafranov shift of the magnetic axis, \(\Delta _s/a\).