Skip to main content Accessibility help


  • Access
  • Cited by 5



MathJax is a JavaScript display engine for mathematics. For more information see
      • Send article to Kindle

        To send this article to your Kindle, first ensure is added to your Approved Personal Document E-mail List under your Personal Document Settings on the Manage Your Content and Devices page of your Amazon account. Then enter the ‘name’ part of your Kindle email address below. Find out more about sending to your Kindle. Find out more about sending to your Kindle.

        Note you can select to send to either the or variations. ‘’ emails are free but can only be sent to your device when it is connected to wi-fi. ‘’ emails can be delivered even when you are not connected to wi-fi, but note that service fees apply.

        Find out more about the Kindle Personal Document Service.

        Analytic study on low- $n$ external ideal infernal modes in tokamaks with large edge pressure gradients
        Available formats

        Send article to Dropbox

        To send this article to your Dropbox account, please select one or more formats and confirm that you agree to abide by our usage policies. If this is the first time you use this feature, you will be asked to authorise Cambridge Core to connect with your <service> account. Find out more about sending content to Dropbox.

        Analytic study on low- $n$ external ideal infernal modes in tokamaks with large edge pressure gradients
        Available formats

        Send article to Google Drive

        To send this article to your Google Drive account, please select one or more formats and confirm that you agree to abide by our usage policies. If this is the first time you use this feature, you will be asked to authorise Cambridge Core to connect with your <service> account. Find out more about sending content to Google Drive.

        Analytic study on low- $n$ external ideal infernal modes in tokamaks with large edge pressure gradients
        Available formats
Export citation


The problem of pressure driven infernal type perturbations near the plasma edge is addressed analytically for a circular limited tokamak configuration which presents an edge flattened safety factor. The plasma is separated from a metallic wall, either ideally conducting or resistive, by a vacuum region. The dispersion relation for such types of instabilities is derived and discussed for two classes of equilibrium profiles for pressure and mass density.

1 Introduction

Attaining high performance tokamak operation regimes is of primary interest for a future fusion reactor. Such a condition is usually achieved in the so-called tokamak H-mode (high confinement) regime. Besides the high confinement time, the H-mode is characterised by large pressure gradients located near the plasma boundary. These gradients favour the development of a particular kind of instability called the edge localised mode (ELM). Such type of perturbations are associated with rapid and violent bursts of energy and particles which can severely damage the plasma facing components with intolerable heat loads. An increasing interest has been therefore devoted to the study of the so-called quiescent high confinement mode (QH-mode). The QH-mode regime (Burrell et al. 2002; Suttrop et al. 2003, 2005) shares with the H-mode a large edge pressure pedestal (edge transport barrier) and high energy confinement time but without the presence of ELMs. These are replaced by a less dangerous continuous magnetohydrodynamic (MHD) activity called edge harmonic oscillation (EHO). The energy loads deposited by EHOs are much lower compared to the ones associated with ELMs. EHOs have low- $n$ toroidal number (e.g. $n=1,2,3$ while ELMs are instead characterised by large $n$ values) and are always experimentally observed during the QH-mode operation (Burrell et al. 2001, 2002, 2005; Suttrop et al. 2004; Solano et al. 2010).

These low- $n$ oscillations have been linked with kink/peeling modes (Connor et al. 1998) which are proposed as a possible candidate for the explanation of the appearance of these perturbations. Linear and nonlinear simulations of QH-mode DIII-D plasma discharges showed that kink/peeling modes are the main unstable modes in the nonlinear phase whose saturation leads to a three-dimensional stationary state with toroidal periodicity characterised by a low- $n$ number (generally $n=1,2$ ) (Liu et al. 2015). These low- $n$ modes are not the linearly most unstable modes in the linear phase, but they are driven by the nonlinear coupling with medium- $n$ ( ${\sim}3,4,5$ ) harmonics (Liu et al. 2015).

Steep edge pressure gradients are associated in the low collisionality regime with a large edge bootstrap current contribution which produces an edge flattening of the safety factor ( $q$ ) profile. Numerical studies of low- $n$ MHD modes in QH regime with a plateau in $q$ near the edge, corresponding to the peak of the bootstrap current, have been found to have infernal-like features (Zheng, Kotschenreuther & Valanju 2013a , b ). Indeed under these conditions, i.e. large pressure gradients and $q$ flat over and extended region, it is likely that infernal type instabilities develop. These instabilities are characterised by a toroidicity induced coupling between neighbouring Fourier harmonics, i.e. a dominant $m_{0}$ mode is coupled with its sidebands $m_{0}\pm 1$ . Three-dimensional free boundary MHD equilibria simulations of Joint European Torus (JET) and Tokamak à Configuration Variable (TCV)-like plasmas, in which a large edge bootstrap current flattens the safety factor, recovered saturated ideal kink/peeling structures (Cooper et al. 2015, 2016a , b ). The computed equilibrium state is characterised by a distorted boundary with a dominant $n=1$ Fourier component where the corrugation is driven mainly by non-axisymmetric components of the parallel current density.

Thus the main aim of this work is to study analytically the stability properties against infernal modes of MHD equilibria which present a local flattening of the safety factor close to the plasma boundary associated with large pressure gradients. We follow the analysis presented in Brunetti et al. (2018) bringing out the work more deeply, and more importantly, looking at extensions which include toroidal mode number analysis and resistive walls. Also we extend the results previously obtained in Brunetti et al. (2018) to different types of more realistic profiles ( $q$ etc.). We assume a vacuum region separating the plasma and the metallic wall. The analysis of the perturbation is split into various regions, each treated separately. The eigensolutions for the Fourier modes in each region are matched with the appropriate choice of the plasma–vacuum and vacuum–wall jump interface conditions. This eventually yields the dispersion relation.

The paper is organised as follows: § 2 describes the geometry and the physical model employed for the plasma modelling. The infernal modes equations are introduced and a discussion on their validity in application to our problem is given. In § 3 the eigensolutions for coupled Fourier harmonics are solved for a particular choice of the safety factor profile and their logarithmic jumps across the transition point between high and low-shear regions are evaluated. In § 4 we solve for the vacuum perturbation and apply the appropriate boundary conditions at the plasma edge and metallic wall interfaces. The solution of the equation for the main eigenmode $m_{0}$ and the matching procedure with the high-shear region eigensolutions is presented in § 5 where two classes of equilibrium profiles for temperature and mass density are considered. The dispersion relation is eventually derived and discussed in § 6 with two different wall boundary conditions. Finally the findings of this work and future outlook are summarised in § 7. The appendices include a description of inertial corrections and residual coupling effects at the resonant surface of the $m_{0}-1$ harmonic.

2 Physical model

We consider a large aspect ratio tokamak of major and minor radii $R_{0}$ and $a$ respectively ( $\unicode[STIX]{x1D700}=a/R_{0}\ll 1$ ) with shifted circular toroidal surfaces. The vacuum region extends for $a<r<b$ and for $r>b+d$ and a metallic wall of thickness $d$ is located at $b<r<b+d$ . The equilibrium assumes a strong toroidal field ( $B_{T}$ ) and a smaller poloidal field ( $B_{P}$ ), i.e. $B_{P}/B_{T}\sim \text{O}(\unicode[STIX]{x1D700})$ with $\unicode[STIX]{x1D6FD}\sim \text{O}(\unicode[STIX]{x1D700}^{2})$ (here $\unicode[STIX]{x1D6FD}=p_{0}/B_{T}^{2}$ where $p_{0}$ is the equilibrium pressure). We use a right-handed straight field line coordinate system $(r,\unicode[STIX]{x1D717},\unicode[STIX]{x1D711})$ where $r$ is a flux label with the dimensions of length, $\unicode[STIX]{x1D717}$ and $\unicode[STIX]{x1D711}$ are the poloidal-like (counter-clockwise) and toroidal angles respectively. The contravariant and covariant basis vectors are denoted hereafter by $\unicode[STIX]{x1D735}C^{i}$ and $\boldsymbol{e}_{C^{i}}$ respectively, with $C^{i}=(r,\unicode[STIX]{x1D717},\unicode[STIX]{x1D711})$ .

The magnetic field in the plasma is represented in terms of flux functions (D’haeseleer et al. 1991):

(2.1) $$\begin{eqnarray}\displaystyle & \displaystyle \boldsymbol{B}=\unicode[STIX]{x1D735}f\times \unicode[STIX]{x1D735}\unicode[STIX]{x1D717}-\unicode[STIX]{x1D735}\unicode[STIX]{x1D713}\times \unicode[STIX]{x1D735}\unicode[STIX]{x1D711}. & \displaystyle\end{eqnarray}$$

The equilibrium fluxes, denoted with $\unicode[STIX]{x1D713}_{0}$ and $f_{0}$ , depend only on $r$ , and $q=f_{0}^{\prime }/\unicode[STIX]{x1D713}_{0}^{\prime }$ with $f_{0}^{\prime }\approx B_{0}r$ which follows from $B_{T}\approx R_{0}B_{0}/R$  (Wesson 2011) ( $B_{0}$ is the magnetic field strength on the axis and the prime denotes the derivative with respect to the radial variable). Normalising $\unicode[STIX]{x1D707}_{0}=1$ ( $\unicode[STIX]{x1D707}_{0}$ is the vacuum permeability), from Ampére’s law $\boldsymbol{J}=\unicode[STIX]{x1D735}\times \boldsymbol{B}$ where $\boldsymbol{J}$ is the current density, we have $J_{0}=J_{0}^{\unicode[STIX]{x1D711}}/B_{0}^{\unicode[STIX]{x1D711}}=1/R_{0}[(r\unicode[STIX]{x1D707})^{\prime }+\unicode[STIX]{x1D707}]$ and $J_{0}^{\prime }=1/(mR_{0})[(r^{2}k_{||,m})^{\prime }/r]^{\prime }$  (Hazeltine & Meiss 1992), where $\unicode[STIX]{x1D707}=1/q$ and $k_{||,m}=m\unicode[STIX]{x1D707}-n$ with $m$ and $n$ integers.

Our stability analysis is based on the following ideal MHD equations (Hazeltine & Meiss 1992):

(2.2) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D70C}[\unicode[STIX]{x2202}_{t}\boldsymbol{v}+\boldsymbol{v}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\boldsymbol{v}]=-\unicode[STIX]{x1D735}p+\boldsymbol{J}\times \boldsymbol{B}, & \displaystyle\end{eqnarray}$$
(2.3) $$\begin{eqnarray}\displaystyle & \displaystyle \boldsymbol{E}+\boldsymbol{v}\times \boldsymbol{B}=0, & \displaystyle\end{eqnarray}$$
(2.4) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x2202}_{t}p+\boldsymbol{v}\boldsymbol{\cdot }\unicode[STIX]{x1D735}p+\unicode[STIX]{x1D6E4}p\unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{v}=0, & \displaystyle\end{eqnarray}$$
(2.5) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x2202}_{t}\unicode[STIX]{x1D70C}+\unicode[STIX]{x1D735}\boldsymbol{\cdot }(\unicode[STIX]{x1D70C}\boldsymbol{v})=0, & \displaystyle\end{eqnarray}$$

where $\unicode[STIX]{x1D6E4}$ is the adiabatic index, $\boldsymbol{E}$ is the electric field and $\boldsymbol{v}$ is the plasma MHD velocity. It is convenient to write $\boldsymbol{v}=\unicode[STIX]{x2202}\unicode[STIX]{x1D743}/\unicode[STIX]{x2202}t$ where $\unicode[STIX]{x1D743}$ is the Lagrangian perturbed fluid displacement. Finally $\unicode[STIX]{x1D70C}$ is the mass density and $p$ the plasma pressure. We stress the point that we consider an ideal plasma so that the development of tearing type instabilities is prevented. By means of (2.3) we have:

(2.6) $$\begin{eqnarray}\displaystyle & \displaystyle \tilde{\boldsymbol{B}}=\unicode[STIX]{x1D735}\times (\unicode[STIX]{x1D743}\times \boldsymbol{B}_{0}), & \displaystyle\end{eqnarray}$$

where $\tilde{\boldsymbol{B}}$ is the perturbed magnetic field and $\boldsymbol{B}_{0}$ the equilibrium magnetic field. From now on equilibrium quantities will be denoted with the subscript 0 and perturbed quantities have a dependence on time and angular variables of the type $\exp [\text{i}(m\unicode[STIX]{x1D717}-n\unicode[STIX]{x1D711})+\unicode[STIX]{x1D6FE}t]$ .

For our analysis we choose a safety factor which has a local edge flattening as shown in figure 1. The gap between the safety factor and the resonance $m_{0}/n$ is denoted by $\unicode[STIX]{x1D6FF}q$ , i.e. $\unicode[STIX]{x1D6FF}q=q-m_{0}/n$ ( $\unicode[STIX]{x1D6FF}q$ can be either positive or negative) with $\unicode[STIX]{x1D6FF}q$ small. It is implicitly assumed that $m_{0}>1$ . In a more realistic geometry with separatrix the safety factor grows to large values, although in a narrow region much smaller than the expected radial extension of the mode. This particular case is outside of the scope of the present work. Nevertheless, motivated by Medvedev et al. (2006), Zheng et al. (2013a , b ), we expect that the main physical ingredients for the development of an edge infernal type MHD instability are captured. According to  Turnbull & Troyon (1989) we assume that the ratio $q_{a}/q_{0}$ is sufficiently large ( $q_{0}$ and $q_{a}$ are the values of the safety factor on the magnetic axis and at the edge respectively) in order to be stable against low- $n$ external kink modes.

Figure 1. Model safety factor with the vacuum region extending from $a$ to $b$ . In the edge low-shear region ( $r_{\ast }<r<a$ ) the $q$ profile is close to $m_{0}/n$ . The position of the $(m_{0}-1)/n$ resonant surface is highlighted.

According to figure 1, we refer to the region which extends from the magnetic axis to $r=r_{\ast }<a$ as the sheared region where the magnetic shear is large ( $q^{\prime }\sim \text{O}(1)$ ) and the pressure gradient sufficiently small. Here mode coupling is prevented and different Fourier harmonics behave independently according to (Mikhailovskii 1998; Brunetti et al. 2014):

(2.7) $$\begin{eqnarray}\displaystyle & \displaystyle {\displaystyle \frac{\text{d}}{\text{d}r}}\left[r^{3}k_{||,m}^{2}{\displaystyle \frac{\text{d}X_{m}}{\text{d}r}}\right]-r(m^{2}-1)k_{||,m}^{2}X_{m}=0. & \displaystyle\end{eqnarray}$$

