2. Heuristic model

Insight into the mechanism and time scales of the vertical instability can be gained with the help of a simple heuristic model, involving currents flowing in three parallel rectilinear wires. This toy model is well known, and a simple derivation of the relevant dispersion relation for the case of an ideal wall can be found, e.g. in Yolbarsop et al. (Reference Yolbarsop, Porcelli and Fitzpatrick2021, Reference Yolbarsop, Porcelli, Liu and Fitzpatrick2022) and Pfefferlé & Bhattacharjee (Reference Pfefferlé and Bhattacharjee2018). This calculation is repeated here for the reader's convenience, together with its extension to the case of a resistive wall.

With reference to figure 1, let $y$ be the vertical direction. The currents flow along the $z$ direction, which mimics the toroidal direction of a tokamak plasma. The two ‘external’ currents, $I_{{\rm Ext}}$, are equal and positive and are fixed at $y = \pm l$, while the ‘plasma’ current, $I_p$, can drift along the vertical direction. A vacuum surrounds the three wires. Clearly, $y=0$ is an unstable equilibrium point for the plasma wire. Also shown in figure 1 are magnetic X-points located at $y=\pm l_y < l$.

Figure 1. Schematic diagram of the heuristic model for the vertical instability.

The equation of motion for the plasma wire is

(2.1)\begin{equation} \mu_m\,\ddot { y } = \frac{4\,{ I }_{ P }\,{ I }_{ {\rm Ext} }}{c^{2}} \, \frac{y}{{ l }^{ 2 }-{ y }^{ 2 }} , \end{equation}

with $\mu _m$ the linear mass density, $c$ the speed of light and an over-dot denotes time derivative (c.g.s. units have been adopted). Neglecting self- and mutual induction currents, $I_P$ and $I_{ {\rm Ext} }$ remain constant as the plasma wire is displaced. For small $y\ll l$, the solution of (2.1) is $y =y_0\,{\rm e}^{\gamma _H\,t}$, where $y_0$ is an initial displacement, and ${ \gamma }_{ H } = ( 1/l ) ( 4I_P I_{ {\rm Ext} } / \mu _m c^{2} )^{ 1/2 }$.

If, instead of a plasma wire, we consider a vertically elongated plasma with uniform current density up to an elliptical magnetic surface with minor semi-axis $a$ and major semi-axis $b$, then, according to the analysis of Porcelli & Yolbarsop (Reference Porcelli and Yolbarsop2019), a relationship is found involving the currents $I_P$ and $I_{{\rm Ext}}$ and the distance $l$: ${ I }_{ {\rm Ext} }/{ I }_{ P } = [(b-a)/(b+a)] [l^{2}/(a^{2} + b^{2})]$. With this expression, ${ \gamma }_{ H }$ turns out to depend only on the plasma current $I_P$ and on the semi-axes $a$ and $b$, but not on $I_{ {\rm Ext} }$ and $l$. Also, $\mu _m$ can be replaced by $\mu _m \rightarrow {\rm \pi}\, a\,b\,{ \varrho }_{ m }$, where ${ \varrho }_{ m }$ is the volume mass density. After straightforward algebra, and taking the limit of small ellipticity, $e_0 \ll 1$, where

(2.2)\begin{equation} e_0 = \frac{b^{2} - a^{2}}{b^{2} + a^{2}}, \end{equation}

the growth rate can be written as

(2.3)\begin{equation} { \gamma }_{ H } = { e_0 }^{ 1/2 } {\tau_A}^{{-}1 }, \end{equation}