It will be shown in § 4 that also the vacuum magnetic perturbation is described by (2.7). Indeed in the vacuum the safety factor is expected to grow as $r^{2}$ so that we can regard also this region as a sheared region.

The low-shear region extends from $r_{\ast }$ to the plasma edge $r=a$ . Here the safety factor is almost flat and close to $m_{0}/n$ . In this region pressure gradients and field line bending weakening allow for mode coupling between neighbouring poloidal Fourier harmonics (Hastie & Hender 1988; Waelbroeck & Hazeltine 1988; Gimblett, Hastie & Hender 1996; Brunetti et al. 2014). The presence of three Fourier modes with mode numbers $(m_{0},n),(m_{0}\pm 1,n)$ is assumed whose relative amplitude is $X_{m_{0}\pm 1,n}\sim \text{O}(\unicode[STIX]{x1D700})X_{m_{0},n}$ , where $X=\unicode[STIX]{x1D743}\boldsymbol{\cdot }\unicode[STIX]{x1D735}r$ . For sake of simplicity hereafter we will drop the subscript $n$ and we introduce the notation $X_{\pm }=X_{m_{0}\pm 1,n}$ . The equations describing the perturbation are (Hastie & Hender 1988; Waelbroeck & Hazeltine 1988; Gimblett et al. 1996; Wahlberg, Graves & Chapman 2013; Brunetti et al. 2014):

(2.8) $$\begin{eqnarray}\displaystyle & \displaystyle {\displaystyle \frac{\text{d}}{\text{d}r}}\left[r^{3}Q{\displaystyle \frac{\text{d}X_{m_{0}}}{\text{d}r}}\right]+r\left[(1-m_{0}^{2})Q+r\left({\displaystyle \frac{A_{2}}{n^{2}}}\right)^{\prime }-{\displaystyle \frac{\unicode[STIX]{x1D6FC}^{2}}{2}}+{\displaystyle \frac{\unicode[STIX]{x1D6FC}r}{R_{0}}}\left({\displaystyle \frac{1}{q^{2}}}-1\right)\right]X_{m_{0}} & \displaystyle \nonumber\\ \displaystyle & \displaystyle +{\displaystyle \frac{\unicode[STIX]{x1D6FC}}{2}}\left[{\displaystyle \frac{r^{-m_{0}}}{1+m_{0}}}(r^{2+m_{0}}X_{+})^{\prime }+{\displaystyle \frac{r^{m_{0}}}{1-m_{0}}}(r^{2-m_{0}}X_{-})^{\prime }\right]=0, & \displaystyle\end{eqnarray}$$
(2.9) $$\begin{eqnarray}\displaystyle & \displaystyle [r^{-1\mp 2m_{0}}(r^{2\pm m_{0}}X_{\pm })^{\prime }]^{\prime }={\displaystyle \frac{1\pm m_{0}}{2}}[\unicode[STIX]{x1D6FC}r^{\mp m_{0}}X_{m_{0}}]^{\prime }, & \displaystyle\end{eqnarray}$$

with $Q=(\unicode[STIX]{x1D6FF}q/q)^{2}+A_{1}/n^{2}$ , $\unicode[STIX]{x1D714}_{A}=B_{0}/(R_{0}\sqrt{\unicode[STIX]{x1D70C}})$ and $\unicode[STIX]{x1D6FC}=-(2R_{0}p_{0}^{\prime }q^{2})/B_{0}^{2}$ . Equilibrium mass density gradients are allowed since they play an important role in the determination of the pressure gradients in the narrow edge region where the perturbation we are interested in is localised. In the limit $\unicode[STIX]{x1D6FE}/[\unicode[STIX]{x1D6E4}p_{0}/(\unicode[STIX]{x1D70C}_{0}R_{0}^{2})]\ll 1$ we approximate (Wahlberg et al. 2013) $A_{1}\simeq (\unicode[STIX]{x1D70C}_{0}\unicode[STIX]{x1D6FE}^{2})/(\bar{\unicode[STIX]{x1D70C}}\unicode[STIX]{x1D714}_{A}^{2})[1+2q^{2}]$ and $A_{2}\approx A_{1}$ ( $\bar{\unicode[STIX]{x1D70C}}$ is the value of the mass density on the magnetic axis).

Equations (2.8) and (2.9) have been derived for $\unicode[STIX]{x1D6FD}\sim \unicode[STIX]{x1D700}^{2}$ . Then the following question arises: are they applicable when the pressure gradient increases, i.e. when $\unicode[STIX]{x1D6FC}\sim 1$ ? We point out that from (2.8), the requirement that the inertial term is of the same order as the coupling term gives $Q\sim (\unicode[STIX]{x1D6FC}\unicode[STIX]{x0394}r/a)^{2}$ where $\unicode[STIX]{x0394}r$ is the width of the low-shear region. In the following analysis we assume that $Q\sim \unicode[STIX]{x1D700}^{2}$ . If $\unicode[STIX]{x1D6FC}\sim 1$ and the pressure gradient is localised within a narrow region, the assumption of concentric circular magnetic surfaces still holds (Greene & Chance 1981; Connor, Ham & Hastie 2016). From the equation for the Shafranov shift $\unicode[STIX]{x1D6E5}_{s}$ (assumed of order $\unicode[STIX]{x1D700}$ ):

(2.10) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6E5}_{s}^{\prime \prime }+3\unicode[STIX]{x1D6E5}_{s}^{\prime }/r=1/R_{0}+\unicode[STIX]{x1D6FC}/r, & \displaystyle\end{eqnarray}$$

we approximate $r\unicode[STIX]{x1D6E5}_{s}^{\prime \prime }\approx \unicode[STIX]{x1D6FC}$ and $T_{0}^{\prime }=-p_{0}^{\prime }R_{0}/B_{0}\sim \unicode[STIX]{x1D6FC}$ ( $T_{0}=B_{0,\unicode[STIX]{x1D711}}$ being the toroidal covariant component of the equilibrium magnetic field (Greene & Weimer 1971; Connor et al. 2016)).

Neglecting effects linear with respect to the magnetic shear (i.e. assume a flat $q$ ) the elements of the metric tensor in our straight field line coordinate system are evaluated up to order $\unicode[STIX]{x1D700}$ giving:

(2.11) $$\begin{eqnarray}\displaystyle & \displaystyle g_{rr}=1+\unicode[STIX]{x1D6FC}^{2}\sin ^{2}\unicode[STIX]{x1D717}-2\unicode[STIX]{x1D6E5}_{s}^{\prime }\cos \unicode[STIX]{x1D717}+4\left({\displaystyle \frac{r}{R_{0}}}-\unicode[STIX]{x1D6E5}_{s}^{\prime }\right)\unicode[STIX]{x1D6FC}\sin ^{2}\unicode[STIX]{x1D717}, & \displaystyle\end{eqnarray}$$
(2.12) $$\begin{eqnarray}\displaystyle & \displaystyle g_{r\unicode[STIX]{x1D717}}=\left(r\unicode[STIX]{x1D6FC}+2{\displaystyle \frac{r^{2}}{R_{0}}}-2r\unicode[STIX]{x1D6E5}_{s}^{\prime }\right)\sin \unicode[STIX]{x1D717}+\left({\displaystyle \frac{r}{R_{0}}}+\unicode[STIX]{x1D6E5}_{s}^{\prime }\right)r\unicode[STIX]{x1D6FC}\sin \unicode[STIX]{x1D717}\cos \unicode[STIX]{x1D717}, & \displaystyle\end{eqnarray}$$
(2.13) $$\begin{eqnarray}\displaystyle & \displaystyle g_{\unicode[STIX]{x1D717}\unicode[STIX]{x1D717}}=r^{2}+{\displaystyle \frac{2r^{3}}{R_{0}}}\cos \unicode[STIX]{x1D717}+2r^{2}\unicode[STIX]{x1D6E5}_{s}^{\prime }\cos \unicode[STIX]{x1D717}, & \displaystyle\end{eqnarray}$$
(2.14a,b ) $$\begin{eqnarray}\displaystyle & \displaystyle g_{\unicode[STIX]{x1D711}\unicode[STIX]{x1D711}}=R_{0}^{2}\left(1+{\displaystyle \frac{2r}{R_{0}}}\cos \unicode[STIX]{x1D717}\right),\quad 1/\sqrt{g}={\displaystyle \frac{1}{rR_{0}}}\left(1-{\displaystyle \frac{2r}{R_{0}}}\cos \unicode[STIX]{x1D717}\right), & \displaystyle\end{eqnarray}$$

where the ratio $g_{\unicode[STIX]{x1D711}\unicode[STIX]{x1D711}}/\sqrt{g}$ depends only on the flux label $r$ having used (2.10) and $\langle R^{2}\rangle ^{\prime }\approx -R_{0}\unicode[STIX]{x1D6FC}$ , where $\langle \boldsymbol{\cdot }\,\rangle :=\,1/2\unicode[STIX]{x03C0}\int _{0}^{2\unicode[STIX]{x03C0}}(\boldsymbol{\cdot })\,\text{d}\unicode[STIX]{x1D717}$ . We see that $g_{rr}$ and $g_{r\unicode[STIX]{x1D717}}$ are modified compared to their standard low $\unicode[STIX]{x1D6FD}$ expressions (Lazzaro & Zanca 2003; Brunetti et al. 2014). The latter are recovered by allowing $\unicode[STIX]{x1D6FC}$ to become small. Note that ellipticity and triangularity have not been included since they do not play a role in the coupling mechanism between $\pm 1$ neighbouring sidebands. This is because of the elongation and triangularity angular dependencies of the type $\cos 2\unicode[STIX]{x1D717}$ and $\sin 2\unicode[STIX]{x1D717}$ respectively (Connor et al. 2016).

For $\unicode[STIX]{x1D6FC}\lesssim 1$ we approximate $g_{rr}\simeq 1$ and at leading order the set of the equations describing the perturbation is still given by (2.8) and (2.9). Conversely, assuming that $\unicode[STIX]{x1D6FC}\gtrsim 1$ in a narrow region of width $\unicode[STIX]{x0394}r$ , we introduce the following orderings:

(2.15a-d ) $$\begin{eqnarray}\displaystyle & \displaystyle r{\displaystyle \frac{\text{d}}{\text{d}r}}\sim {\displaystyle \frac{r}{\unicode[STIX]{x0394}r}}\sim {\displaystyle \frac{1}{\unicode[STIX]{x1D700}}},\quad Q\sim \unicode[STIX]{x1D700}^{2},\quad m_{0}\sim \unicode[STIX]{x1D6FC}\sim 1,\quad {\displaystyle \frac{X_{\pm }}{X_{m_{0}}}}\sim \unicode[STIX]{x1D700}. & \displaystyle\end{eqnarray}$$

By means of (2.15) and employing the large $\unicode[STIX]{x1D6FC}$ metric tensor, at leading order we can derive the following coupled equations (we approximate $r_{\ast }\sim a$ ):

(2.16) $$\begin{eqnarray}\displaystyle & \displaystyle a^{2}[QX_{m_{0}}^{\prime }]^{\prime }-{\displaystyle \frac{\unicode[STIX]{x1D6FC}^{2}}{2}}X_{m_{0}}+{\displaystyle \frac{\unicode[STIX]{x1D6FC}}{2}}\left[{\displaystyle \frac{aX_{+}^{\prime }}{1+m_{0}}}+{\displaystyle \frac{aX_{-}^{\prime }}{1-m_{0}}}\right]=0, & \displaystyle\end{eqnarray}$$
(2.17) $$\begin{eqnarray}\displaystyle & \displaystyle aX_{\pm }^{\prime \prime }={\displaystyle \frac{1\pm m_{0}}{2}}[\unicode[STIX]{x1D6FC}X_{m_{0}}]^{\prime }. & \displaystyle\end{eqnarray}$$

The two simpler equations above can be used as an alternative to (2.8) and (2.9) in order to make analytically treatable more complex/realistic profiles. This will be discussed in § 5.2. We arrive to the same set of equations also when we have $1/m\sim \unicode[STIX]{x0394}r/r\sim \unicode[STIX]{x1D700}^{1/2}$ (in this case $X_{\pm }/X_{m_{0}}\sim \unicode[STIX]{x1D700}^{1/2}$ ). Pushing further the ordering of $m$ , i.e. $1/m\sim \unicode[STIX]{x1D700}$ , we enter in a ballooning-like configuration, where the weight of all the harmonics is approximately the same. This problem is not addressed and is beyond the scope of the paper.

Therefore equations (2.7), (2.8) and (2.9) (or alternatively (2.16) and (2.17)) form the basis for the analysis developed in the next sections. The solutions to this set of equations must fulfil the boundary conditions imposed at the plasma-‘outer world’ interface at $r=a$ . These boundary conditions are obtained by solving the perturbation in the vacuum region. Hence the solution to the problem proceeds as follows: first the equations for the perturbations are solved in the sheared and vacuum regions (with the appropriate boundary conditions at the magnetic axis and at the external metallic shell). Then equations (2.8) and (2.9) (or (2.16) and (2.17)) are solved in the low-shear region. Matching of the eigensolutions across the transition points yields eventually the dispersion relation.

3 Sheared region eigenfunctions

The aim of this section is to derive the shape of the main mode $m_{0}$ and the logarithmic derivative of the neighbouring $X_{\pm }$ . These quantities are required for the derivation of the dispersion relation when matching of the eigenfunctions at $r_{\ast }$ is required.

The profile of the rotational transform in $0<r<r_{\ast }$ is taken of the form (Mikhailovskii 1998):

(3.1) $$\begin{eqnarray}\displaystyle & \displaystyle (m_{0}-1)\unicode[STIX]{x1D707}-n=S[1-(r/r_{s})^{\unicode[STIX]{x1D706}}], & \displaystyle\end{eqnarray}$$

where $r_{s}$ is the position of the resonant surface of the mode $(m_{0}-1,n)$ (cf. figure 1). This choice for the safety factor allows the eigenproblem to be analytically solvable in terms of the hypergeometric equation. By requiring that $\unicode[STIX]{x1D707}(r_{\ast })\approx n/m_{0}$ we obtain $S=(n/m_{0})/[(r_{\ast }/r_{s})^{\unicode[STIX]{x1D706}}-1]$ . Note that in this approximation $\unicode[STIX]{x1D6FF}q$ corrections appearing in the logarithmic derivatives of the eigenfunctions are neglected. We shall analyse each harmonic separately. Note that an alternative model (though less realistic) for the safety factor which provides simpler expressions for the logarithmic jumps of the $m_{0}\pm 1$ sidebands has been used in Brunetti et al. (2018).

3.1 Main harmonic $m_{0}$

We multiply (2.7) by $X_{m_{0}}$ and integrate between 0 and $r_{\ast }$  (Hazeltine & Meiss 1992; Mikhailovskii 1998) yielding:

(3.2) $$\begin{eqnarray}\displaystyle & \displaystyle r^{3}k_{||,m_{0}}^{2}X_{m_{0}}\left.{\displaystyle \frac{\text{d}X_{m_{0}}}{\text{d}r}}\right|_{0}^{r_{\ast }}-\int _{0}^{r_{\ast }}\,\text{d}r\left[r^{3}k_{||,m_{0}}^{2}\left|{\displaystyle \frac{\text{d}X_{m_{0}}}{\text{d}r}}\right|^{2}+r(m_{0}^{2}-1)k_{m_{0},m_{0}}^{2}|X_{m_{0}}|^{2}\right]=0, & \displaystyle\end{eqnarray}$$

where the boundary condition $X_{m_{0}}(0)=0$ is assumed. Since $[r^{3}k_{||,m_{0}}^{2}X_{m}(\text{d}X_{m}/\text{d}r)]_{r_{\ast }}\sim k_{||,m_{0}}(r_{\ast })^{2}\ll 1$ (with the reasonable assumption that the derivatives are not diverging), we are left with an integral of positive terms which must be zero. The only possibility is that the function under the sign of integration is vanishing, this is:

(3.3) $$\begin{eqnarray}\displaystyle & \displaystyle X_{m_{0}}=0. & \displaystyle\end{eqnarray}$$

This is in accordance with numerical calculations presented e.g. in Zheng et al. (2013a , b ). Note that the equation above still holds if the $q=m_{0}/n$ surface is crossed, viz. $\unicode[STIX]{x1D6FF}q>0$ , if the resonant point of the main $m_{0}$ mode (to whom we refer with $r_{m}$ ) is not far from $r_{\ast }$ (i.e. $r_{m}\approx r_{\ast }$ ) (Gimblett et al. 1996).

3.2 Lower sideband ( $m_{-}=m_{0}-1$ )

With the choice of the rotational transform given by (3.1), the expressions for the sidebands are readily obtained. Introducing the variable $z=(r/r_{s})^{\unicode[STIX]{x1D706}}$ , when $z<1$ , the solution of (2.7) for the lower sideband fulfilling the boundary condition on the magnetic axis ( $X_{-}(0)<\infty$ ) is (Kuvshinov & Mikhailovskii 1991; Mikhailovskii 1998; Brunetti et al. 2014):

(3.4) $$\begin{eqnarray}\displaystyle & \displaystyle X_{-}=A_{{<}}z^{(m_{-}-1)/\unicode[STIX]{x1D706}}(1-z)^{-1}F(\unicode[STIX]{x1D702},\unicode[STIX]{x1D701};\unicode[STIX]{x1D702}+\unicode[STIX]{x1D701}+1;z), & \displaystyle\end{eqnarray}$$

where $F$ is the hypergeometric function (see Abramowitz & Stegun 1968, p. 555) and we defined $\unicode[STIX]{x1D702}=(m_{-}-\bar{m}_{-})/\unicode[STIX]{x1D706}$ , $\unicode[STIX]{x1D701}=(m_{-}+\bar{m}_{-})/\unicode[STIX]{x1D706}$ , $\bar{m}_{-}=\sqrt{m_{-}^{2}+2\unicode[STIX]{x1D706}+\unicode[STIX]{x1D706}^{2}}$ . Conversely, when $z>1$ , the solution of (3.3) reads (Kuvshinov & Mikhailovskii 1991; Mikhailovskii 1998; Brunetti et al. 2014):

(3.5) $$\begin{eqnarray}\displaystyle & \displaystyle X_{-}={\displaystyle \frac{z^{-(1+\bar{m}_{-})/\unicode[STIX]{x1D706}}}{z-1}}(A_{{>}}^{\ast }F(\unicode[STIX]{x1D701},-\unicode[STIX]{x1D702};1+\unicode[STIX]{x1D701}-\unicode[STIX]{x1D702};1/z)+B_{{>}}^{\ast }z^{\unicode[STIX]{x1D701}-\unicode[STIX]{x1D702}}F(-\unicode[STIX]{x1D701},\unicode[STIX]{x1D702};1+\unicode[STIX]{x1D702}-\unicode[STIX]{x1D701};1/z)), & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

where $A_{{>}}^{\ast }=-A_{{>}}(\unicode[STIX]{x1D702}\unicode[STIX]{x1D701}\unicode[STIX]{x1D6E4}(-\unicode[STIX]{x1D702})\unicode[STIX]{x1D6E4}(\unicode[STIX]{x1D701}))/(\unicode[STIX]{x1D6E4}(1-\unicode[STIX]{x1D702}+\unicode[STIX]{x1D701}))$ and $B_{{>}}^{\ast }=-B_{{>}}(\unicode[STIX]{x1D702}\unicode[STIX]{x1D701}\unicode[STIX]{x1D6E4}(\unicode[STIX]{x1D702})\unicode[STIX]{x1D6E4}(-\unicode[STIX]{x1D701}))/$ $(\unicode[STIX]{x1D6E4}(1+\unicode[STIX]{x1D702}-\unicode[STIX]{x1D701}))$ .

The asymptotic behaviour of (3.4) and (3.5) close to $r_{\ast }$ , depends on the ratio $B_{{>}}/A_{{>}}$ . Neglecting for sake of simplicity inertial corrections, $B_{{>}}/A_{{>}}$ is found by imposing that for an ideal mode the displacement is finite at its own rational surface (modifications of the ratio $B_{{>}}/A_{{>}}$ due to inertial corrections or residual coupling effects are discussed in appendix A). Thus we immediately obtain $B_{{>}}/A_{{>}}\rightarrow -1$ and $A_{{<}}=0$ .

In the limit $(r_{s}/r_{\ast })^{\unicode[STIX]{x1D706}}\ll 1$ , far from the resonant surface ( $z\gg 1$ ), we have:

(3.6) $$\begin{eqnarray}\displaystyle & \displaystyle X_{-}\sim \left({\displaystyle \frac{r}{r_{\ast }}}\right)^{-(\bar{m}_{-}+\ell )}+C_{0}\left({\displaystyle \frac{r}{r_{\ast }}}\right)^{\bar{m}_{-}-\ell }, & \displaystyle\end{eqnarray}$$

with $C_{0}=-(\unicode[STIX]{x1D6E4}(\unicode[STIX]{x1D702})\unicode[STIX]{x1D6E4}(-\unicode[STIX]{x1D701})\unicode[STIX]{x1D6E4}(1-\unicode[STIX]{x1D702}+\unicode[STIX]{x1D701}))/(\unicode[STIX]{x1D6E4}(-\unicode[STIX]{x1D702})\unicode[STIX]{x1D6E4}(\unicode[STIX]{x1D701})\unicode[STIX]{x1D6E4}(1+\unicode[STIX]{x1D702}-\unicode[STIX]{x1D701}))\left(r_{\ast }/r_{s}\right)^{2\bar{m}_{-}}$ and $\ell =\unicode[STIX]{x1D706}+1$ . Thus by means of (3.6) it is straightforward to obtain:

(3.7) $$\begin{eqnarray}\displaystyle & \displaystyle \mathbb{C}_{-}={\displaystyle \frac{r_{\ast }X_{-}^{\prime }(r_{\ast })}{X_{-}(r_{\ast })}}=(\bar{m}_{-}){\displaystyle \frac{C_{0}-1}{C_{0}+1}}-\ell . & \displaystyle\end{eqnarray}$$

When $(r_{s}/r_{\ast })^{\unicode[STIX]{x1D706}}\sim 1$ it is possible to show that $\mathbb{C}_{-}\sim \text{O}(z-1)$ . However since generally speaking the ratio $(r_{s}/r_{\ast })^{\unicode[STIX]{x1D706}}$ is neither small nor close to unity, in the numerical computation of the dispersion relation we shall use the complete expression obtained from (3.5).

3.3 Upper sideband ( $m_{+}=m_{0}+1$ )

It is easy to show that with the safety factor given by (3.1) the resonant surface of the mode $(m_{0}+1,n)$ occurs at $r_{+}=r_{s}[1+2m_{0}[(r_{\ast }/r_{s})^{\unicode[STIX]{x1D706}}-1]/(m_{0}+1)]^{1/\unicode[STIX]{x1D706}}$ . Hence the expression for the upper sideband regular on the axis ( $X_{+}<\infty$ ) reads:

(3.8) $$\begin{eqnarray}\displaystyle & \displaystyle X_{+}\sim \tilde{z}^{(m_{+}-1)/\unicode[STIX]{x1D706}}(1-\tilde{z})^{-1}F(\tilde{\unicode[STIX]{x1D702}},\tilde{\unicode[STIX]{x1D701}};\tilde{\unicode[STIX]{x1D702}}+\tilde{\unicode[STIX]{x1D701}}+1;\tilde{z}), & \displaystyle\end{eqnarray}$$

where $\tilde{z}=(r/r_{+})^{\unicode[STIX]{x1D706}}$ , $\tilde{\unicode[STIX]{x1D702}}=(m_{+}-\bar{m}_{+})/\unicode[STIX]{x1D706}$ , $\tilde{\unicode[STIX]{x1D701}}=(m_{+}+\bar{m}_{+})/\unicode[STIX]{x1D706}$ , $\bar{m}_{+}=\sqrt{m_{+}^{2}+2\unicode[STIX]{x1D706}+\unicode[STIX]{x1D706}^{2}}$ . Thus if $(r_{\ast }/r_{s})^{\unicode[STIX]{x1D706}}\gg 1$ and $m_{0}\gg 1$ then $r_{+}/r_{\ast }\sim 2^{1/\unicode[STIX]{x1D706}}$ so that from (3.8) we obtain:

(3.9) $$\begin{eqnarray}\displaystyle & \displaystyle \mathbb{C}_{+}={\displaystyle \frac{r_{\ast }X_{+}^{\prime }(r_{\ast })}{X_{+}(r_{\ast })}}\approx (m_{+}-1)+\unicode[STIX]{x1D706}+{\displaystyle \frac{\unicode[STIX]{x1D706}\tilde{\unicode[STIX]{x1D702}}\tilde{\unicode[STIX]{x1D701}}}{1+\tilde{\unicode[STIX]{x1D702}}+\tilde{\unicode[STIX]{x1D701}}}}{\displaystyle \frac{F(\tilde{\unicode[STIX]{x1D702}}+1,\tilde{\unicode[STIX]{x1D701}}+1;\tilde{\unicode[STIX]{x1D702}}+\tilde{\unicode[STIX]{x1D701}}+2;{\textstyle \frac{1}{2}})}{2F(\tilde{\unicode[STIX]{x1D702}},\tilde{\unicode[STIX]{x1D701}};\tilde{\unicode[STIX]{x1D702}}+\tilde{\unicode[STIX]{x1D701}}+1;\displaystyle {\textstyle \frac{1}{2}})}}. & \displaystyle\end{eqnarray}$$

Conversely in the limit $(r_{\ast }/r_{s})^{\unicode[STIX]{x1D706}}\sim 1+\unicode[STIX]{x1D6FF}$ with $\unicode[STIX]{x1D6FF}\ll 1$ we have $(r_{+}/r_{\ast })^{\unicode[STIX]{x1D706}}\sim 1+\unicode[STIX]{x1D6FF}_{\ast }$ with $\unicode[STIX]{x1D6FF}_{\ast }=(m_{0}-1)/(m_{0}+1)\unicode[STIX]{x1D6FF}$ , which yields:

(3.10) $$\begin{eqnarray}\displaystyle & \displaystyle \mathbb{C}_{+}\approx {\displaystyle \frac{\unicode[STIX]{x1D706}}{\unicode[STIX]{x1D6FF}_{\ast }}}+m_{+}-1-\unicode[STIX]{x1D706}-\unicode[STIX]{x1D706}\tilde{\unicode[STIX]{x1D702}}\tilde{\unicode[STIX]{x1D701}}[2\unicode[STIX]{x1D6FE}_{E}+\unicode[STIX]{x1D6F9}(\tilde{\unicode[STIX]{x1D702}}+1)+\unicode[STIX]{x1D6F9}(\tilde{\unicode[STIX]{x1D701}}+1)+\ln \unicode[STIX]{x1D6FF}_{\ast }]. & \displaystyle\end{eqnarray}$$

In order to avoid approximations which could introduce unphysical behaviours, as in the case for the lower sideband, the full expression for $X_{+}$ (i.e. (3.8)) is preferred when the growth rate is computed numerically.

In the next sections the quantities $\mathbb{C}_{\pm }$ will be matched to the low-shear region solutions and then used for the evaluation of the dispersion relation.

4 Vacuum region

In the vacuum region the magnetic field must fulfil the conditions $\unicode[STIX]{x1D735}\times \boldsymbol{B}=\unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{B}=0$ . Hence we write $\boldsymbol{B}=\unicode[STIX]{x1D735}\unicode[STIX]{x1D712}$ with the constraint $\unicode[STIX]{x1D735}^{2}\unicode[STIX]{x1D712}=0$  (Wesson 1978). In large aspect ratio and under the assumption $m>n$ for each Fourier harmonic, the equation determining the $m$ th component of the vacuum perturbation reads:

(4.1) $$\begin{eqnarray}\displaystyle & \displaystyle r(r\unicode[STIX]{x1D712}_{m}^{\prime })^{\prime }-m^{2}\unicode[STIX]{x1D712}_{m}=0. & \displaystyle\end{eqnarray}$$

Recalling that $\tilde{B}_{m}^{r}=\unicode[STIX]{x1D712}_{m}^{\prime }$ , this gives $\tilde{B}_{m}^{r}\sim (r/b)^{-m-1}-D(r/b)^{m-1}$ for $a<r<b$ where $D$ is a constant determined by the boundary conditions at the plasma–metallic wall interface. In the region $r>b+d$ the solution of the vacuum perturbation is written in terms of modified Bessel functions (Lashmore-Davies 2001). However for $r\gtrsim b+d$ the behaviour of the radial perturbed magnetic field is well approximated by $\tilde{B}_{m}^{r}\sim (r/b)^{-m-1}$  (Mikhailovskii 1998; Lashmore-Davies 2001).

By means of Faraday’s law $\boldsymbol{E}=\unicode[STIX]{x1D702}\boldsymbol{J}$ , and assuming that the radial derivatives of the perturbation are dominant, the equation for the perturbation within the wall is:

(4.2) $$\begin{eqnarray}\displaystyle & \displaystyle {\displaystyle \frac{\text{d}^{2}\tilde{B}^{r}}{\text{d}r^{2}}}=\unicode[STIX]{x1D70E}\unicode[STIX]{x1D6FE}\tilde{B}^{r}. & \displaystyle\end{eqnarray}$$

In the thin wall approximation ( $d/b\ll 1$ ) (Gimblett 1986) we integrate (4.2) across the wall, so that $[r(\tilde{B}_{m}^{r})^{\prime }/\tilde{B}_{m}^{r}]_{b}^{b+d}=(\unicode[STIX]{x1D6FE}/\unicode[STIX]{x1D714}_{A})S_{w}(d/b)$ with $S_{w}=\unicode[STIX]{x1D70E}b^{2}\unicode[STIX]{x1D714}_{A}$ . Here $\unicode[STIX]{x1D714}_{A}$ is computed in the plasma core. Requiring continuity of $\tilde{B}_{m}^{r}$ across the thin wall and employing the wall jump condition gives $D=F/(F+2m)$ where $F=(\unicode[STIX]{x1D6FE}/\unicode[STIX]{x1D714}_{A})S_{w}(d/b)$ . We stress the point that the constant $D$ depends on $m$ , thus it varies for different Fourier modes. Note that if we consider the case of a perfectly conducting wall (i.e. $\unicode[STIX]{x1D702}\rightarrow 0$ , $S_{w}\rightarrow \infty$ ) we must have $D=1$ . Hence $\tilde{B}_{m}^{r}=0$ for any $m$ which corresponds to the condition $\hat{\boldsymbol{n}}\boldsymbol{\cdot }\boldsymbol{B}=0$ where $\hat{\boldsymbol{n}}$ is the normal vector pointing outward the $r=a$ surface (we identify $\hat{\boldsymbol{n}}\equiv \unicode[STIX]{x1D735}r$ ) (Bernstein et al. 1958; Wesson 1978; Freidberg 1987).

The jump condition at the plasma–vacuum interface is (Bernstein et al. 1958; Freidberg 1987):

(4.3) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x27E6}\hat{\boldsymbol{n}}\boldsymbol{\cdot }\boldsymbol{B}\unicode[STIX]{x27E7}_{a}=0, & \displaystyle\end{eqnarray}$$

where $\unicode[STIX]{x27E6}\cdot \unicode[STIX]{x27E7}_{a}=(\cdot )_{a+\unicode[STIX]{x1D6FF}}-(\cdot )_{a-\unicode[STIX]{x1D6FF}}$ with $\unicode[STIX]{x1D6FF}\rightarrow 0$ . By means of (2.6) we can extend the definition of perturbed displacement $X$ outside the plasma, i.e. we write $\tilde{B}^{r}$ in terms of $X$ . This yields $\tilde{B}_{m}^{r}\sim k_{||,m}X_{m}$ so that since $\tilde{B}_{m}^{r}$ must be continuous (cf. (4.3)), $X_{m}$ also is continuous since $k_{||}$ is continuous by hypothesis. It follows that the vacuum perturbation fulfils (2.7). Thus it is easy to see that in the vacuum for a generic Fourier mode $m$ we have (Wesson 1978; Mikhailovskii 1998):

(4.4) $$\begin{eqnarray}\displaystyle & \displaystyle \left.{\displaystyle \frac{rX_{m}^{\prime }}{X_{m}}}\right|_{a+\unicode[STIX]{x1D6FF}}={\displaystyle \frac{2\unicode[STIX]{x1D707}_{\ast }}{\unicode[STIX]{x1D707}_{\ast }-n/m}}-{\displaystyle \frac{m+1+\displaystyle {\displaystyle \frac{F}{F+2m}}(m-1)\left(a/b\right)^{2m}}{1-\displaystyle {\displaystyle \frac{F}{F+2m}}\left(a/b\right)^{2m}}}, & \displaystyle\end{eqnarray}$$

where in the vacuum $\unicode[STIX]{x1D707}=\unicode[STIX]{x1D707}_{\ast }a^{2}/r^{2}$ (i.e. vanishing toroidal current) (Mikhailovskii 1998). Let us introduce the notation $\mathbb{B}_{\pm }=(rX_{m_{0}\pm 1}^{\prime }/X_{m_{0}\pm 1})|_{a+\unicode[STIX]{x1D6FF}}$ . We point out that, as calculated in § 3, because we approximate $\unicode[STIX]{x1D707}_{\ast }\approx n/m_{0}$ , $\unicode[STIX]{x1D6FF}q$ corrections $\mathbb{B}_{\pm }$ are not taken into account.

For the Fourier mode $m_{0}$ , since $\unicode[STIX]{x1D707}_{\ast }-n/m_{0}\sim \unicode[STIX]{x1D6FF}q\ll 1$ we set $X_{m}(a)\approx 0$ . This can be deduced by the fact that close to $a$ the eigenfunction behaves as $X_{m_{0}}(r)\sim (X_{m_{0}}(a)/\unicode[STIX]{x1D6FF}q)[1-(a/b)^{2m_{0}}D]$ , where obviously $D$ is computed for the mode $m_{0}$ . According to (4.3) the displacement is continuous across $a$ , thus in order to prevent $X_{m_{0}}$ to become arbitrarily large when $\unicode[STIX]{x1D6FF}q\ll 1$ , in general we must set $X_{m_{0}}(a)=0$ . In the case of an ideal or nearly ideal wall we note that in the vacuum region the equation describing the perturbation is (2.7) (having introduced the quantity $X_{m}$ also in the vacuum). Thus if we multiply (2.7) by $X_{m_{0}}$ and integrate by parts from $a$ to $b$ we obtain (3.2) with the replacements $r_{\ast }\rightarrow a$ and $a\rightarrow b$ . Since $k_{||,m_{0}}^{2}(a)\ll 1$ and $X_{m_{0}}(b)=0$ eventually to leading order we get $X_{m_{0}}=0$ . We shall still approximate $X_{m_{0}}(a)\approx 0$ when the resonance $q=m_{0}/n$ is in the vacuum gap if the resonant point is sufficiently close to $a$ .

As in the treatment of the sheared region, the logarithmic derivatives of the sideband harmonics are required for matching with the solutions in the low-shear region when the dispersion relation is derived. This will be shown in the next section.

5 Low-shear region matching and dispersion relation

In this section we solve for the main mode and the sideband harmonics deriving eventually the dispersion relation. In this region of width $a-r_{\ast }=\unicode[STIX]{x1D6E5}\sim \unicode[STIX]{x1D700}\ll 1$ , large pressure gradients drive large edge bootstrap current contributions which in turn flatten the safety factor. There are several choices for the equilibrium profiles for $p_{0}$ and $\unicode[STIX]{x1D70C}_{0}$ (e.g. linear, piecewise continuous polynomials etc.) which allow for an exact treatment of the perturbation. Our analysis concentrates first on a case in which the profiles are step-like and then on a more realistic tanh type profiles. The step-like case has been already treated in Brunetti et al. (2018) and it turns out to be useful when toroidal rotation is considered (this is not discussed here). In the next subsection we limit to summarise the findings and to show the mathematical procedure for the derivation of the dispersion relation. Part of these techniques are employed when the tanh model is considered, though a slightly different approach must be used.

5.1 Step-like functions

In order to mimic the abrupt decrease of the pressure and mass density in the flat- $q$ region (Burrell et al. 2001, 2005), following Wahlberg et al. (2013) we introduce step-like profiles:

(5.1) $$\begin{eqnarray}\displaystyle & \displaystyle p_{0}/p_{\ast }\sim \unicode[STIX]{x1D70C}/\unicode[STIX]{x1D70C}_{0}=\unicode[STIX]{x1D6E9}(r_{p}-r), & \displaystyle\end{eqnarray}$$

where $\unicode[STIX]{x1D6E9}(x)$ is the Heaviside step function of argument $x$ , $r_{p}=(a+r_{\ast })/2$ and $p_{\ast }$ is the value of plasma pressure at $r_{\ast }$ . The mass density is assumed constant in $0<r<r_{\ast }$ .

By integrating once (2.9) we obtain:

(5.2) $$\begin{eqnarray}\displaystyle & \displaystyle (r^{2\pm m_{0}}X_{m_{0}\pm 1})^{\prime }=r^{1\pm 2m_{0}}L_{\pm }+{\displaystyle \frac{1\pm m_{0}}{2}}\unicode[STIX]{x1D6FC}r^{1\pm m_{0}}X_{m_{0}}. & \displaystyle\end{eqnarray}$$

Integration of (2.9) across $r_{p}$ gives $\unicode[STIX]{x27E6}(r^{2\pm m_{0}}X_{\pm })^{\prime }\unicode[STIX]{x27E7}_{r_{p}}=0$ implying that in (5.2) we must have $L_{\pm }(r<r_{p})=L_{\pm }(r>r_{p})$ , i.e. the constant $L_{+}$ (or $L_{-}$ ) on the left and on the right of $r_{p}$ must be the same. Thus plugging (5.2) into (2.8), under the assumption $(1/q^{2}-1)\approx -1$ , yields:

(5.3) $$\begin{eqnarray}\displaystyle \left[r^{3}QX_{m_{0}}^{\prime }\right]^{\prime }+r\left[(1-m_{0}^{2})Q+\hat{\unicode[STIX]{x1D6FE}}^{2}{\displaystyle \frac{r\unicode[STIX]{x1D70C}_{0}^{\prime }}{\bar{\unicode[STIX]{x1D70C}}}}-{\displaystyle \frac{\unicode[STIX]{x1D6FC}r}{R_{0}}}\right]X_{m_{0}}+{\displaystyle \frac{\unicode[STIX]{x1D6FC}}{2}}\left[{\displaystyle \frac{r^{1+m_{0}}L_{+}}{1+m_{0}}}+{\displaystyle \frac{r^{1-m_{0}}L_{-}}{1-m_{0}}}\right]=0, & & \displaystyle\end{eqnarray}$$

where $\hat{\unicode[STIX]{x1D6FE}}^{2}=\unicode[STIX]{x1D6FE}^{2}(1+2q^{2})/(n\unicode[STIX]{x1D714}_{A})^{2}$ .

The three harmonics must be supplied with appropriate boundary conditions. The solution of the main harmonic $m_{0}$ in the sheared and vacuum regions (cf. §§ 3.1 and 4) provides the constraint:

(5.4) $$\begin{eqnarray}\displaystyle & \displaystyle X_{m_{0}}(r_{\ast })=X_{m_{0}}(a)=0. & \displaystyle\end{eqnarray}$$

The boundary conditions for the sideband harmonics are obtained by integrating (2.9) across $r_{\ast }$ and $a$ . Since the pressure and its gradient at these points are vanishingly small, this yields:

(5.5a,b ) $$\begin{eqnarray}\displaystyle & \displaystyle \left.{\displaystyle \frac{rX_{\pm }^{\prime }(r)}{X_{\pm }(r)}}\right|_{r\rightarrow r_{\ast }^{+}}=\mathbb{C}_{\pm },\quad \left.{\displaystyle \frac{rX_{\pm }^{\prime }(r)}{X_{\pm }(r)}}\right|_{r\rightarrow a^{-}}=\mathbb{B}_{\pm }, & \displaystyle\end{eqnarray}$$

where the quantities $\mathbb{C}_{\pm }$ and $\mathbb{B}_{\pm }$ have been evaluated in §§ 3 and 4 respectively.

Before solving for the main harmonic we first determine the constants $L_{\pm }$ . These are obtained by evaluating (5.2) at $r_{\ast }$ and $a$ and using the constraint (5.4). Since the pressure gradient and $X_{m_{0}}$ are vanishing at the plasma boundary, we obtain $X_{\pm }(r_{\ast })=r_{\ast }^{\pm m_{0}}L_{\pm }/(2\pm m_{0}+\mathbb{C}_{\pm })$ and $X_{\pm }(a)=a^{\pm m_{0}}L_{\pm }/(2\pm m_{0}+\mathbb{B}_{\pm })$ . By means of (5.1), integration of (5.2) from $r_{\ast }$ to $a$ finally gives $(r_{p}^{\pm m_{0}}L_{\pm }/1\pm m_{0})=X_{0}(\hat{\unicode[STIX]{x1D6FD}}/\unicode[STIX]{x1D700}_{p})\unicode[STIX]{x1D6EC}^{(\pm )}$ where $\hat{\unicode[STIX]{x1D6FD}}=2p_{\ast }q^{2}/B_{0}^{2}$ , $\unicode[STIX]{x1D700}_{p}=r_{p}/R_{0}$ and (Gimblett et al. 1996):

(5.6) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6EC}^{(\pm )}={\displaystyle \frac{(r_{p}/r_{\ast })^{2\pm 2m_{0}}(1\pm m_{0})[2\pm m_{0}+\mathbb{C}_{\pm }][2\pm m_{0}+\mathbb{B}_{\pm }]}{(\mathbb{C}_{\pm }\mp m_{0})[2\pm m_{0}+\mathbb{B}_{\pm }]-(\mathbb{B}_{\pm }\mp m_{0})[2\pm m_{0}+\mathbb{C}_{\pm }]\left({\displaystyle \frac{a}{r_{\ast }}}\right)^{2\pm 2m_{0}}}}. & \displaystyle\end{eqnarray}$$

Integration of (5.2) and (5.3) across $r_{p}$ with the profiles given in (5.1), shows that the singularities arising from $p_{0}^{\prime }$ and $\unicode[STIX]{x1D70C}_{0}^{\prime }$ produce discontinuities at this point of $X_{\pm }$ , $X_{\pm }^{\prime }$ and $X_{m_{0}}^{\prime }$ while $X_{m_{0}}$ remains continuous. By means of (5.1), the equation for $X_{m_{0}}$ is greatly simplified:

(5.7) $$\begin{eqnarray}\displaystyle & \displaystyle [r^{3}X_{m_{0}}^{\prime }]^{\prime }+r(1-m_{0}^{2})X_{m_{0}}=0. & \displaystyle\end{eqnarray}$$

The equation above is solved separately for $r<r_{p}$ and $r>r_{p}$ imposing continuity at $r_{p}$ and the constraints (5.4), giving:

(5.8) $$\begin{eqnarray}\displaystyle X_{m_{0}}=X_{0}\times \left\{\begin{array}{@{}ll@{}}\displaystyle {\displaystyle \frac{\left({\displaystyle \frac{r}{r_{\ast }}}\right)^{m_{0}-1}-\left({\displaystyle \frac{r}{r_{\ast }}}\right)^{-m_{0}-1}}{\left({\displaystyle \frac{r_{p}}{r_{\ast }}}\right)^{m_{0}-1}-\left({\displaystyle \frac{r_{p}}{r_{\ast }}}\right)^{-m_{0}-1}}}, & r<r_{p},\\[31.20007pt] \displaystyle {\displaystyle \frac{\left({\displaystyle \frac{r}{a}}\right)^{m_{0}-1}-\left({\displaystyle \frac{r}{a}}\right)^{-m_{0}-1}}{\left({\displaystyle \frac{r_{p}}{a}}\right)^{m_{0}-1}-\left({\displaystyle \frac{r_{p}}{a}}\right)^{-m_{0}-1}}}, & r>r_{p}.\end{array}\right. & & \displaystyle\end{eqnarray}$$

Following Wahlberg et al. (2013), the dispersion relation is obtained by integrating (5.3) across $r_{p}$ . This gives:

(5.9) $$\begin{eqnarray}\displaystyle & \displaystyle r_{p}\unicode[STIX]{x27E6}QX_{m_{0}}^{\prime }\unicode[STIX]{x27E7}_{r_{p}}-[\hat{\unicode[STIX]{x1D6FE}}^{2}+\hat{\unicode[STIX]{x1D6FD}}]X_{0}+\left({\displaystyle \frac{\hat{\unicode[STIX]{x1D6FD}}}{\unicode[STIX]{x1D700}_{p}}}\right)^{2}[\unicode[STIX]{x1D6EC}^{(+)}+\unicode[STIX]{x1D6EC}^{(-)}]X_{0}/2=0. & \displaystyle\end{eqnarray}$$

Let us introduce $g_{\pm }=(rX_{m_{0}}^{\prime }/X_{m_{0}})|_{r_{p}^{\pm }}$ . By means of (5.8), it is easy to show that if $r_{p}\approx a$ and $m_{0}\sim 1$ then $g_{-}\approx -g_{+}\approx 2a/\unicode[STIX]{x1D6E5}$ . Hence (5.9) becomes:

(5.10) $$\begin{eqnarray}\displaystyle & \displaystyle {\displaystyle \frac{\hat{\unicode[STIX]{x1D6FE}}^{2}}{2}}+\left({\displaystyle \frac{\unicode[STIX]{x1D6FF}q}{q}}\right)^{2}={\displaystyle \frac{\unicode[STIX]{x1D6E5}}{4a}}\left[{\displaystyle \frac{\hat{\unicode[STIX]{x1D6FD}}^{2}}{2\unicode[STIX]{x1D700}_{p}^{2}}}(\unicode[STIX]{x1D6EC}^{(+)}+\unicode[STIX]{x1D6EC}^{(-)})-\hat{\unicode[STIX]{x1D6FD}}\right]. & \displaystyle\end{eqnarray}$$

The last term on the right-hand side of (5.10) corresponds to the Mercier term in (2.8) and has a weak stabilising influence. It is clear that the instability drive are the pressure gradient and the field line bending weakening. This dispersion relation has been discussed in detail in Brunetti et al. (2018). We shall now use this result as a reference for the more realistic tanh case presented in the next section.

5.2 Tanh model

In this section we model $p_{0}$ and $\unicode[STIX]{x1D70C}_{0}$ with more realistic smooth profiles. It is important to note that since the step and the corresponding discontinuities at $r_{p}$ are lost, in order to derive the dispersion relation we must adopt a different procedure compared to the one employed in the previous case.

The total pressure is written as $p_{0}=n_{0}(T_{i}+T_{e})$ where $T_{s}$ is the temperature of the $s$ species and $n_{0}$ the numerical density (quasineutrality is imposed). Assuming that the electron temperature profile is proportional to the density profile, i.e. $T_{e}\approx T_{e0}n_{0}/\bar{n}$ ( $\bar{n}$ is the value of the numerical density in the core), and $T_{i}\gg T_{e0}$ with $T_{i}^{\prime }\sim 0$ we eventually obtain $p_{0}\sim n_{0}T_{i}$ and $p_{0}^{\prime }\sim n_{0}^{\prime }T_{i}$ . Hence pressure and mass density have approximately the same shape (Zheng et al. 2013a , b ) in qualitative accordance with experimental data (Burrell et al. 2001, 2005). Thus a more realistic choice for the profiles of $p_{0}$ and $\unicode[STIX]{x1D70C}_{0}$ is the following:

(5.11) $$\begin{eqnarray}\displaystyle & \displaystyle p_{0}/p_{\ast }\sim \unicode[STIX]{x1D70C}_{0}/\bar{\unicode[STIX]{x1D70C}}={\textstyle \frac{1}{2}}[1-\tanh [(r-r_{p})/\unicode[STIX]{x1D6FF}]], & \displaystyle\end{eqnarray}$$

where $\unicode[STIX]{x1D6FF}/\unicode[STIX]{x1D6E5}\lesssim 1$ and $r_{p}$ has been already defined in the previous section (we recall that $\unicode[STIX]{x1D6E5}$ is the width of the low-shear region). Note that for $\unicode[STIX]{x1D6FF}\rightarrow 0$ we shall recover the step-like case (this will be indeed shown later).

In this region we adopt the ordering presented in (2.15), and therefore we employ (2.16) and (2.17). As in the previous section, we integrate once (2.17) and we plug the result into (2.16) obtaining:

(5.12) $$\begin{eqnarray}\displaystyle & \displaystyle a^{2}[QX_{m_{0}}]^{\prime }+{\displaystyle \frac{\unicode[STIX]{x1D6FC}}{2}}\left[{\displaystyle \frac{\hat{L}_{+}}{1+m_{0}}}+{\displaystyle \frac{\hat{L}_{-}}{1-m_{0}}}\right]=0, & \displaystyle\end{eqnarray}$$
(5.13) $$\begin{eqnarray}\displaystyle & \displaystyle aX_{\pm }^{\prime }=\hat{L}_{\pm }+{\displaystyle \frac{1\pm m_{0}}{2}}\unicode[STIX]{x1D6FC}X_{m_{0}}, & \displaystyle\end{eqnarray}$$

where explicitly $Q=Z[1-c\tanh x/\unicode[STIX]{x1D6FF}]$ with $Z=(\unicode[STIX]{x1D6FF}q/q)^{2}+\hat{\unicode[STIX]{x1D6FE}}^{2}/2$ and $c=\hat{\unicode[STIX]{x1D6FE}}^{2}/2Z$ . The solution of (5.12) is easily obtained and it reads:

(5.14) $$\begin{eqnarray}\displaystyle & \displaystyle X_{m_{0}}=c_{1}+\left(c_{2}+{\displaystyle \frac{H\unicode[STIX]{x1D6FF}}{c}}\right){\displaystyle \frac{x}{\unicode[STIX]{x1D6FF}}}+c_{2}\left\{c\ln \left[\cosh \left({\displaystyle \frac{x}{\unicode[STIX]{x1D6FF}}}\right)-c\sinh \left({\displaystyle \frac{x}{\unicode[STIX]{x1D6FF}}}\right)\right]\right\}, & \displaystyle\end{eqnarray}$$

where $x=r-r_{p}$ and $H=(q^{2}p_{\ast }/2)/(B_{0}^{2}Z\unicode[STIX]{x1D700})[(\hat{L}_{+}/a)/(1+m_{0})+(\hat{L}_{-}/a)/(1-m_{0})]$ .

The constants $c_{1,2}$ of the main harmonic are determined by imposing the boundary conditions (5.4) yielding ( $h=\unicode[STIX]{x1D6E5}/2$ ):

(5.15a,b ) $$\begin{eqnarray}\displaystyle & \displaystyle c_{1}={\displaystyle \frac{hH\ln [1+(1-c^{2})\sinh ^{2}(h/\unicode[STIX]{x1D6FF})]/2}{h/\unicode[STIX]{x1D6FF}-c\tanh ^{-1}[c\tanh (h/\unicode[STIX]{x1D6FF})]}},\quad c_{2}=-{\displaystyle \frac{hH/c}{h/\unicode[STIX]{x1D6FF}-\tanh ^{-1}[c\tanh (h/\unicode[STIX]{x1D6FF})]}}. & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

By integrating (2.17) across $a$ and $r_{\ast }$ under the assumption that the pressure gradient is not too large at these points, employing (5.4) we find that the sidebands fulfil (5.5).

As shown in the previous section, to compute the constants $\hat{L}_{\pm }$ equation (5.13) is first evaluated both at $r_{\ast }$ and $a$ and then integrated from $r_{\ast }$ to $a$ yielding:

(5.16) $$\begin{eqnarray}\displaystyle & \displaystyle \hat{L}_{\pm }={\displaystyle \frac{\mathbb{B}_{\pm }\mathbb{C}_{\pm }(1\pm m_{0})/(2a)}{[\mathbb{C}_{\pm }-\mathbb{B}_{\pm }-\mathbb{B}_{\pm }\mathbb{C}_{\pm }\unicode[STIX]{x1D6E5}/a]}}\int _{r_{\ast }}^{a}\unicode[STIX]{x1D6FC}X_{m_{0}}\,\text{d}r. & \displaystyle\end{eqnarray}$$

We immediately note that $\hat{L}_{\pm }(X_{m_{0}})$ are linear in $X_{m_{0}}$ , i.e. for $k$ constant $k\hat{L}_{\pm }(X_{m_{0}})=\hat{L}_{\pm }(kX_{m_{0}})$ . The expressions for $\hat{L}_{\pm }$ are inserted into the equation for the main harmonic (5.12) and the result is divided by the real number $\int _{r_{\ast }}^{a}\unicode[STIX]{x1D6FC}X_{m_{0}}\,\text{d}r$ . The solution to this equation is still given by (5.14) with an appropriate rescaling of the constants $c_{1,2}$ (i.e. the term $H$ ). Hence this is equivalent to solve (5.12) imposing the integral condition (Wahlberg & Graves 2007):

(5.17) $$\begin{eqnarray}\displaystyle & \displaystyle \int _{r_{\ast }}^{a}\unicode[STIX]{x1D6FC}X_{m_{0}}\,\text{d}r=1. & \displaystyle\end{eqnarray}$$

Equation (5.17) provides the dispersion relation with $\hat{L}_{\pm }=(\mathbb{B}_{\pm }\mathbb{C}_{\pm }(1\pm m_{0})/(2a))/$ $[\mathbb{C}_{\pm }-\mathbb{B}_{\pm }-\mathbb{B}_{\pm }\mathbb{C}_{\pm }\unicode[STIX]{x1D6E5}/a]$ . Plugging (5.14) and (5.15) into (5.17), some straightforward algebra produces the following result:

(5.18) $$\begin{eqnarray}\displaystyle & \displaystyle {\displaystyle \frac{\hat{\unicode[STIX]{x1D6FE}}^{2}}{2}}+\left({\displaystyle \frac{\unicode[STIX]{x1D6FF}q}{q}}\right)^{2}=\left({\displaystyle \frac{\hat{\unicode[STIX]{x1D6FD}}}{2\unicode[STIX]{x1D700}}}\right)^{2}\left[{\displaystyle \frac{\hat{L}_{+}}{1+m_{0}}}+{\displaystyle \frac{\hat{L}_{-}}{1-m_{0}}}\right]A_{\unicode[STIX]{x1D6FE}}, & \displaystyle\end{eqnarray}$$

where $\hat{\unicode[STIX]{x1D6FD}}$ has been introduced in the previous section and $A_{\unicode[STIX]{x1D6FE}}=h/c$ $[(ch/\unicode[STIX]{x1D6FF}-\tanh ^{-1}(c\tanh (h/\unicode[STIX]{x1D6FF})))/(h/\unicode[STIX]{x1D6FF}-c\tanh ^{-1}(c\tanh (h/\unicode[STIX]{x1D6FF})))]$ . We shall now analyse more in detail the dispersion relation equation (5.18) and the associated eigenfunctions.

6 Analysis of the dispersion relation of the tanh model

The main difficulty arising in the solution of (5.18) is the dependence upon the growth rate embedded in the function $A_{\unicode[STIX]{x1D6FE}}$ , i.e. inside the coefficient $c$ . Generally speaking the solution of such a dispersion relation must be tackled numerically. Nonetheless an explicit dependence of $\unicode[STIX]{x1D6FE}$ on the plasma parameters can be obtained for some special cases.

Let us consider first the case of a perfectly conducting wall ( $\unicode[STIX]{x1D70E}\rightarrow \infty$ ). If the metallic wall is directly interfaced with the plasma then $\mathbb{B}_{\pm }\rightarrow \infty$ and because $\mathbb{C}_{\pm }>0$ we immediately have $\hat{L}_{\pm }\lessgtr 0$ . Therefore, since $A_{\unicode[STIX]{x1D6FE}}>0$ , the right-hand side of (5.18) is negative indicating stability. Thus with ideal wall boundary conditions, a vacuum region is necessary for the instability to develop. This also suggests a threshold in the wall position for the growth rate in accordance with the results presented in Zheng et al. (2013a , b ).

Let us now examine the behaviour of $A_{\unicode[STIX]{x1D6FE}}$ near marginal stability ( $c\ll 1$ ) and in the strong instability region ( $c\sim 1$ ). By taking the corresponding limits we obtain at leading order:

(6.1a,b ) $$\begin{eqnarray}\displaystyle & \displaystyle A_{\unicode[STIX]{x1D6FE}}\sim h\left(1-{\displaystyle \frac{\tanh (h/\unicode[STIX]{x1D6FF})}{h/\unicode[STIX]{x1D6FF}}}\right),\quad c\ll 1\;(2\unicode[STIX]{x1D6FF}q^{2}/(q\hat{\unicode[STIX]{x1D6FE}})^{2}\gg 1), & \displaystyle\end{eqnarray}$$
(6.2a,b ) $$\begin{eqnarray}\displaystyle & \displaystyle A_{\unicode[STIX]{x1D6FE}}\rightarrow h\left(1-{\displaystyle \frac{4h/\unicode[STIX]{x1D6FF}}{2h/\unicode[STIX]{x1D6FF}+\sinh (2h/\unicode[STIX]{x1D6FF})}}\right),\quad c\sim 1\;(2\unicode[STIX]{x1D6FF}q^{2}/(q\hat{\unicode[STIX]{x1D6FE}})^{2}\ll 1). & \displaystyle\end{eqnarray}$$

Thus the $\unicode[STIX]{x1D6FE}$ dependence in $A_{\unicode[STIX]{x1D6FE}}$ is removed and an explicit expression for the growth rate can be obtained. The resulting growth rates with respect to the plasma parameters are shown in figure 2.

Figure 2. Ideal wall boundary condition growth rates for the main mode $m_{0}/n=4/1$ with $q(r_{\ast })=m_{0}/n-\unicode[STIX]{x1D6FF}q$ , $\unicode[STIX]{x1D700}=1/3$ , $r_{s}=0.8$ , $\unicode[STIX]{x1D706}=2$ , $r_{\ast }/a=0.95$ and $h/\unicode[STIX]{x1D6FF}=3$ with respect to $\unicode[STIX]{x1D6FF}q$ (a), $\hat{\unicode[STIX]{x1D6FD}}$ (b) and wall position (c). Note that $\hat{\unicode[STIX]{x1D6FD}}$ can be moderately large due to the factor $q^{2}$ . The asymptotic behaviour obtained from formulae (6.1) and (6.2) is shown in (a).

It is worth noticing that particular attention has to be devoted to the expansion of $A_{\unicode[STIX]{x1D6FE}}$ near $c=1$ . Indeed a series expansion $A_{\unicode[STIX]{x1D6FE}}$ in $(c-1)$ cannot be performed if $h/\unicode[STIX]{x1D6FF}$ becomes sufficiently large. Therefore we decided to retain only the first term which provides the correct asymptotic behaviour for $h/\unicode[STIX]{x1D6FF}\rightarrow \infty$ . In such a case ( $\unicode[STIX]{x1D6FF}\rightarrow 0$ ) it is easy to see that both for $c\ll 1$ and $c\sim 1$ the step-like dispersion relation (cf. (5.10)) is formally recovered substituting $\unicode[STIX]{x1D6EC}^{(\pm )}\rightarrow \hat{L}_{\pm }$ and without the Mercier contribution (this has been dropped due to the ordering (2.15)). The shape of the eigenfunctions also qualitatively reduces to the step-like case as shown in figure 3. Differences remain due to the fact that in the system of equations (5.12) and (5.13) we dropped terms which conversely are retained in (5.2) and (5.3). In addition we point out that in the tanh model, contrarily to the step case, the shape of the eigenfunctions depends on the value of $\unicode[STIX]{x1D6FF}q$ (and so on the corresponding growth rate). We note that the radial structure of the main and sidebands harmonics closely resembles what has been found numerically in three-dimensional equilibria simulations (Cooper et al. 2015, 2016a , b ) and in MHD stability calculations (Medvedev et al. 2006; Zheng et al. 2013a , b ).

Figure 3. Plot of the eigenfunctions normalised with respect to the maximum of the main mode $m_{0}=4,n=1$ ( $X_{M}$ ) for the tanh (a)–(b) and the step-like model (c) with an ideal wall boundary condition with $\unicode[STIX]{x1D700}=1/3$ , $\hat{\unicode[STIX]{x1D6FD}}=3\,\%$ , $r_{s}=0.8$ , $\unicode[STIX]{x1D706}=2$ , $r_{\ast }/a=0.95$ and $q_{\ast }\approx 4$ . In (a) and (b) the eigenfunctions are computed for $\unicode[STIX]{x1D6FF}q=0.03$ corresponding to $\unicode[STIX]{x1D6FE}/\unicode[STIX]{x1D714}_{A}\approx 1.4\times 10^{-3}$ . The dashed vertical line indicates the position of the low-shear region middle point $r_{p}$ . We set in (a) $h/\unicode[STIX]{x1D6FF}=3$ while in (b) $h/\unicode[STIX]{x1D6FF}=20$ . Qualitatively the behaviour of the step model is recovered for large $h/\unicode[STIX]{x1D6FF}$ . Note that the maximum of the main harmonic ( $m_{0}=4$ ) in the tanh model is slightly shifted to the right with respect to $r_{p}$ .

We shall now consider the case with a resistive wall. For sake of simplicity in the following analysis we set $\unicode[STIX]{x1D714}_{A}=1$ . We assume that the conductivity is sufficiently large so that we can expand $\mathbb{B}_{\pm }$ in $1/\unicode[STIX]{x1D70E}$ . This yields $\mathbb{B}_{\pm }=\mathbb{B}_{\pm }|_{\unicode[STIX]{x1D70E}=\infty }+(\mathbb{B}_{\pm }^{w}/(S_{w}\unicode[STIX]{x1D6FE}(d/b)))$ with:

(6.3) $$\begin{eqnarray}\displaystyle & \displaystyle \mathbb{B}_{\pm }^{w}={\displaystyle \frac{4(m_{0}\pm 1)^{2}(a/b)^{2m_{0}\pm 2}}{1-(a/b)^{2m_{0}\pm 2}}}>0. & \displaystyle\end{eqnarray}$$

Inserting (6.3) into the expressions for $\hat{L}_{\pm }$ , to leading order we obtain:

(6.4) $$\begin{eqnarray}\displaystyle & \displaystyle {\displaystyle \frac{\hat{L}_{\pm }}{1\pm m_{0}}}\approx {\displaystyle \frac{\mathbb{B}_{\pm }^{0}\mathbb{C}_{\pm }/(2a)}{[\mathbb{C}_{\pm }-\mathbb{B}_{\pm }^{0}-\mathbb{B}_{\pm }^{0}\mathbb{C}_{\pm }\unicode[STIX]{x1D6E5}/a]}}+{\displaystyle \frac{1}{S_{w}\unicode[STIX]{x1D6FE}}}\left({\displaystyle \frac{(b/d)\mathbb{B}_{\pm }^{w}\mathbb{C}_{\pm }^{2}/(2a)}{[\mathbb{C}_{\pm }-\mathbb{B}_{\pm }^{0}-\mathbb{B}_{\pm }^{0}\mathbb{C}_{\pm }\unicode[STIX]{x1D6E5}/a]^{2}}}\right), & \displaystyle\end{eqnarray}$$

where $\mathbb{B}_{\pm }^{0}=\mathbb{B}_{\pm }|_{\unicode[STIX]{x1D70E}=\infty }$ . Let us formally write the expression above as $(\hat{L}_{\pm }/1\pm m_{0})\approx \mathscr{D}_{\pm }+\mathscr{R}_{\pm }/(S_{w}\unicode[STIX]{x1D6FE})$ . Hence the equation (5.18) can be cast in the following manner:

(6.5) $$\begin{eqnarray}\displaystyle & \displaystyle \hat{\unicode[STIX]{x1D6FE}}\left[\hat{\unicode[STIX]{x1D6FE}}^{2}-\hat{\unicode[STIX]{x1D6FE}}_{I}^{2}\right]=\unicode[STIX]{x1D6FE}_{w}^{3}, & \displaystyle\end{eqnarray}$$

where $\hat{\unicode[STIX]{x1D6FE}}_{I}^{2}=2(\hat{\unicode[STIX]{x1D6FD}}/2\unicode[STIX]{x1D700})^{2}h[\mathscr{D}_{+}+\mathscr{D}_{-}]$ is the growth rate with the ideal wall and:

(6.6) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6FE}_{w}^{3}=2{\displaystyle \frac{\sqrt{1+2q^{2}}}{nS_{w}}}\left({\displaystyle \frac{\hat{\unicode[STIX]{x1D6FD}}}{2\unicode[STIX]{x1D700}}}\right)^{2}h[\mathscr{R}_{+}+\mathscr{R}_{-}]>0, & \displaystyle\end{eqnarray}$$