where ${\tau _A}^{-1 } = {{ B_P }^{\prime }} / { ( 4{\rm \pi} \,{ \varrho }_{ m } ) }^{ 1/2 }$ is the inverse Alfvén time, and ${ B_P }^{'}$ is the radial derivative of the poloidal magnetic field on the magnetic axis. Note that, for a circular plasma cross-section, $e_0 = 0$ and the growth rate $\gamma _H$ vanish. Considering typical values of present-day tokamak experiments, $\gamma _H^{-1}$ is indeed a swift growth time of the order of a few microseconds.

A perfectly conducting wall close to the plasma can provide passive feedback stabilization of the vertical instability. The stabilization mechanism is explained as follows. When the plasma is displaced from its equilibrium position, image currents are induced at the wall. The sign of these currents is such that the corresponding net force opposes the motion of the plasma wire. In the heuristic model, we can model this effect by assuming two currents of opposite sign, $\pm \delta I$, proportional to the displacement of the plasma column, and localized at $y=\pm l$. Then,

(2.4)\begin{equation} L \delta \dot{I } = L D{ I }_{ {\rm Ext} }\dot { y }/l , \end{equation}

where $L$ is an effective inductance and $D$ is a dimensionless proportionality constant, which, in the case of a tokamak plasma, can be determined in terms of the wall geometry (see § 4 here below). After straightforward algebra, and taking the limit $y\ll l$, the equation of motion modified by the feedback currents becomes

(2.5)\begin{equation} \mu_m\, { c }^{ 2 } \,\ddot { y } \approx 4{ I }_{ P }{ I }_{ {\rm Ext}}\,(1-D) { y }/{ { l }^{ 2 } }. \end{equation}

Thus, the vertical instability is suppressed when $D > 1$, which, in the actual tokamak case, translates into a criterion related to the proximity of the wall to the plasma column. As we shall see in § 4, for the special case where the plasma boundary and the wall are modelled by confocal ellipses, the marginal stability criterion corresponds to the magnetic X-points lying on the wall, while feedback stabilization, $D>1$, requires the X-points to lie beyond the wall. In the latter case, (2.5) implies that the plasma oscillates vertically with frequency

(2.6)\begin{equation} \omega ={\pm}{ \omega }_{ H }={\pm}\sqrt { D-1 } \,{ \gamma }_{ H }. \end{equation}

These two solutions with frequency $\omega = \pm \omega _H$ are a first indication of what we dub the VDOM.

For the more realistic case of a resistive wall, (2.4) is modified into

(2.7)\begin{equation} R\delta I + L \delta \dot{I } = L D{ I }_{ {\rm Ext} }\dot { y }/l . \end{equation}

After straightforward algebra, and again in the limit $y \ll l$, we find

(2.8)\begin{equation} \mu_m c^{2}\dddot{y} = 4 I_P I_{ {\rm Ext} }\dot{y}/ l^{2} - 4 I_P \dot{\delta I}/l . \end{equation}

Eliminating $\dot {\delta I}$, we obtain the third-order differential equation for $y$

(2.9)\begin{equation} \dddot{y} + \frac{1}{\tau_R}\ddot{y} + \omega_H^{2}\dot{y}-\frac{1}{\tau_R}\frac{\omega_H^{2}}{D-1}y = 0, \end{equation}

where $\tau _R = L/R$ is the resistive-wall penetration time.

We restrict the discussion that follows to the limit $D>1$. Thus, searching for solutions of the type $y(t)\sim y_0{\rm e}^{\gamma t}$, a cubic dispersion relation for complex $\gamma$ is obtained

(2.10)\begin{equation} \gamma^{3} + \frac{1}{\tau_R}\gamma^{2} + \omega_H^{2}\gamma-\frac{1}{\tau_R}\frac{\omega_H^{2}}{D-1} = 0, \end{equation}

where $\omega _H$ is defined in (2.6). In the limit of small wall resistivity, $\omega _H\tau _R \gg 1$, the two VDOM oscillatory roots of the ideal wall case turn out to be damped; setting $\gamma = -{\rm i}\omega$, the two damped modes have a complex frequency

(2.11)\begin{equation} \omega ={\pm} \omega_H -\frac{{\rm i}D}{2(D-1)\tau_R}. \end{equation}

On the other hand, the third root corresponds to a purely growing resistive instability with a growth rate

(2.12)\begin{equation} \gamma = \frac{1}{(D-1)\tau_R}\equiv\gamma_R . \end{equation}

For this non-rotating, unstable mode, the inverse growth rate scales with $\tau _R$, which, for typical tokamak parameters, gives rise to a characteristic growth time of the order of a few milliseconds. It is then possible to counter the growth of the resistive-wall vertical instability by means of active feedback control associated with currents flowing outside the tokamak vacuum chamber (Lazarus et al. Reference Lazarus, Lister and Neilson1990; Lister et al. Reference Lister, Lazarus, Kellman, Moret, Ferron, Helton, Lao, Leuer, Strait and Taylor1990; Albanese et al. Reference Albanese, Mattei and Villone2004; Portone Reference Portone2017). We point out that, as long as the VDOM are damped, and as long as their oscillation frequency is well separated from the bandwidth of the vertical stabilization closed loop, they are not expected to have an impact on the active feedback stabilization system applied to the non-rotating, vertical resistive-wall mode.

In § 4, a more general dispersion relation, valid for arbitrary values of $e_0$, is obtained based on a complete solution of the ideal-MHD normal-mode problem applied to $n=0$ perturbations. In the limit of small ellipticity, the general dispersion relation reduces to the cubic equation (2.10), confirming in this way the validity of the heuristic model.

3. Vertical rigid-shift mode structure

The analysis of this section follows that already published by us in Porcelli et al. (Reference Porcelli, Yolbarsop, Barberis and Fitzpatrick2021) and Yolbarsop et al. (Reference Yolbarsop, Porcelli and Fitzpatrick2021, Reference Yolbarsop, Porcelli, Liu and Fitzpatrick2022), which we briefly summarize here for the reader's convenience. An adequate model for the normal-mode analysis of the vertical instability is the well-known reduced ideal-MHD model (Strauss Reference Strauss1976). For axisymmetric modes, toroidal effects do not play an important role and, therefore, the straight-tokamak approximation is adopted. The magnetic field is ${\boldsymbol {B} }={\boldsymbol {e} }_{ z }\times \boldsymbol {\nabla } \psi +{ B }_{ z }\,{\boldsymbol {e} }_{ z }$, where ${\boldsymbol {e}}_{ z }$ is the unit vector along the ignorable $z$-direction, which mimics the toroidal coordinate, and $B_z$ is constant. The plasma flow is $\boldsymbol {v} = {\boldsymbol {e }}_{z} \times \boldsymbol {\nabla } \varphi +v_{ z }\,{\boldsymbol {e} }_{ z }$. In the standard low-$\beta$ limit for a tokamak plasma, the fields $B_z$ and $v_z$ decouple from the fields $\psi$ and $\varphi$. Therefore, the magnetic flux function, $\psi$, and the streamfunction, $\varphi$, obey the model equations (Strauss Reference Strauss1976)

(3.1)\begin{gather} \frac { \partial \psi }{ \partial t } +\left[ \varphi ,\psi \right] =0, \end{gather}

(3.2)\begin{gather}\frac { \partial }{ \partial t } \,\boldsymbol{\nabla}\,\boldsymbol{\cdot}\, \left( \varrho \,\boldsymbol{\nabla} \varphi \right) + [\varrho,(\boldsymbol{\nabla}\varphi)^{2}] + [\varphi,\varrho]U +\left[ \varphi ,U \right] =\left[ \psi ,J \right]. \end{gather}

In these equations, all quantities are dimensionless, brackets are defined as $[ \chi,\eta ] ={\boldsymbol {e} }_{ z }\,\boldsymbol {\cdot }\,\boldsymbol {\nabla } \chi \times \boldsymbol {\nabla } \eta$, $J={\nabla }^{2}\psi$ is the normalized current density and $U={ \nabla }^{2}\varphi$ is the normalized flow vorticity. Space and time are normalized as $\hat {r } =r /{ r }_{ 0 }$, where ${ r }_{ 0 }=a\,b/[ ( { a }^{ 2 }+{ b }^{ 2 } ) /2 ]^{ 1/ 2 }$ is a convenient equilibrium scale length, and $\hat {t} = t / { \tau }_{ A }$, where ${ \tau }_{A}$ is the poloidal Alfvén time defined below (2.3). The dimensionless fields are normalized as $\hat { \psi } =\psi /( { { B }_{ p } }^{\prime } { r }_{ 0 }^{ 2 } )$, $\hat {\varphi } = ( \tau _A/{ r_{ 0 }}^{2}) \varphi$; the plasma density is normalized to its on-axis value, $\hat {\varrho } =\varrho _m /{ \varrho }_{ m0 }$, and the current density is $\hat {J } = (4{\rm \pi} /c{ { B }_{ p } }^{\prime }){ J }_{ z }$. In order to simplify the notation, over-hats are actually dropped in (3.1) and (3.2), and in the following.

At equilibrium, fields are stationary, and, by assumption, the equilibrium plasma velocity is negligibly small. The current density, $J_{ {\rm eq} }$, is assumed to be uniform up to an elliptical boundary with minor semi-axis $a$ and major semi-axis $b$, and to drop to zero beyond that boundary, where vacuum is assumed. Clearly, the elliptical boundary must be a magnetic flux surface, which lies necessarily within the region bounded by the magnetic separatrix. In elliptical coordinates $( \mu,\theta )$, where $x= A \sinh ( \mu ) \cos ( \theta )$ and $y = A \cosh ( \mu ) \sin ( \theta )$, with $A = \sqrt { { b }^{ 2 } - { a }^{ 2 } }$, the elliptical boundary corresponds to $\mu = { \mu }_{ b }$, such that $a = A\sinh {\mu _b}$ and $b = A\cosh {\mu _b}$. The equilibrium current density is ${ J }_{ {\rm eq} }( \mu ) = 2 H( { \mu }_{ b } - \mu )$, where $H(x)$ is the Heaviside unit step function.

According to the analysis of Porcelli & Yolbarsop (Reference Porcelli and Yolbarsop2019), in the plasma region inside the elliptical boundary, where $\mu < { \mu }_{ b }$ and $\psi = { \psi }_{ {\rm eq} }^{ - }$, the solution of ${ \nabla }^{ 2 } { \psi }_{ {\rm eq} }^{ - }=2$ that is well behaved on the magnetic axis and that reduces to a constant on the elliptical boundary is best written in terms of Cartesian components

(3.3)\begin{equation} { \psi }_{ {\rm eq} }^{ - }( x , y ) = \frac{1}{2} \left( \frac{ { x }^{ 2 } }{ { b }^{ 2 }} + \frac{ { y }^{ 2 } }{ { a }^{ 2 } } \right) . \end{equation}

In the vacuum region outside the elliptical boundary, where $\mu > { \mu }_{ b }$, the equilibrium flux $\psi _{ {\rm eq} } = { \psi }_{ {\rm eq} }^{ + }$ satisfies ${ \nabla }^{ 2 } { \psi }_{ m }^{ + }=0$. Here, and in the following, the superscripts ‘$-$’ and ‘$+$’ indicate the plasma and vacuum regions, respectively. As we assume no equilibrium current sheets, $\psi _{{\rm eq}}$ and its derivative along the normal to the boundary must be continuous across the boundary. The relevant analytic solution is

(3.4)\begin{equation} { \psi }_{ {\rm eq} }^{ + }( \mu , \theta ) = \frac{1}{2} + { \alpha }^{ 2 } \left\{ \mu - { \mu }_{ b } + \frac{ { e }_{0} }{2} \sinh \left[ 2( \mu - { \mu }_{ b })\right] \cos ( 2 \theta)\right\} \end{equation}

with $\alpha ^{2} = ab/r_0^{2}$ and ${ e }_{ 0 }$ the ellipticity parameter defined in (2.2). Magnetic flux surfaces $\psi _{ eq }(\mu, \theta ) = {\rm {const}}.$ exhibit a magnetic separatrix at $\psi _{ {\rm eq} }(\mu, \theta ) = \psi _X = { \mu }_{ b } \, { \alpha }^{ 2 }$, with X-points located at $\mu = \mu _X =2\mu _b$ and $\theta = \theta _X = {\rm \pi}/2 \pm n{\rm \pi}$.

The equilibrium plasma density profile is also assumed to be uniform up to the elliptical boundary $\mu = \mu _b$, i.e. $\varrho _{ {\rm eq} } = H( { \mu }_{ b } - \mu )$. Thus, this model equilibrium, where the plasma terminates well before the magnetic separatrix, and a wall is present, resembles more closely a limiter tokamak scenario, with the equilibrium magnetic X-points that may lie either outside the plasma containment chamber, or in the vacuum region inside the chamber.

Two main reasons justify the choice of this particular equilibrium. First, we intend to carry out a fully analytic treatment of the normal-mode problem, which can be done if the starting point is a relatively simple equilibrium. The second is that, based on previous studies of the vertical stability problem, effects associated with gradients of the equilibrium plasma current density are not expected to play an important role. Furthermore, the equilibrium adopted in this article is the same as that used in Laval et al. (Reference Laval, Pellat and Soule1974). On the other hand, equilibrium flows are neglected in this article and ought to be considered in future studies.

For stability considerations, set $\psi ( \mu,\theta, t) = \psi _{{\rm eq}}( \mu, \theta ) + \tilde {\psi }( \mu, \theta )\,{\rm e}^{\gamma \,t}$ and $\varphi ( \mu,\theta, t ) = \tilde {\phi }( \mu, \theta )\,{\rm e}^{\gamma t}$, where the over-tilde denotes small perturbed quantities and $\gamma = -{\rm i}\omega$. The linearized versions of (3.1) and (3.2) are

(3.5)\begin{gather} \gamma \,\tilde{\psi}+\left[ \tilde {\varphi },\psi_{ {\rm eq} } \right] =0, \end{gather}

(3.6)\begin{gather}\gamma \,\boldsymbol{\nabla}\,\boldsymbol{\cdot}\, ( \varrho_{ {\rm eq} } \boldsymbol{\nabla} \tilde{ \varphi } ) =\left[ \tilde{\psi}, J_{ {\rm eq} } \right] + \left[ \psi_{ {\rm eq} } , \tilde{ J } \right]. \end{gather}

In the region $\mu < { \mu }_{ b }$, the streamfunction corresponding to a rigid vertical shift is represented in elliptical coordinates by

(3.7)\begin{equation} \tilde{ \varphi }( \mu , \theta ) = \gamma \, \xi \, a \frac{\sinh{\mu}}{\sinh{\mu_b}} \cos{\theta}, \end{equation}

where $\xi$ is the vertical displacement of the plasma column. From the flux freezing condition (3.5) we obtain the corresponding perturbed magnetic flux

(3.8)\begin{equation} \tilde{ \psi }^{ - }( \mu , \theta ) ={-} \frac{\xi}{b}\frac{ \cosh{ \mu } }{ \cosh{ { \mu }_{ b } } } \sin{ \theta }. \end{equation}

Since $\tilde { U }$, $\tilde { J }$, $\boldsymbol {\nabla } \varrho _{ {\rm eq} }$ and $\boldsymbol {\nabla } J_{ {\rm eq} }$ all vanish inside the elliptical boundary, (3.6) is trivially satisfied. Note that in elliptical coordinates, ${\nabla }^{ 2 }\chi = h^{ -2 }(\partial ^{2} \chi /\partial \mu ^{2} + \partial ^{2} \chi /\partial \theta ^{2} )$, where $h = 1 / | \boldsymbol {\nabla } \mu | = 1 / | \boldsymbol {\nabla } \theta |$ is a scale factor, with $h^{2} = A^{2}(\cosh {2\mu }+\cos {2\theta })/2$.

When an ideal or resistive wall is present, the rigid-shift solutions (3.7) and (3.8) for the streamfunction and the perturbed flux in the plasma region retain their validity, while, in the vacuum region, the perturbed flux satisfies $\nabla ^{2}\tilde {\psi }^{+} =0$, whose solution can be best represented as

(3.9)\begin{equation} \tilde{\psi}^{+}(\mu , \theta ) ={-}\frac{\xi_\infty}{b}\exp{[-(\mu -\mu_b)]}\sin{\theta} + \frac{\xi_{{\rm ext}}}{b}\frac{\cosh{\mu}}{\cosh{\mu_b}}\sin{\theta}. \end{equation}

In this expression, $\xi _\infty$ is the amplitude of the rigid vertical displacement in the limit where the wall is moved to infinity, and the term proportional to $\xi _{{\rm ext}}$ represents the contribution to the perturbed flux due to the image currents that form on the wall when this is at a finite distance from the plasma boundary. Continuity of flux at the plasma boundary requires that $\xi = \xi _\infty - \xi _{{\rm ext}}$, and so the actual vertical displacement $\xi$ is reduced, as compared with the no-wall case, by the amount $\xi _{{\rm ext}}$. A perturbed current sheet forms at the plasma boundary

(3.10)\begin{equation} \tilde{ J }( \mu ,\theta ) = \tilde{ j }_b ( \theta ) \delta ( \mu - { \mu }_{ b } ) = \frac{ 1 }{ { h } ^ { 2 } } \left. \left( \frac{ \partial \tilde{ \psi }^{ + }}{ \partial \mu } - \frac{ \partial \tilde{ \psi }^{ - }}{ \partial \mu } \right) \right| _{ { \mu }_{ b } } \delta ( \mu - { \mu }_{ b } ), \end{equation}

where $\delta (x)$ is the Dirac delta function. A staightforward calculation yields the elliptical-angle modulation of the current sheet

(3.11)\begin{equation} \tilde{ j }_b ( \theta ) = \frac{ 2 ( a + b ) }{ b^{2} ( { a }^{ 2 } + { b }^{ 2 } ) } \frac{ { \xi }_{ \infty } \sin{ \theta } }{ 1 + { e }_{ 0 } \cos{ 2 \theta } } . \end{equation}

Note that $\tilde { j }_b ( \theta )$ depends only on $\xi _\infty$ and not on $\xi _{{\rm ext}}$, as the current sheet at the plasma boundary does not depend on the wall image currents.

4. Dispersion relation and quadratic forms

In this section, with the help of quadratic forms, we derive a dispersion relation for $n=0$ vertical modes, which depends on geometrical parameters, $a$ and $b$, and on a function, $D_w(\gamma )$, determined by the geometry and the resistivity of the wall, which becomes independent of $\gamma$ in the ideal wall limit.

Let us introduce the auxiliary streamfunction $\tilde {\varphi }^{{\dagger} }$, such that it equals the complex-conjugate streamfunction in the volume occupied by the plasma (up to $\mu = \mu _b + \epsilon$, i.e. including the perturbed current sheet at the plasma boundary), that we denote by the symbol $\varOmega$, while $\tilde {\varphi }^{{\dagger} }$ equals zero in the vacuum region. We multiply the perturbed plasma equation of motion (3.6) by $\tilde {\varphi }^{{\dagger} }/2\gamma ^{*}$ and integrate it over the whole volume extending to infinity. After standard manipulations, we obtain

(4.1)\begin{equation} -\gamma^{2}\tfrac{1}{2}\int_\varOmega {\rm d}^{3}\,x\rho_{{\rm eq}}\boldsymbol{\xi}\,\boldsymbol{\cdot}\,\boldsymbol{\xi}^{*} ={-}\tfrac{1}{2} \int_\varOmega {\rm d}^{3}\,x\boldsymbol{\xi}^{*} \,\boldsymbol{\cdot}\, \left[ (\tilde{\boldsymbol{J}} \times {\boldsymbol{B}}_{{\rm eq}}) + ({\boldsymbol{J}}_{{\rm eq}}\times \tilde{\boldsymbol{B}})\right], \end{equation}

where $\boldsymbol {\xi } ={\boldsymbol {e}}_z\times \boldsymbol {\nabla } \tilde {\varphi }/\gamma$ is the displacement vector. Thus, the dispersion relation can be written as

(4.2)\begin{equation} -\gamma^{2} = \delta W/\delta I, \end{equation}

where

(4.3)\begin{equation} \delta I = \tfrac{1}{2} \int_\varOmega \rho_{{\rm eq}}\boldsymbol{\xi}\,\cdot\,\boldsymbol{\xi}^{*}\, {\rm d}^{3}x, \end{equation}

and

(4.4)\begin{equation} \delta W ={-} \tfrac{1}{2} \int_\varOmega\boldsymbol{\xi}^{*} \,\boldsymbol{\cdot}\, \boldsymbol{F}(\boldsymbol{\xi}) \,{\rm d}^{3}x, \end{equation}

with $\boldsymbol {F}(\boldsymbol {\xi }) = [ (\tilde {\boldsymbol {J}}\times {\boldsymbol {B}}_{{\rm eq}}) + ({\boldsymbol {J}}_{{\rm eq}}\times \tilde {\boldsymbol {B}})]$ the force density operator in the limit where $\beta$ effects are negligibly small. Note that, even though the potential energy integral is extended to the volume occupied by the plasma, wall effects are included in $\delta W$ through the perturbed magnetic flux, which depends also on the contribution due to image currents that flow on the wall as the plasma boundary is displaced. As it will be shown in the following, $\delta W$ is a real quantity, independent of the mode frequency, in the limit of an ideal wall, but it becomes frequency dependent when the wall is resistive.

Straightforward algebra, using ${\rm d}^{3}x = h^{2}\,{\rm d}\theta \,{\rm d}\mu \,{\rm d} z$, leads to

(4.5)\begin{equation} \delta I = \frac{1}{2} \int_\varOmega {\rm d}^{3}\,x \rho_{{\rm eq}}\frac{\boldsymbol{\nabla}\tilde{\varphi}\,\boldsymbol{\cdot}\,\boldsymbol{\nabla} \tilde{\varphi}^{*}}{|\gamma|^{2}} = \frac{\rm \pi}{2}abL_z \xi ^{2} , \end{equation}

where $L_z$ is the length of the straight tokamak. Without loss of generality, we can take $\xi$ to be real and to correspond to the amplitude of the vertical displacement introduced in the previous section. The perturbed energy integral is

(4.6)\begin{equation} \delta W = \frac{1}{2}\int_\varOmega {\rm d}^{3}\,x \left({\boldsymbol{e}}_z\times\frac{\boldsymbol{\nabla}\tilde{\varphi}^{*}}{\gamma^{*}}\right) \,\boldsymbol{\cdot}\,(\tilde{J}\boldsymbol{\nabla}\psi_{{\rm eq}}+J_{{\rm eq}}\boldsymbol{\nabla}\tilde{\psi}) . \end{equation}

The last term can be further manipulated, using the results for the perturbed streamfunction and the perturbed magnetic flux obtained in § 3. $\delta W$ is written as the sum of two terms. The first term is

(4.7)\begin{align} & \frac{1}{2} \int_{\varOmega} {\rm d}^{3}\,x \tilde{J} \boldsymbol{\nabla}\psi_{{\rm eq}}\,\boldsymbol{\cdot}\, \left({\boldsymbol{e}}_z \times \frac{\boldsymbol{\nabla} \varphi^{*}}{\gamma^{*}}\right)\nonumber\\ & \quad = \frac{L_z}{2}\int_0^{2{\rm \pi}}{\rm d}\theta\int_0^{\mu_b+\epsilon} {\rm d}\mu \,h^{2}\tilde{ j }_b ( \theta ) \delta ( \mu - { \mu }_{ b } ) \boldsymbol{\nabla}\psi_{{\rm eq}}\,\boldsymbol{\cdot}\,\left({\boldsymbol{e}}_z\times\frac{\boldsymbol{\nabla} \varphi^{{{\dagger}}}}{\gamma^{*}}\right) = \frac{\rm \pi}{2}\frac{a+b}{b^{3}}L_z \xi_\infty\xi. \end{align}

The second term can be expressed as

(4.8)\begin{align} & \frac{1}{2} \int_\varOmega {\rm d}^{3}\,x J_{{\rm eq}} \boldsymbol{\nabla}\tilde{\psi}\,\boldsymbol{\cdot}\, \left({\boldsymbol{e}}_z \times \frac{\boldsymbol{\nabla} \varphi^{*}}{\gamma^{*}}\right) \nonumber\\ & \quad = \frac{L_z}{2}\int_0^{2{\rm \pi}}{\rm d}\theta\int_0^{\mu_b+\epsilon} {\rm d}\mu \,h^{2} J_{{\rm eq}} \boldsymbol{\nabla}\tilde{\psi}^{-}\,\boldsymbol{\cdot}\, \left({\boldsymbol{e}}_z \times \frac{\boldsymbol{\nabla} \varphi^{*}}{\gamma^{*}}\right) ={-}\frac{\rm \pi}{2}\frac{a^{2}+b^{2}}{ab^{3}}L_z\xi^{2} . \end{align}

It is convenient to introduce the quantity

(4.9)\begin{equation} D_w = \frac{\xi_{{\rm ext}}}{\hat{e}_0\xi_{\infty}}, \end{equation}

where $\hat {e}_0=e_0b/(a+b)$. Then, combining (4.7), (4.8) and (4.9), remembering that $\xi = \xi _{\infty }-\xi _{ext}$, we obtain

(4.10)\begin{equation} \delta W ={-}\frac{\rm \pi}{2}L_z\frac{1-a/b}{ab}\frac{1-D_w}{1-\hat{e}_0 D_w}\xi^{2} . \end{equation}

Thus, (4.2) yields the dispersion relation

(4.11)\begin{equation} \gamma^{2} = \frac{1-a/b}{a^{2}b^{2}}\frac{1-D_w}{1-\hat{e}_0 D_w}. \end{equation}

This dispersion relation is ‘general’, in the sense that it can be applied to the three cases of interest, i.e. the no-wall limit, the ideal wall case and the resistive-wall case, as wall effects are included in (4.11) through the single stability function, $D_w$. In (4.11), parameters $a$ and $b$ are normalized to the scale length $r_0$ defined below (3.2). Reintroducing dimensions, the dispersion relation (4.11) takes the form

(4.12)\begin{equation} (\gamma\tau_A )^{2} =\frac{r_0^{4}}{a^{2}b^{2}} \left(1-\frac{a}{b}\right)\frac{1-D_w}{1-\hat{e}_0D_w}. \end{equation}

Thus, the mode's stability depends on the quantity $D_w$, which is determined once the location of the wall and its nature, whether ideal or resistive, is decided. As shown in §§ 4.1 and 4.2, if the wall is perfectly conducting, $D_w$ is a real quantity, independent of $\gamma$, that reduces to the parameter $D$ introduced in the heuristic model of S 2. In this case, the relevant dispersion relation is quadratic in $\gamma$, and the sign of $\gamma ^{2}$, which depends on the sign of $1-D$, determines the stability of vertical displacements in the ideal-MHD limit. However, if a resistive wall is considered, then $D_w = D_w(\gamma )$ and the dispersion relation becomes cubic in $\gamma$, as discussed in § 4.3.

4.1. Wall at infinity

The first case we investigate is the no-wall limit, or equivalently, the case where the wall is moved to infinity. The analysis of § 3 shows that, in this case, $\xi _{{\rm ext}} = 0$, hence $D_w=0$. The vertical mode is unstable, with $\gamma ^{2}$ a positive real value (for $b>a$), which reduces to zero in the circular limit, $a=b$. The mode growth rate is

(4.13)\begin{equation} {\gamma} = \frac{ r_{0}^{2} }{ a b } \left( 1- \frac{ a }{ b } \right)^{1/2} \tau_{ A }^{ {-}1 }\equiv\gamma_\infty, \end{equation}

and $\gamma ^{2}$ can also be written in terms of the ellipticity parameter,

(4.14)\begin{equation} \gamma^{2} = ( 1 - { e }_{ 0 } )( 1 + { e }_{ 0 } - \sqrt{ 1 - { e }_{ 0 }^{ 2 } }\, ) \tau_A^{ {-}2 }\equiv\gamma_\infty^{2}. \end{equation}

We remark that this result is valid for arbitrary values of $e_0$ in the interval $0 \leqslant e_0 \leqslant 1$. In the limit of small ellipticity, $\gamma _\infty$ reduces to $\gamma _H$ in (2.3) and the mode growth rate agrees perfectly well with the heuristic result of § 2.

4.2. Passive feedback stabilization: ideal wall

Assume that the wall is represented by an elliptical coordinate surface, $\mu = \mu _w$, confocal to the elliptical plasma boundary at $\mu = \mu _b$, with $\mu _w \geqslant \mu _b$. If the wall is ideal, the perturbed flux in (3.9) must vanish at $\mu = \mu _w$, i.e. $\tilde {\psi }^{+}(\mu _w , \theta ) = 0$. It follows that

(4.15)\begin{equation} \frac{\xi_{{\rm ext}}}{\xi_{\infty}} = \frac{\exp({-(\mu_w-\mu_b)})\cosh{\mu_b}}{\cosh{\mu_w}} = \frac{1+\exp{(2\mu_b)}}{1+\exp{(2\mu_w)}}. \end{equation}

Now, let us use

(4.16a,b)\begin{gather} a_{ w} = A \sinh \mu_{ w}, \quad b_{ w} = A \cosh \mu_{ w}, \end{gather}

(4.17a,b)\begin{gather}a = A \sinh \mu_{ b}, \quad b = A \cosh \mu_{ b}, \end{gather}

where $A = \sqrt {b_w^{2} - a_w^{2}} = \sqrt {b^{2} - a^{2}}$, with $b_w$ and $a_w$ the major and minor semi-axes of the elliptical wall, respectively. It follows that

(4.18)\begin{gather} {\rm e}^{2\mu_b} = \frac{b + a}{b - a}, \end{gather}

(4.19)\begin{gather}{\rm e}^{2\mu_w} = \frac{b_w + a_w}{b_w-a_w}. \end{gather}

Using (4.15)–(4.19) and the definition of $D_w$ in (4.9), we obtain

(4.20)\begin{equation} D_w = \frac { b^{2}+a^{2} }{ (b-a)^{2} } \frac { { b }_{ w }-{ a }_{ w } }{ b_w } \equiv D. \end{equation}

Note that, as the wall is moved further away from the plasma boundary, it becomes more and more ‘circular’, due to the assumption of confocality with the elliptical plasma boundary. Therefore, as $b_w/b\to \infty$, $b_w\to a_w$ and $D\to 0$, as it should. We also observe that, apart from the extra term $\hat {e}_0D$ in the denominator of (4.11), the dispersion relation (4.11) with $D_w = D$ agrees very well with the heuristic dispersion relation obtained in § 2. We can conclude that, if $D>1$, $\gamma ^{2}$ is negative, and thus the vertical mode oscillates with a real frequency

(4.21)\begin{equation} \omega ={\pm}\omega_0 ={\pm}\left[\frac{D-1}{1-\hat{e}_0D}\right]^{1/2}\gamma_\infty , \end{equation}

which reduces to $\omega _H$ defined in (2.6) in the limit of small ellipticity. The dispersion relation is quadratic and the only solutions that can be extracted, in the limit $D>1$, are the VDOM solutions. Thus, an ideal elliptical wall, confocal with the elliptical plasma boundary, can provide passive feedback stabilization of the vertical instability when $D_w = D>1$.

The term $\hat {e}_0D$ can be neglected if it is small as compared with unity, but in fact it plays an important role if the wall coincides with the plasma boundary. Indeed, $D$ reaches its maximum value in this limit, where $\mu _w = \mu _b$, $D=D_{\max }=(a+b)/(e_0b) = 1/\hat {e}_0$ and the denominator in (4.21) vanishes. Therefore, in this special limit, the oscillation frequency goes to infinity, but also, the amplitude of the displacement, $\xi = \xi _\infty - \xi _{{\rm ext}}$, goes to zero, as $\xi _{{\rm ext}}\to \xi _\infty$. As a numerical example, take $a_w = a$, $b_w = b$ and $b/a = 1.8$. In this case, $D_{\max }\approx 2.9$. As the wall is placed further away from the plasma boundary, the value of $D$ decreases monotonically, and a purely oscillatory solution is found for as long as $D$ remains larger than unity.

It can be easily checked, using (4.15)–(4.19), that the marginally stable value, $D=1$, is obtained when $\mu _w = 2\mu _b$, which corresponds to the elliptical wall intercepting the X-points. Thus, our model indicates that values of $D<1$, for which no passive feedback stabilization is possible, are found when the X-points lie inside the volume bounded by the wall. We point out that this result is in perfect agreement with that obtained by Laval et al. (Reference Laval, Pellat and Soule1974) on the basis of the ideal-MHD energy principle. However, this result holds true for the special case of confocal wall and plasma boundaries. For different wall shapes that fit more tightly with the plasma boundary, passive feedback stabilization due to the ideal wall does not necessarily require that the X-points lie beyond the wall. Furthermore, in real tokamak plasmas, additional coils and/or short-circuited plates located inside the vacuum chamber can contribute to passive feedback stabilization (Albanese et al. Reference Albanese, Mattei and Villone2004; Portone Reference Portone2017).

Finally, we point out that, as shown in Yolbarsop et al. (Reference Yolbarsop, Porcelli and Fitzpatrick2021, Reference Yolbarsop, Porcelli, Liu and Fitzpatrick2022), if the plasma extends to the magnetic separatrix, as is normally the case for a divertor tokamak plasma, axisymmetric plasma currents associated with vertical displacement perturbations are driven in the vicinity of the magnetic X-points. As a consequence, these currents can also contribute to passive feedback stabilization of the vertical instability.

4.3. Passive feedback stabilization: resistive wall

Lastly, we consider the more realistic case of a resistive wall. In this situation the perturbed magnetic flux is not zero at the wall, but can diffuse across on the resistive-wall time scale. We identify three different regions: inside the elliptical boundary, $\mu \leqslant \mu _b$, where the perturbed flux $\tilde {\psi }$ equals $\tilde {\psi }^{-}$ given by (3.8); the vacuum region between the plasma boundary and the wall, $\mu _b \leqslant \mu \leqslant \mu _w$, where $\tilde {\psi }$ equals $\tilde {\psi }^{+}$ given by (3.9); and the vacuum outside the wall, where $\tilde {\psi }=\tilde {\psi }_{{\rm out}}$ is the solution of $\nabla ^{2}\tilde {\psi }_{{\rm out}} = 0$ that decays to zero at infinity. Focusing on vertical displacements with elliptical mode number $m=1$, the relevant solution is

(4.22)\begin{equation} \tilde{\psi}_{{\rm out}}(\mu,\theta) = \psi_{0}\exp({-(\mu-\mu_w)})\sin{\theta} \quad\text{for } \mu \geqslant \mu_w+\delta_w , \end{equation}

where $\delta _{w}$ is a small parameter representing the average width of the thin wall.

Two conditions at the wall determine the parameter $\psi _{0}$ and $D_w$; the latter will turn out to be a function of complex $\gamma$. The first condition is continuity of flux at the wall, $\tilde {\psi }^{+}(\mu _w,\theta ) = \tilde {\psi }_{{\rm out}}(\mu _w,\theta )$, which gives

(4.23)\begin{equation} \psi_{0} ={-}\frac{\xi_{\infty}}{b}\exp({-(\mu_w-\mu_b)}) + \frac{\xi_{{\rm ext}}}{b}\frac{\cosh{\mu_w}}{\cosh{\mu_b}} . \end{equation}

The second condition involves the current flowing along the wall. Consider the resistive Ohm's law for the perturbed magnetic flux within the wall, $\partial {\psi }_{w}/\partial t = (\eta c^{2}/4{\rm \pi} )\nabla ^{2}{\psi }_{w}$, where $\eta$ is the wall resistivity. After proper normalization, we obtain

(4.24)\begin{equation} \nabla^{2}\tilde{\psi}_w = \frac{\gamma}{\varepsilon_{\eta}}\tilde{\psi}_w, \end{equation}

where $\varepsilon _{\eta } = (\eta \tau _A c^{2})/(4{\rm \pi} r_0^{2})$ is the inverse of the relevant dimensionless Lundquist number. Since the wall is relatively thin, we can neglect the dependence on $\mu$ of the perturbed flux at the right-hand side of (4.24). This approximation is similar to the standard constant-$\psi$ approximation used in magnetic reconnection theory. Set $\tilde {\psi } = \hat {\psi }(\mu )\sin {\theta }$. We assume the ordering $\partial _\mu ^{2}\hat {\psi } \sim \hat {\psi }/\mu _w\delta _w \gg \hat {\psi }$. Then, (4.24) can be approximated as

(4.25)\begin{equation} \frac{{\rm d}^{2}\hat{\psi}}{{\rm d}\mu^{2}} = h^{2}\frac{\gamma}{\varepsilon_{\eta}}\psi_0. \end{equation}

Integrating across the thin wall,

(4.26)\begin{equation} \int_{\mu_w}^{\mu_w+\delta_{w}} \frac{{\rm d}^{2}\hat{\psi}}{{\rm d}\mu^{2}} {\rm d}\mu = \int_{\mu_w}^{\mu_w+\delta_{w}} h^{2}\frac{\gamma}{\varepsilon_{\eta}} \psi_0 \,{\rm d}\mu . \end{equation}

The scale factor $h$ is defined below (3.8). Strictly speaking, $h$ depends on both $\mu$ and $\theta$. However, in the limit of small wall ellipticity, which is assumed in the following, the $\theta$ modulation of $h$ along the wall is small, so that we can approximate $h(\mu _w , \theta ) \approx b_w$. Carrying out the integration, we find

(4.27)\begin{equation} \left.\left( \frac{{\rm d} \hat{\psi}_{{\rm out}}}{{\rm d} \mu} - \frac{{\rm d} \hat{\psi}^{+}}{{\rm d} \mu} \right) \right|_{\mu_w} = \frac{\gamma \delta_w b_w^{2}}{\varepsilon_{\eta}} \psi_0, \end{equation}

or equivalently

(4.28)\begin{equation} {-}\frac{\xi_\infty}{b}\exp({-(\mu_w - \mu_b)}) - \frac{\xi_{ext}}{b}\frac{\sinh{\mu_w}}{\cosh{\mu_b}} = \left( 1 + \frac{\gamma \delta_w b_w^{2}}{\varepsilon_{\eta}} \right)\psi_0. \end{equation}

Equations (4.23) and (4.28) are two equations that can be used to determine the two parameters, $\psi _0\propto \xi _\infty$, and $\xi _{ext}/\xi _\infty \propto D_w$, cf. (4.9). After straightforward algebra, we obtain

(4.29)\begin{equation} D_w(\gamma) = D\frac{\gamma \tau_{\eta}}{1+\gamma \tau_{\eta} }, \end{equation}

and

(4.30)\begin{equation} \psi_0 ={-}\frac{\xi_{\infty}}{b}\frac{b_w}{(a+b)}\frac{e_0D}{(1 + \gamma \tau_{\eta})} . \end{equation}

In these equations, $D$ was defined in (4.15), and $\tau _{\eta } = [b_w^{3}/(a_w+b_w)]\delta _{w} / \varepsilon _{\eta }$ is the resistive-wall time (normalized to the relevant Alfvén time). Notice that the stability parameter $D_w$ depends on $\gamma$ due to finite wall resistivity. In the ideal wall limit, $\tau _\eta \rightarrow \infty$, $D_w(\gamma )\to D$ reduces to the parameter $D$ discussed in § 4.2. Also, (4.30) shows that the perturbed flux at the wall goes to zero in the ideal wall limit, as it should.

Substituting $D_w(\gamma )$ in (4.11), one obtains the dispersion relation for vertical displacements in the case of feedback stabilization via a resistive wall

(4.31)\begin{equation} \gamma^{3} + \gamma^{2}\frac{1}{\tau_{\eta}} \frac{1}{1-\hat{e}_0D} + \gamma \omega_0^{2} + \omega_0^{2}\frac{1}{\tau_{\eta}} \frac{1}{1-D} = 0, \end{equation}

with $\omega _{0}$ defined in (4.21). It is reassuring that (4.31) agrees very well with the dispersion relation found heuristically in § 2, cf. (2.10). This is particularly true in the small ellipticity limit, where the term $\hat {e}_0D$ can be neglected and the frequency $\omega _0$ reduces to $\omega _H$, cf. (2.6).

The three roots of this cubic dispersion relation are easily determined in the realistic limit $\omega _0 \tau _{\eta } \gg 1$. After straightforward algebra, for $D>1$, one obtains two damped oscillatory roots and an unstable growing solution ($\gamma = -{\rm i}\omega$)

(4.32)\begin{gather} \omega ={\pm} \omega_0 - {\rm i} \frac{1}{2\tau_{\eta}}\frac{ D(1-\hat{e}_0)}{(D-1)(1-\hat{e}_0D)} , \end{gather}

(4.33)\begin{gather}\gamma = \frac{1}{(D-1)\tau_{\eta}}. \end{gather}

The two oscillatory, weakly damped roots of (4.32) correspond to the two VDOM that oscillate with a frequency slightly below the poloidal Alfvén frequency, as already found in § 4.2. The effect of a resistive wall is to introduce a small damping rate of the order of the inverse resistive-wall time. The third root, (4.33), corresponds to an unstable mode, growing on the resistive-wall time scale. Typical values of the resistive wall time in tokamak devices are of the order of a few milliseconds. As we have already pointed out in § 2, due to its slow growth, this non-rotating, unstable mode is normally suppressed by means of active feedback control systems (Albanese et al. Reference Albanese, Mattei and Villone2004).