having assumed for sake of simplicity that $h/\unicode[STIX]{x1D6FF}\gg 1$ , i.e. $A_{\unicode[STIX]{x1D6FE}}\rightarrow h$ .

Suppose that $\hat{\unicode[STIX]{x1D6FE}}_{I}^{2}>0$ . Therefore if $\unicode[STIX]{x1D6FE}_{w}/\hat{\unicode[STIX]{x1D6FE}}_{I}\ll 1$ , the growth rate can be written as $\hat{\unicode[STIX]{x1D6FE}}=\hat{\unicode[STIX]{x1D6FE}}_{I}+(\unicode[STIX]{x1D6FE}_{w}^{3}/2\hat{\unicode[STIX]{x1D6FE}}_{I}^{2})$ . Conversely if $\unicode[STIX]{x1D6FE}_{w}/\hat{\unicode[STIX]{x1D6FE}}_{I}\gg 1$ , the growth rate is $\hat{\unicode[STIX]{x1D6FE}}=\unicode[STIX]{x1D6FE}_{w}+(\hat{\unicode[STIX]{x1D6FE}}_{I}^{2}/3\unicode[STIX]{x1D6FE}_{w})$ . In both cases the destabilising effect of the resistive wall is evident. In the opposite case, when $\hat{\unicode[STIX]{x1D6FE}}_{I}^{2}<0$ and $\unicode[STIX]{x1D6FE}_{w}$ are sufficiently small (i.e. we are in the stability region with the ideal wall) instability is still possible with growth rate $\hat{\unicode[STIX]{x1D6FE}}=\unicode[STIX]{x1D6FE}_{w}^{3}/|\hat{\unicode[STIX]{x1D6FE}}_{I}^{2}|$ . The behaviour of $\unicode[STIX]{x1D6FE}_{w}$ with respect to $n$ is shown in figure 4.

Figure 4. Plot of $\unicode[STIX]{x1D6FE}_{w}/\unicode[STIX]{x1D714}_{A}$ with respect to the toroidal mode number $n$ with $\unicode[STIX]{x1D700}=1/3$ , $b/a=1.3$ , $d/a=10^{-2}$ $\hat{\unicode[STIX]{x1D6FD}}=5\,\%$ , $r_{s}=0.8$ , $\unicode[STIX]{x1D706}=2$ , $r_{\ast }/a=0.95$ , $S_{w}=10^{5}$ and $q_{\ast }=4$ (here we have dropped the normalisation $\unicode[STIX]{x1D714}_{A}=1$ ). Note that $\unicode[STIX]{x1D6FE}_{w}/\unicode[STIX]{x1D714}_{A}$ peaks for small values of $n$ .

Although toroidal rotation effects are not considered in the present work, based on the results shown in Brunetti et al. (2018) we expect that the rotation induced Doppler shift should be included in the low-shear inertial contribution with a modification of the dispersion relation ( $\unicode[STIX]{x1D6FA}$ is the rotation frequency evaluated at $r_{\ast }$ ):

(6.7) $$\begin{eqnarray}\displaystyle & \displaystyle \hat{\unicode[STIX]{x1D6FE}}\left[{\displaystyle \frac{(\unicode[STIX]{x1D6FE}-\text{i}n\unicode[STIX]{x1D6FA})^{2}(1+2q^{2})}{n^{2}}}-\hat{\unicode[STIX]{x1D6FE}}_{I}^{2}\right]=\unicode[STIX]{x1D6FE}_{w}^{3}. & \displaystyle\end{eqnarray}$$

It is easy to see that for sufficiently large toroidal rotation, the resistive wall effects become negligible.

7 Conclusions

In this work we addressed the linear stability properties of infernal type perturbations with a safety factor flat and close to a rational number ( $q\approx m_{0}/n$ ) near the edge. This perturbation presents a main ( $m_{0},n$ ) Fourier harmonic, resonant in the flat $q$ region, which is coupled with its neighbouring sidebands ( $m_{0}\pm 1,n$ ). Although it is well known that plasma resistivity can grow to large values near the edge (Turnbull et al. 2016), the standard ideal MHD model has been adopted. This is motivated by the fact that the $m$ values resonating with the edge $q$ must be moderately large, and thus they are expected to be stable against standard tearing perturbations (Hegna & Callen 1994).

An exact analytic treatment of the magnetic perturbation in the inner region of large magnetic shear has been possible with a careful choice, although rather general, of the safety factor. In the region of the local edge flattening of the safety factor, the equilibrium profiles for pressure and mass density have been modelled with two simple classes of analytic functions. In one case (already discussed in Brunetti et al. (2018)) the step-like model has been employed. In the second case we introduced a layer ordering in the low-shear region resulting in a simplification of the coupled equations. This allowed us to employ a more realistic tanh model. A vacuum region with a parabolic safety factor ( $q\sim r^{2}$ ) has been included with either ideally conducting or resistive wall boundary conditions (note that in the ideally conducting wall case, the presence of a vacuum region separating plasma and wall is necessary for the instability to develop). Each region has been analysed separately and the solutions eventually matched across the transitions points. The advantage of this choice for the equilibrium profiles, i.e. safety factor and either step or tanh-like pressure and mass density, consists in an exact mathematical description of both the eigenfunctions and the associated dispersion relation. We point out that additional corrections such as inertia and residual coupling associated with the behaviour of the sidebands inside the sheared region have been neglected (a brief description of their effects is given in the appendices). The dispersion relation which has been derived is sufficiently simple to show the explicit dependence of the growth rate upon the engineering plasma parameters, i.e. plasma current profile and $\unicode[STIX]{x1D6FD}$ .

Similarly to the peeling–ballooning (PB) modes, the drive of the instability is the combined effect of pressure gradients and the edge current (field line bending weakening) which translates in a toroidal coupling of neighbouring poloidal Fourier harmonics. The structure of edge infernal modes resembles the one of PB perturbations (Snyder et al. 2004), consisting of a non-vanishing peeling component accompanied by an inner located bell-shaped displacement. Hence we can regard these perturbations as very low- $n$ PB modes. Such a model seems to be appropriate for the description of the edge harmonic oscillations. Indeed the shape of the eigenfunctions calculated analytically exhibits similarities with both numerical results (Medvedev et al. 2006; Zheng et al. 2013a , b ) and with experimental findings (Chen et al. 2016) (in particular the bell shape and the radial localisation are correctly recovered by the model we proposed). We also note that experimental measurement of density and temperature fluctuations show even parity (in the radial direction) (Burrell et al. 2001; Greenfield et al. 2001), excluding the tearing character of the perturbation and thus supporting even more the ideal plasma approximation.

It is worth emphasising the facts that this work is based on very simple ‘model’ profiles and is developed in the linear framework. Nevertheless, although the observed perturbations are in their nonlinear stage, we shall expect that the linear structural properties (shape, parity, location of (main) peak) may persist from the linear into the nonlinear stage.

It is found that in the case of ideal wall boundary conditions the growth rates become larger with increasing $n$  (Zheng et al. 2013b ; Liu et al. 2015). This might be relevant for QH regimes with zero injected torque (Burrell et al. 2016; Chen et al. 2017). In such a case a weakly quasi-coherent MHD (broadband) activity with a cessation of the EHO is observed. This could be linked with the nonlinear interplay of different $n$ values (also moderately large) which are all allowed to grow in the linear phase (Liu et al. 2015, 2018). Consequently in the nonlinear stage they can compete with each other without favouring a particular $n$ to be dominant (Liu et al. 2018). The inclusion of a finite wall conductivity has a destabilising effect and near the ideal wall marginal stability boundary the resistive wall growth rate peaks for small values of the toroidal (poloidal) mode numbers.

We point out that from the results obtained in Brunetti et al. (2018) with step-like profiles and an ideal wall, a subsonic toroidal rotation Doppler shifts the eigenfrequency. Although in the present work effects of plasma rotation have been neglected, we may argue that also with the tanh model a purely toroidal rotation should enter through a Doppler shift of the eigenfrequencies. This predicts a rotating mode with frequency $n\unicode[STIX]{x1D6FA}_{\ast }$ ( $\unicode[STIX]{x1D6FA}_{\ast }$ is the toroidal rotation at the pedestal top) in accordance with experimental data (see e.g. Burrell et al. 2001; Greenfield et al. 2001; Solano et al. 2010). We may also expect resistive wall effects to become negligible for sufficiently large toroidal flows (Nave & Wesson 1990).

In addition, as shown in Brunetti et al. (2018) in the ideal wall case, the mode stability properties are not altered by a subsonic toroidal flow. Therefore we shall infer that the toroidal rotation alone is not sufficient to explain why in the experiments only a single $n$ dominant mode is observed (see e.g. (Chen et al. 2016)). Moreover experimental findings presented in Garofalo et al. (2011), Burrell et al. (2012) show that the key ingredient for accessing the QH mode is the shear associated with the $\boldsymbol{E}\times \boldsymbol{B}$ drift. This induces us to conclude that diamagnetic drifts and sheared poloidal flows with the associated electric field effects should play a fundamental role in the accessibility of the QH operational space. This thus will require further analytic investigations.

Finally we argue that the effects of the separatrix may alter the jump conditions of the magnetic perturbation, resulting in a modification of the dispersion relation. Shaping effects also (mainly elongation) could alter the mode dynamics. It is envisaged that a numerical approach is required to address such a difficult problem.

Appendix A. Effects of inertial corrections on the lower sideband

The quantity $B_{2}/A_{2}$ which appears in (3.5), depends on the growth rate $\unicode[STIX]{x1D6FE}$ as inertial corrections may become important when the resonant surface of the lower sideband is approached. We expand (3.4) and (3.5) about $r_{s}$ giving (we use the same notation as in Mikhailovskii (1998), Brunetti, Lazzaro & Nowak (2017)):

(A 1a,b ) $$\begin{eqnarray}\displaystyle & \displaystyle X_{-}(z<1)\sim {\displaystyle \frac{1}{1-z}}+\unicode[STIX]{x1D702}\unicode[STIX]{x1D701}\ln (1-z)+\unicode[STIX]{x1D6E5}_{c},\quad X_{-}(z>1)\sim {\displaystyle \frac{1}{z-1}}-\unicode[STIX]{x1D702}\unicode[STIX]{x1D701}\ln (z-1)+\unicode[STIX]{x1D6E5}_{p}, & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

where $\unicode[STIX]{x1D6E5}_{c}=[\unicode[STIX]{x1D702}\unicode[STIX]{x1D701}(2\unicode[STIX]{x1D6FE}_{E}-1+\unicode[STIX]{x1D6F9}(\unicode[STIX]{x1D701}+1)+\unicode[STIX]{x1D6F9}(\unicode[STIX]{x1D702}+1))-m_{-}-1/\unicode[STIX]{x1D706}]$ ( $\unicode[STIX]{x1D6FE}_{E}$ is the Euler–Mascheroni constant and $\unicode[STIX]{x1D6F9}$ the digamma function (see Abramowitz & Stegun 1968, p. 253)), $\unicode[STIX]{x1D6E5}_{p}=d_{+}+((B_{{>}}/A_{{>}})/(1+B_{{>}}/A_{{>}}))(d_{-}-d_{+})$ and:

(A 2) $$\begin{eqnarray}\displaystyle & \displaystyle d_{\pm }=\unicode[STIX]{x1D702}\unicode[STIX]{x1D701}[1-2\unicode[STIX]{x1D6FE}_{E}-\unicode[STIX]{x1D6F9}(\pm \unicode[STIX]{x1D701}+1)-\unicode[STIX]{x1D6F9}(\mp \unicode[STIX]{x1D702}+1)]\mp {\displaystyle \frac{\bar{m}_{-}\pm 1}{\unicode[STIX]{x1D706}}}. & \displaystyle\end{eqnarray}$$

Equations (A 1) contain a logarithmic term which is generally neglected when the inertial layer is studied in slab geometry. Following the analysis presented in Mikhailovskii (1998), we allow for inertial effects in (2.7) close to $r_{s}$ so that ( $^{\prime }\equiv \text{d}/\text{d}x$ ):

(A 3) $$\begin{eqnarray}\displaystyle & \displaystyle [(1+x)^{2/\unicode[STIX]{x1D706}+1}(x^{2}+\bar{\unicode[STIX]{x1D6FE}}^{2})X^{\prime }]^{\prime }-(1+x)^{2/\unicode[STIX]{x1D706}-1}{\displaystyle \frac{m_{-}^{2}-1}{\unicode[STIX]{x1D706}^{2}}}[x^{2}+\bar{\unicode[STIX]{x1D6FE}}^{2}]X=0, & \displaystyle\end{eqnarray}$$

where $z=1+x$ and $\bar{\unicode[STIX]{x1D6FE}}^{2}=\unicode[STIX]{x1D6FE}^{2}/(S\unicode[STIX]{x1D714}_{A})^{2}(1+2q^{2})$ . Let us introduce the following ordering: $x\sim \bar{\unicode[STIX]{x1D6FE}}\sim \unicode[STIX]{x1D6FF},\text{d}/\text{d}x\sim 1/\unicode[STIX]{x1D6FF}$ with $\unicode[STIX]{x1D6FF}\ll 1$ and $(m_{-}^{2}-1/\unicode[STIX]{x1D706}^{2})\sim S\sim 1$ . We write the solution as $X=X_{0}+\unicode[STIX]{x1D6FF}X_{1}+\text{O}(\unicode[STIX]{x1D6FF}^{2})$ and then expand the equation above in a series of $\unicode[STIX]{x1D6FF}$ . Eventually we have:

(A 4) $$\begin{eqnarray}\displaystyle & \displaystyle X=c_{0}-c_{1}/\bar{\unicode[STIX]{x1D6FE}}\arctan (x/\bar{\unicode[STIX]{x1D6FE}})\pm c_{3}\ln (x^{2}+\bar{\unicode[STIX]{x1D6FE}}^{2}), & \displaystyle\end{eqnarray}$$

where $\pm$ is for $x\lessgtr 0$ . As $x/\bar{\unicode[STIX]{x1D6FE}}\rightarrow \infty$ the logarithmic terms in (A 4) are immediately matched with (A 1) by choosing $c_{3}=\unicode[STIX]{x1D702}\unicode[STIX]{x1D701}/2$ , while the non-logarithmic part of (A 4) behaves as:

(A 5) $$\begin{eqnarray}\displaystyle & \displaystyle X\sim c_{0}\pm {\displaystyle \frac{c_{1}}{|x|}}\left(1-{\displaystyle \frac{\unicode[STIX]{x03C0}}{2}}|x|/\bar{\unicode[STIX]{x1D6FE}}\right), & \displaystyle\end{eqnarray}$$

where $\pm$ is for $x\lessgtr 0$ as before. Matching (A 5) with (A 1) gives $\unicode[STIX]{x03C0}/\hat{\unicode[STIX]{x1D6FE}}+\unicode[STIX]{x1D6E5}_{c}+\unicode[STIX]{x1D6E5}_{p}=0$ and eventually we obtain:

(A 6) $$\begin{eqnarray}\displaystyle & \displaystyle B_{{>}}/A_{{>}}=-{\displaystyle \frac{1-[\unicode[STIX]{x1D702}\unicode[STIX]{x1D701}\cot (\unicode[STIX]{x03C0}\unicode[STIX]{x1D702})]\bar{\unicode[STIX]{x1D6FE}}}{1-[\unicode[STIX]{x1D702}\unicode[STIX]{x1D701}\cot (\unicode[STIX]{x03C0}\unicode[STIX]{x1D701})]\bar{\unicode[STIX]{x1D6FE}}}}\approx -1+\bar{\unicode[STIX]{x1D6FE}}\unicode[STIX]{x1D702}\unicode[STIX]{x1D701}[\cot (\unicode[STIX]{x03C0}\unicode[STIX]{x1D702})-\cot (\unicode[STIX]{x03C0}\unicode[STIX]{x1D701})], & \displaystyle\end{eqnarray}$$

where the last approximation holds if $\hat{\unicode[STIX]{x1D6FE}}\ll 1$ . For sufficiently small growth rates we have $B_{{>}}/A_{{>}}=-1$ .

Appendix B. Allowance for residual coupling effects at the resonant layer of the lower sideband

We assume that a small coupling contribution remains near the resonant point $r_{s}$ of the lower harmonic $X_{-}$ . This corresponds to the fact that the main harmonic $m_{0}$ and pressure gradients can be non-exactly vanishing for $r<r_{\ast }$ . For sake of simplicity we write $m=m_{0}-1$ . Thus we modify the (2.7) by adding the constant $U$ which represents the coupling corrections similar to the ones appearing in (2.8) and (2.9):

(B 1) $$\begin{eqnarray}\displaystyle & \displaystyle {\displaystyle \frac{\text{d}}{\text{d}r}}\left[r^{3}Q{\displaystyle \frac{\text{d}X_{-}}{\text{d}r}}\right]-r(m^{2}-1)QX_{-}+U=0, & \displaystyle\end{eqnarray}$$

where here $Q=k_{||}^{2}+\unicode[STIX]{x1D6FE}^{2}(1+2q^{2})\unicode[STIX]{x1D714}_{A}^{2}$ . We assume that $U$ is small far from $r_{s}$ , but sufficiently large close to this point.

To simplify the analysis we employ a safety factor which is constant $1/q=\unicode[STIX]{x1D707}_{0}$ for $0<r<r_{0}$ and $1/q=\unicode[STIX]{x1D707}_{\ast }\sim n/m_{0}<\unicode[STIX]{x1D707}_{0}$ for $r_{\ast }<r<a$ , while in the region $r_{0}<r<r_{\ast }$ behaves as (Brunetti et al. 2018):

(B 2) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D707}\approx n/m_{0}(r_{\ast }/r)^{2}. & \displaystyle\end{eqnarray}$$

This corresponds to a current profile of the form (the profiles for $q$ and $J_{0}$ are shown in figure 5):

(B 3) $$\begin{eqnarray}\displaystyle \displaystyle R_{0}J_{0}=\left\{\begin{array}{@{}ll@{}}2\unicode[STIX]{x1D707}_{0}, & 0<r<r_{0}\\ 0, & r_{0}<r<r_{\ast },\\ 2\unicode[STIX]{x1D707}_{\ast }, & r_{\ast }<r<a,\end{array}\right. & & \displaystyle\end{eqnarray}$$

Figure 5. Step current model: the $q$ profile is assumed parabolic (i.e. vanishing current) for $r_{0}<r<r_{\ast }$ .

In the region $r_{0}<r<r_{\ast }$ , sufficiently far from $r_{s}$ we neglect the terms proportional to $\unicode[STIX]{x1D6FE}$ , hence the solution of (B 1) reads:

(B 4) $$\begin{eqnarray}\displaystyle & \displaystyle X_{-}=T_{\lessgtr }{\displaystyle \frac{(r/r_{s})^{m-1}}{k_{||}}}+K_{\lessgtr }{\displaystyle \frac{(r/r_{s})^{-m-1}}{k_{||}}}-{\displaystyle \frac{U}{2m}}\left[{\displaystyle \frac{r^{m-1}}{k_{||}}}\int {\displaystyle \frac{r^{-m-1}}{k_{||}}}\,\text{d}r-{\displaystyle \frac{r^{-m-1}}{k_{||}}}\int {\displaystyle \frac{r^{m-1}}{k_{||}}}\,\text{d}r\right], & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

where $K_{\lessgtr }$ and $T_{\lessgtr }$ are constants and the integrals in the expression above can be easily represented in terms of the Gauss hypergeometric function. Note that a more careful treatment is needed when $m$ is even.

When the resonant surface is approached we have:

(B 5) $$\begin{eqnarray}\displaystyle \displaystyle X_{-}(r\lessgtr r_{s}) & = & \displaystyle \pm {\displaystyle \frac{(K_{\lessgtr }+T_{\lessgtr })}{2n|\hat{x}|}}-{\displaystyle \frac{K_{\lessgtr }(1/2-m)+T_{\lessgtr }(1/2+m)}{2n}}-{\displaystyle \frac{U}{4n^{2}r_{s}}}\ln |\hat{x}|\nonumber\\ \displaystyle \displaystyle & & \displaystyle \mp {\displaystyle \frac{U}{4m^{2}n^{2}r_{s}|\hat{x}|}}-{\displaystyle \frac{U}{4n^{2}r_{s}}}\left(\unicode[STIX]{x1D6FE}_{E}-1+\ln 2+\unicode[STIX]{x1D6F9}\left({\displaystyle \frac{m}{2}}\right)+{\displaystyle \frac{1}{m}}-{\displaystyle \frac{1}{2m^{2}}}\right),\end{eqnarray}$$

where $\hat{x}=(r-r_{s})/r_{s}$ . Far from $r_{s}$ , neglecting the term proportional to $U$ , $K_{{<}}$ is expressed in terms of $T_{{<}}$ by imposing smoothness at $r_{0}$ (if $r_{0}\ll r_{\ast }$ and $m_{0}$ sufficiently large we can approximate $K_{{<}}\approx 0$ ). If the inertial layer is thin, we may assume that the poloidal flux $\unicode[STIX]{x1D713}_{-}\sim \hat{x}X_{-}$ is continuous across $r_{s}$ . This implies that $(K_{{<}}+T_{{<}})=(K_{{>}}+T_{{>}})$ . The part of (B 5) which does not contain the logarithm can be written as:

(B 6) $$\begin{eqnarray}\displaystyle & \displaystyle X_{-}(r\lessgtr r_{s})\sim {\displaystyle \frac{1}{|\hat{x}|}}\left(1+\unicode[STIX]{x1D6E5}_{\mp }|\hat{x}|\right), & \displaystyle\end{eqnarray}$$

where we defined:

(B 7) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6E5}_{\pm }=\pm {\displaystyle \frac{{\displaystyle \frac{1}{2}}+m+\left({\displaystyle \frac{1/2-m}{1/2+m}}\right){\displaystyle \frac{K_{\lessgtr }}{T_{\lessgtr }}}+{\displaystyle \frac{U/T_{\lessgtr }}{2m^{2}nr_{s}}}\mathscr{H}}{1+{\displaystyle \frac{K_{\lessgtr }}{T_{\lessgtr }}}-{\displaystyle \frac{U/T_{\lessgtr }}{2m^{2}nr_{s}}}}}, & \displaystyle\end{eqnarray}$$

with $\mathscr{H}=m^{2}[\unicode[STIX]{x1D6FE}_{E}-1+\ln 2+\unicode[STIX]{x1D6F9}(m/2)]+m-1/2$ .

In the inertial layer we approximate (B 1) as:

(B 8) $$\begin{eqnarray}\displaystyle & \displaystyle {\displaystyle \frac{\text{d}}{\text{d}\hat{x}}}\left[(\hat{x}^{2}+\bar{\unicode[STIX]{x1D6FE}})^{2}{\displaystyle \frac{\text{d}X_{-}}{\text{d}\hat{x}}}\right]+U_{1}=0, & \displaystyle\end{eqnarray}$$

where $\bar{\unicode[STIX]{x1D6FE}}^{2}=\unicode[STIX]{x1D6FE}^{2}(1+2q^{2})/(\unicode[STIX]{x1D714}_{A}mr_{s}\unicode[STIX]{x1D707}_{s}^{\prime })^{2}$ and $U_{1}=(1/r_{s})U/(mr_{s}\unicode[STIX]{x1D707}_{s}^{\prime })^{2}$ . The solution of (B 8) is $X_{-}=c_{0}+c_{1}\arctan (\hat{x}/\bar{\unicode[STIX]{x1D6FE}})-(U_{1}/2)\ln (\hat{x}^{2}+\bar{\unicode[STIX]{x1D6FE}}^{2})$ , where $d_{0,1}$ are constants of integration. The logarithmic terms in (B 5) and in the solution of (B 8) are automatically matched. If $\hat{x}/\bar{\unicode[STIX]{x1D6FE}}\rightarrow \infty$ , the asymptotic behaviour of the non-logarithmic part of the layer solution behaves as (A 5). Matching gives $\unicode[STIX]{x03C0}/\hat{\unicode[STIX]{x1D6FE}}+\unicode[STIX]{x1D6E5}_{+}+\unicode[STIX]{x1D6E5}_{-}=0$ (see previous section).


Abramowitz, M. & Stegun, I. 1968 Handbook of Mathematical Functions With Formulas, Graphs, and Mathematical Tables. Dover.
Bernstein, I. B., Frieman, E. A., Kruskal, M. D. & Kulsrud, R. M. 1958 An energy principle for hydromagnetic stability problems. Proc. R. Soc. Lond. A 244, 17.
Brunetti, D., Graves, J. P., Cooper, W. A. & Wahlberg, C. 2014 Fast growing resistive two fluid instabilities in hybrid-like tokamak configuration. Plasma Phys. Control. Fusion 56, 075025.
Brunetti, D., Graves, J. P., Lazzaro, E., Mariani, A., Nowak, S., Cooper, W. A. & Wahlberg, C. 2018 Analytic stability criteria for edge MHD oscillations in high performance ELM free tokamak regimes. Nucl. Fusion 58, 014002.
Brunetti, D., Lazzaro, E. & Nowak, S. 2017 Ideal and resistive magnetohydrodynamic instabilities in cylindrical geometry with a sheared flow along the axis. Plasma Phys. Control. Fusion 59, 055012.
Burrell, K. H., Austin, M. E., Brennan, D. P., DeBoo, J. C., Doyle, E. J., Fenzi, C., Fuchs, C., Gohil, P., Greenfield, C. M., Groebner, R. J. et al. 2001 Quiescent double barrier high-confinement mode plasmas in the DIII-D tokamak. Phys. Plasmas 8, 2153.
Burrell, K. H., Austin, M. E., Brennan, D. P., DeBoo, J. C., Doyle, E. J., Gohil, P., Greenfield, C. M., Groebner, R. J., Lao, L. L., Luce, T. C. et al. 2002 Quiescent H-mode plasmas in the DIII-D tokamak. Plasma Phys. Control. Fusion 44, A253.
Burrell, K. H., Barada, K., Chen, X., Garofalo, A. M., Groebner, R. J., Muscatello, C. M., Osborne, T. H., Petty, C. C., Rhodes, T. L., Snyder, P. B. et al. 2016 Discovery of stationary operation of quiescent H-mode plasmas with net-zero neutral beam injection torque and high energy confinement on DIII-D. Phys. Plasmas 23, 056103.
Burrell, K. H., Garofalo, A. M., Solomon, W. M., Fenstermacher, M. E., Osborne, T. H., Park, J.-K., Schaffer, M. J. & Snyder, P. B. 2012 Reactor-relevant quiescent H-mode operation using torque from non-axisymmetric, non-resonant magnetic fields. Phys. Plasmas 19, 056117.
Burrell, K. H., West, W. P., Doyle, E. J., Austin, M. E., Casper, T. A., Gohil, P., Greenfield, C. M., Groebner, R. J., Hyatt, A. W., Jayakumar, R. J. et al. 2005 Advances in understanding quiescent H-mode plasmas in DIII-D. Phys. Plasmas 12, 056121.
Chen, X., Burrell, K., Ferraro, N. M., Osborne, T. H., Austin, M. E., Garofalo, A. M., Groebner, R. J., Kramer, G. J. Jr, Luhmann, N. C., McKee, G. R. et al. 2016 Rotational shear effects on edge harmonic oscillations in DIII-D quiescent H-mode discharges. Nucl. Fusion 56, 076011.
Chen, X., Burrell, K. H., Osborne, T. H., Solomon, W. M., Barada, K., Garofalo, A. M., Groebner, R. J., Luhmann, N. C., McKee, G. R., Muscatello, C. M. et al. 2017 Stationary QH-mode plasmas with high and wide pedestal at low rotation on DIII-D. Nucl. Fusion 57, 022007.
Connor, J. W., Ham, C. J. & Hastie, R. J. 2016 The effect of plasma beta on high-n ballooning stability at low magnetic shear. Plasma Phys. Control. Fusion 58, 085002.
Connor, J. W., Hastie, R. J., Wilson, H. R. & Miller, R. L. 1998 Magnetohydrodynamic stability of tokamak edge plasmas. Phys. Plasmas 7, 2687.
Cooper, W. A., Brunetti, D., Duval, B. P., Faustin, J. M., Graves, J. P., Kleiner, A., Patten, H., Pfefferlé, D., Porte, L., Raghunathan, M. et al. 2016a Saturated ideal kink/peeling formations described as three-dimensional magnetohydrodynamic tokamak equilibrium states. Phys. Plasmas 23, 040701.
Cooper, W. A., Graves, J. P., Duval, B. P., Porte, L., Reimerdes, H., Sauter, O. & Tran, T.-M. 2015 A 3-D MHD equilibrium description of nonlinearly saturated ideal external kink/peeling structures in tokamaks. J. Plasma Phys. 81, 515810605.
Cooper, W. A., Graves, J. P., Duval, B. P., Sauter, O., Faustin, J. M., Kleiner, A., Lanthaler, S., Patten, H., Raghunathan, M., Tran, T.-M. et al. 2016b Three-dimensional magnetohydrodynamic equilibrium of quiescent H-modes in tokamak systems. Plasma Phys. Control. Fusion 58, 064002.
D’haeseleer, W. D., Hitchon, W. G., Callen, J. D. & Shohet, J. L. 1991 Flux Coordinates and Magnetic Field Structure. Springer.
Freidberg, J. P. 1987 Ideal Magnetohydrodynamics. Plenum Press.
Garofalo, A. M., Solomon, W. M., Park, J.-K., Burrell, K. H., DeBoo, J. C., Lanctot, M., McKee, G. R., Reimerdes, H., Schmitz, L., Schaffer, M. et al. 2011 Advances towards QH-mode viability for ELM-stable operation in ITER. Nucl. Fusion 51, 083018.
Gimblett, C. G. 1986 On free boundary instabilities induced by a resistive wall. Nucl. Fusion 26, 617.
Gimblett, C. G., Hastie, R. J. & Hender, T. C. 1996 An analytic study of the magnetohydrodynamic stability of inverse shear profiles. Phys. Plasmas 3, 3369.
Greene, J. M. & Chance, M. S. 1981 The second region of stability against ballooning modes. Nucl. Fusion 21, 453.
Greene, J. M. & Weimer, J. L. J. K. E. 1971 Tokamak equilibrium. Phys. Fluids 14, 671.
Greenfield, C. M., Burrell, K. H., DeBoo, J. C., Doyle, E. J., Stallard, B. W., Synakowski, E. J., Fenzi, C., Gohil, P., Groebner, R. J., Lao, L. L. et al. 2001 Quiescent double barrier regime in the DIII-D Tokamak. Phys. Rev. Lett. 86, 4544.
Hastie, R. J. & Hender, T. C. 1988 Toroidal internal kink stability in tokamaks with ultra flat q profiles. Nucl. Fusion 28, 585.
Hazeltine, R. D. & Meiss, J. D. 1992 Plasma Confinement. Addison-Wesley Publishing Company.
Hegna, C. C. & Callen, J. D. 1994 Stability of tearing modes in tokamak plasmas. Phys. Plasmas 1, 2308.
Kuvshinov, B. N. & Mikhailovskii, A. B. 1991 Analytic theory of internal kink and tearing modes in a tokamak. Sov. J. Plasma Phys. 16, 639.
Lashmore-Davies, C. N. 2001 The resistive wall instability and critical flow velocity. Phys. Plasmas 8, 151.
Lazzaro, E. & Zanca, P. 2003 Effects of a nonresonant $m=0$ , $n=1$ perturbation driven by coupling to the $m=2$ , $n=1$ mode in an elongated tokamak. Phys. Plasmas 10, 2399.
Liu, F., Huijsmans, G. T. A., Loarte, A., Garofalo, A. M., Solomon, W. M., Hoelzl, M., Nkonga, B., Pamela, S., Becoulet, M., Orain, F. et al. 2018 Nonlinear MHD simulations of QH-mode DIII-D plasmas and implications for ITER high Q scenarios. Plasma Phys. Control. Fusion 60, 014039.
Liu, F., Huijsmans, G. T. A., Loarte, A., Garofalo, A. M., Solomon, W. M., Snyder, P. B., Hoelzl, M. & Zeng, L. 2015 Nonlinear MHD simulations of Quiescent H-mode plasmas in DIII-D. Nucl. Fusion 55, 113002.
Medvedev, S. Y., Martynov, A. A., Martin, Y. R., Sauter, O. & Villard, L. 2006 Edge kink/ballooning mode stability in tokamaks with separatrix. Plasma Phys. Control. Fusion 48, 927.
Mikhailovskii, A. B. 1998 Instabilities in a Confined Plasma. IOP.
Nave, M. F. F. & Wesson, J. A. 1990 Mode locking in tokamaks. Nucl. Fusion 30, 2575.
Snyder, P. B., Wilson, H. R., Ferron, J. R., Lao, L. L., Leonard, A. W., Mossessian, D., Murakami, M., Osborne, T. H., Turnbull, A. D. & Xu, X. Q. 2004 ELMs and constraints on the H-mode pedestal: peeling-ballooning stability calculation and comparison with experiment. Nucl. Fusion 44, 320.
Solano, E. R., Lomas, P. J., Alper, B., Xu, G. S., Andrew, Y., Arnoux, G., Boboc, A., Barrera, L., Belo, P., Beurskens, M. N. A. et al. 2010 Observation of confined current ribbon in JET plasmas. Phys. Rev. Lett. 104, 185003.
Suttrop, W., Conway, G. D., Fattorini, L., Horton, L. D., Kurki-Suonio, T., Maggi, C. F., Maraschek, M., Meister, H., Neu, R., Putterich, T. et al. 2004 Study of quiescent H-mode plasmas in ASDEX Upgrade. Plasma Phys. Control. Fusion 46, A151.
Suttrop, W., Hynönen, V., Kurki-Suonio, T., Lang, P., Maraschek, M., Neu, R., Stäbler, A., Conway, G., Hacquin, S., Kempenaars, M. et al. 2005 Studies of the ‘Quiescent H-mode’ regime in ASDEX Upgrade and JET. Nucl. Fusion 45, 721.
Suttrop, W., Maraschek, M., Conway, G. D., Fahrbach, H.-U., Haas, G., Horton, L. D., Kurki-Suonio, T., Lasnier, C. J., Leonard, A. W., Maggi, C. F. et al. 2003 ELM-free stationary H-mode plasmas in the ASDEX Upgrade tokamak. Plasma Phys. Control. Fusion 45, 1399.
Turnbull, A. D., Hanson, J. M., Turco, F., Ferraro, N. M., Lanctot, M. J., Lao, L. L., Strait, E. J., Piovesan, P. & Martin, P. 2016 The external kink mode in diverted tokamaks. J. Plasma Phys. 82, 515820301.
Turnbull, A. D. & Troyon, F. 1989 Toroidal effects on current driven modes in tokamaks. Nucl. Fusion 29, 1887.
Waelbroeck, F. L. & Hazeltine, R. D. 1988 Stability of low-shear tokamaks. Phys. Fluids 31, 1217.
Wahlberg, C. & Graves, J. P. 2007 Stability analysis of internal ideal modes in low-shear tokamaks. Phys. Plasmas 14, 110703.
Wahlberg, C., Graves, J. P. & Chapman, I. T. 2013 Analysis of global hydromagnetic instabilities driven by strongly sheared toroidal flows in tokamak plasmas. Plasma Phys. Control. Fusion 55, 105004.
Wesson, J. A. 1978 Hydromagnetic stability of tokamaks. Nucl. Fusion 18, 87.
Wesson, J. A. 2011 Tokamaks. Oxford University Press.
Zheng, L. J., Kotschenreuther, M. T. & Valanju, P. 2013a Behavior of $n=1$ magnetohydrodynamic modes of infernal type at high-mode pedestal with plasma rotation. Phys. Plasmas 20, 012501.
Zheng, L. J., Kotschenreuther, M. T. & Valanju, P. 2013b Low- $n$ magnetohydrodynamic edge instabilities in quiescent H-mode plasmas with a safety- factor plateau. Nucl. Fusion 53, 063009.