Theory of Edge Localized Mode Suppression by Static Resonant Magnetic Perturbations in the DIII-D Tokamak

TL;DR

This paper develops an analytic model to quickly simulate and understand how static resonant magnetic perturbations suppress edge localized modes in the DIII-D tokamak, building on previous nonlinear MHD simulations.

Contribution

It simplifies existing nonlinear MHD models into an analytic form enabling rapid simulation of RMP-induced ELM suppression in DIII-D.

Findings

Analytic model reproduces key features of ELM suppression

Model predicts ELM suppression thresholds accurately

Provides physical insight into RMP penetration mechanism

Abstract

The plasma response to an externally generated, static, n=2, resonant magnetic perturbation (RMP) in the pedestal region of DIII-D discharge #158115 is investigated. In this particular discharge, the resonant amplitudes of the RMP are modulated in a cycloidal manner at a frequency of 1 Hz. Adopting the plausible hypothesis that mode penetration at the top of the pedestal is a necessary and sufficient condition for the RMP-induced suppression of edge localized modes (ELMs), recent cylindrical, nonlinear, reduced-magnetohydrodynamical (MHD) simulations performed by Hu, Nazikian, et al., (2019) can account, in a quantitative fashion, for the density-pump out and RMP-induced ELM suppression threshold observed in DIII-D discharge #158115. The primary aim of this paper is to employ analytic theory to further simplify the model of Hu, Nazikian, et al., in such a manner that a complete simulation of RMP-induced ELM-suppression in a DIII-D H-mode discharge can be performed in a matter of minutes of real time. A secondary aim is to gain a more exact understanding of the physical mechanism that underlies RMP-induced ELM suppression in the DIII-D tokamak.

Tables

Table 1. Table 1: Measured and estimated physics parameters at two rational surfaces in the pedestal of DIII-D discharge #158115 (see Figure 2 of [ 15 ] ). m 𝑚 m is the poloidal mode number, n 𝑛 n the toroidal mode number, B ϕ subscript 𝐵 italic-ϕ B_{\phi} the toroidal magnetic field (T), R 0 subscript 𝑅 0 R_{0} the major radius (m), a 𝑎 a the minor radius (m), n e subscript 𝑛 𝑒 n_{e} the electron number density ( 10  19 ​ m − 3 superscript 10 19 superscript m 3 10^{\,19}\,{\rm m}^{-3} ), T e subscript 𝑇 𝑒 T_{e} the electron temperature (keV), T i subscript 𝑇 𝑖 T_{i} the ion temperature (keV), η i = d ​ ln ⁡ T i / d ​ ln ⁡ n e subscript 𝜂 𝑖 𝑑 subscript 𝑇 𝑖 𝑑 subscript 𝑛 𝑒 \eta_{i}=d\ln T_{i}/d\ln n_{e} , Z eff subscript 𝑍 eff Z_{\rm eff} the conventional measure of impurity content, χ ⟂ subscript 𝜒 perpendicular-to \chi_{\perp} the perpendicular momentum diffusivity ( m  2 ​ s − 1 superscript m 2 superscript s 1 {\rm m}^{\,2}\,{\rm s}^{-1} ), ω E = E r / ( R 0 ​ B θ ) subscript 𝜔 𝐸 subscript 𝐸 𝑟 subscript 𝑅 0 subscript 𝐵 𝜃 \omega_{E}=E_{r}/(R_{0}\,B_{\theta}) the 𝐄 × 𝐁 𝐄 𝐁 {\bf E}\times{\bf B} frequency ( krad ​ s − 1 krad superscript s 1 {\rm krad}\,{\rm s}^{-1} ), ω ∗ e = ( d ​ p e / d ​ r ) / ( e ​ n e ​ R 0 ​ B θ ) subscript 𝜔 ∗ absent 𝑒 𝑑 subscript 𝑝 𝑒 𝑑 𝑟 𝑒 subscript 𝑛 𝑒 subscript 𝑅 0 subscript 𝐵 𝜃 \omega_{\ast\,e}=(dp_{e}/dr)/(e\,n_{e}\,R_{0}\,B_{\theta}) the electron diamagnetic frequency ( krad ​ s − 1 krad superscript s 1 {\rm krad}\,{\rm s}^{-1} ), r ^ s subscript ^ 𝑟 𝑠 \hat{r}_{s} the rational surface radius normalized to the plasma minor radius, s 𝑠 s the magnetic shear, M i subscript 𝑀 𝑖 M_{i} the majority ion mass number, Δ ≡ Δ ′ ​ r s / ( 2 ​ m ) Δ superscript Δ ′ subscript 𝑟 𝑠 2 𝑚 {\mit\Delta}\equiv{\mit\Delta}^{\prime}\,r_{s}/(2\,m) , and 𝒜 𝒜 {\cal A} is the amplification factor.

$m$	$n$	$B_{\varphi }$	$R_0$	$a$	$n_e$	$T_e$	$T_i$	${\eta }_i$	$Z_{\mathrm{eff}}$
-8	$2$	$-1.94$	1.75	0.93	2.8	1.4	1.4	1.9	2.5
-11	$2$	$-1.94$	1.75	0.93	0.75	0.12	0.12	1.8	2.5
$m$	$n$	${\chi }_{⟂}$	${\omega }_E$	${\omega }_{\ast e}$	${\hat{r}}_s$	$s$	$M_i$	$\Delta$	$𝒜$
-8	$2$	1.0	$-21.3$	$-21.5$	0.853	2.3	2.0	1.0	1.0
-11	$2$	1.0	$-9.8$	$-12.8$	0.974	11.2	2.0	1.0	1.0

Table 2. Table 2: Input parameters for the analytic, cylindrical, single-helicity, four-field, linear, resonant plasma response model of [ 20 ] at two rational surfaces in the pedestal of DIII-D discharge #158115. m 𝑚 m is the poloidal mode number, n 𝑛 n the toroidal mode number, c β = β subscript 𝑐 𝛽 𝛽 c_{\beta}=\sqrt{\beta} (where β 𝛽 \beta is the usual dimensionless measure of plasma pressure), D = S  1 / 3 ​ ρ s / r s 𝐷 superscript 𝑆 1 3 subscript 𝜌 𝑠 subscript 𝑟 𝑠 D=S^{\,1/3}\,\rho_{s}/r_{s} (where S = τ R / τ H 𝑆 subscript 𝜏 𝑅 subscript 𝜏 𝐻 S=\tau_{R}/\tau_{H} ), P = τ R / τ V 𝑃 subscript 𝜏 𝑅 subscript 𝜏 𝑉 P=\tau_{R}/\tau_{V} , and Q = S  1 / 3 ​ | ω 0 | ​ τ H / 2 𝑄 superscript 𝑆 1 3 subscript 𝜔 0 subscript 𝜏 𝐻 2 Q=S^{\,1/3}\,|\omega_{0}|\,\tau_{H}/2 , with ω 0 = − n ​ ( ω E + ω ∗ e ) subscript 𝜔 0 𝑛 subscript 𝜔 𝐸 subscript 𝜔 ∗ absent 𝑒 \omega_{0}=-n\,(\omega_{E}+\omega_{\ast\,e}) .

$m$	$n$	$c_{\beta }$	$D$	$P$	$Q$
-8	$2$	$6.48\times {10}^{-2}$	$1.54\times {10}^{+0}$	$7.16\times {10}^{+0}$	$1.27\times {10}^{+0}$
-11	$2$	$9.81\times {10}^{-3}$	$3.43\times {10}^{-1}$	$4.63\times {10}^{-1}$	$6.16\times {10}^{-2}$

Table 3. Table 3: Input parameters for the linear response model at two rational surfaces in the pedestal of DIII-D discharge #158115. m 𝑚 m is the poloidal mode number, n 𝑛 n the toroidal mode number, r ^ s subscript ^ 𝑟 𝑠 \hat{r}_{s} the rational surface radius normalized to the plasma minor radius, ϵ ≡ ( ϵ a / q s )  2 italic-ϵ superscript subscript italic-ϵ 𝑎 subscript 𝑞 𝑠 2 \epsilon\equiv(\epsilon_{a}/q_{s})^{\,2} , ν θ subscript 𝜈 𝜃 \nu_{\theta} the dimensionless poloidal flow damping parameter, ν ϕ subscript 𝜈 italic-ϕ \nu_{\phi} the dimensionless toroidal flow-damping parameter, ν μ subscript 𝜈 𝜇 \nu_{\mu} the dimensionless perpendicular viscosity parameter, L 𝐿 L the dimensionless locking parameter, and ω ^ 0 subscript ^ 𝜔 0 \hat{\omega}_{0} the normalized natural frequency. The latter quantity is calculated assuming that ω 0 = − n ​ ( ω E + ω ∗ e ) subscript 𝜔 0 𝑛 subscript 𝜔 𝐸 subscript 𝜔 ∗ absent 𝑒 \omega_{0}=-n\,(\omega_{E}+\omega_{\ast\,e}) .

$m$	$n$	${\hat{r}}_s$	$ϵ$	${\nu }_{\theta }$	${\nu }_{\varphi }$	${\nu }_{\mu }$	$L$	${\hat{\omega }}_0$
-8	$2$	$0.853$	$1.77\times {10}^{-2}$	$4.35\times {10}^{+2}$	$0.0$	$2.16\times {10}^{-3}$	$9.41\times {10}^{-3}$	$1.60\times {10}^{+2}$
-11	$2$	$0.974$	$9.34\times {10}^{-3}$	$6.50\times {10}^{+1}$	$0.0$	$1.32\times {10}^{-4}$	$1.60\times {10}^{-2}$	$5.16\times {10}^{+0}$

Table 4. Table 4: Important physical parameters at two rational surfaces in the pedestal of DIII-D discharge #158115. m 𝑚 m is the poloidal mode number, n 𝑛 n the toroidal mode number, ν ∗ e subscript 𝜈 ∗ absent 𝑒 \nu_{\ast\,e} the dimensionless electron collisionality parameter, ν ∗ i subscript 𝜈 ∗ absent 𝑖 \nu_{\ast\,i} the dimensionless ion collisionality parameter, τ SC subscript 𝜏 SC \tau_{\rm SC} the semi-collisional reconnection timescale (s), δ SC subscript 𝛿 SC \delta_{\rm SC} the semi-collisional layer width (m), ρ s subscript 𝜌 𝑠 \rho_{s} the ion sound radius (m), and b f / b v subscript 𝑏 𝑓 subscript 𝑏 𝑣 b_{f}/b_{v} the ratio of the normalized radial magnetic field to the vacuum radial magnetic field at the rational surface.

$m$	$n$	${\nu }_{\ast e}$	${\nu }_{\ast i}$	${\tau }_{\mathrm{SC}}$	${\delta }_{\mathrm{SC}}$	${\rho }_s$	$b_f/b_v$
-8	$2$	$9.26\times {10}^{-2}$	$2.97\times {10}^{-2}$	$1.86\times {10}^{-3}$	$4.14\times {10}^{-3}$	$2.79\times {10}^{-3}$	$1.45\times {10}^{+1}$
-11	$2$	$3.80\times {10}^{+0}$	$1.22\times {10}^{+0}$	$1.14\times {10}^{-4}$	$5.75\times {10}^{-3}$	$8.16\times {10}^{-4}$	$1.77\times {10}^{+0}$

Table 5. Table 5: Input parameters for the nonlinear response model at two rational surfaces in the pedestal of DIII-D discharge #158115. m 𝑚 m is the poloidal mode number, n 𝑛 n the toroidal mode number, r ^ s subscript ^ 𝑟 𝑠 \hat{r}_{s} the rational surface radius normalized to the plasma minor radius, ϵ ≡ ( ϵ a / q s )  2 italic-ϵ superscript subscript italic-ϵ 𝑎 subscript 𝑞 𝑠 2 \epsilon\equiv(\epsilon_{a}/q_{s})^{\,2} , ν θ subscript 𝜈 𝜃 \nu_{\theta} the dimensionless poloidal flow damping parameter, ν ϕ subscript 𝜈 italic-ϕ \nu_{\phi} the dimensionless toroidal flow-damping parameter, ν μ subscript 𝜈 𝜇 \nu_{\mu} the dimensionless perpendicular viscosity parameter, L 𝐿 L the dimensionless locking parameter, and ω ^ 0 subscript ^ 𝜔 0 \hat{\omega}_{0} the normalized natural frequency. The latter quantity is calculated assuming that ω 0 subscript 𝜔 0 \omega_{0} is given by ( 134 )–( 136 ).

$m$	$n$	${\hat{r}}_s$	$ϵ$	${\nu }_{\theta }$	${\nu }_{\varphi }$	${\nu }_{\mu }$	$L$	${\hat{\omega }}_0$
-8	$2$	$0.853$	$1.77\times {10}^{-2}$	$4.35\times {10}^{+2}$	$0.0$	$2.16\times {10}^{-3}$	$9.41\times {10}^{-3}$	$6.97\times {10}^{+1}$
-11	$2$	$0.974$	$9.34\times {10}^{-3}$	$6.50\times {10}^{+1}$	$0.0$	$1.32\times {10}^{-4}$	$1.60\times {10}^{-2}$	$-9.72\times {10}^{-1}$

Equations

$${\mu }_0{j}_{\varphi }(r)=\frac{1}{r}\frac{d(r{B}_{\theta })}{dr}.$$

$$q(r)=\frac{r{B}_{\varphi }}{R_0{B}_{\theta }},$$

$$\delta \mathbf{B}$$

$${\mu }_0\delta \mathbf{j}$$

$$\psi (r,\theta ,\varphi ,t)=\hat{\psi }(r,t)\exp [\mathrm{i}(m\theta -n\varphi )].$$

$${\text{mit}\Delta }^{\prime }={\left[\frac{d\ln \hat{\psi }}{dr}\right]}_{r_s-}^{r_s+},$$

$$\frac{{\delta }_{\mathrm{SC}}}{r_s}{\tau }_R\left(\frac{d}{dt}+\mathrm{i}\omega \right){\text{mit}\Psi }_s={\text{mit}\Delta }^{\prime }{r}_s{\text{mit}\Psi }_s+2|m|\mathcal{A}{\text{mit}\Psi }_v.$$

$${\delta }_{\mathrm{SC}}=\pi \frac{|n{\omega }_{\ast e}{|}^{1/2}{\tau }_H}{({\rho }_s/r_s){\tau }_R^{1/2}}{r}_s$$

$${\tau }_H$$

$${\tau }_R$$

$${\omega }_{\ast e}$$

$${\rho }_s$$

$$\omega (t)=m{\text{mit}\Omega }_{\theta }(r_s,t)-n{\text{mit}\Omega }_{\varphi }(r_s,t),$$

$${\tau }_V=\frac{r_s^2\rho (r_s)}{\mu (r_s)},$$

$$\frac{{\delta }_{\mathrm{SC}}}{r_s}{\tau }_R\frac{d}{dt}\left(\frac{W}{r_s}\right)=\frac{1}{2}\frac{W}{r_s}\left[{\text{mit}\Delta }^{\prime }{r}_s+2|m|\mathcal{A}{\left(\frac{W_v}{W}\right)}^2\cos \phi \right],$$

$$\frac{{\delta }_{\mathrm{SC}}}{r_s}{\tau }_R\left(\frac{d\phi }{dt}-\omega \right)=-2|m|\mathcal{A}{\left(\frac{W_v}{W}\right)}^2\sin \phi .$$

$$W(t)=4{\left(\frac{|{\text{mit}\Psi }_s|}{s{r}_s|B_{\theta }(r_s)|}\right)}^{1/2}{r}_s$$

$$W_v(t)=4{\left(\frac{|{\text{mit}\Psi }_v|}{s{r}_s|B_{\theta }(r_s)|}\right)}^{1/2}{r}_s$$

$$\phi (t)={\phi }_s(t)-{\phi }_v$$

$$\mathcal{I}{\tau }_R\frac{d}{dt}\left(\frac{W}{r_s}\right)={\text{mit}\Delta }^{\prime }{r}_s+2|m|\mathcal{A}{\left(\frac{W_v}{W}\right)}^2\cos \phi ,$$

$$\frac{d{\phi }_s}{dt}-\omega =0,$$

$$T_{\theta \mathrm{EM}}$$

$$T_{\varphi \mathrm{EM}}$$

$${\text{mit}\Omega }_{\theta }(r,t)$$

$${\text{mit}\Omega }_{\varphi }(r,t)$$

$$4{\pi }^2{R}_0\left[(1+2q^2)\rho r^3\frac{\partial {\text{mit}\Delta \Omega }_{\theta }}{\partial t}-\frac{\partial }{\partial r}\left(\mu r^3\frac{\partial {\text{mit}\Delta \Omega }_{\theta }}{\partial r}\right)+\rho r^3\frac{{\text{mit}\Delta \Omega }_{\theta }}{{\tau }_{\theta }}\right]$$

$$=T_{\theta \mathrm{EM}}\delta (r-r_s),$$

$$4{\pi }^2{R}_0^3\left[\rho r\frac{\partial {\text{mit}\Delta \Omega }_{\varphi }}{\partial t}-\frac{\partial }{\partial r}\left(\mu r\frac{\partial {\text{mit}\Delta \Omega }_{\varphi }}{\partial r}\right)+\rho r\frac{{\text{mit}\Delta \Omega }_{\varphi }}{{\tau }_{\varphi }}\right]=T_{\varphi \mathrm{EM}}\delta (r-r_s),$$

$$\frac{\partial {\text{mit}\Delta \Omega }_{\theta }(0,t)}{\partial r}$$

$${\text{mit}\Delta \Omega }_{\theta }(a,t)$$

$$\omega (t)={\omega }_0+m{\text{mit}\Delta \Omega }_{\theta }(r_s,t)-n{\text{mit}\Delta \Omega }_{\varphi }(r_s,t),$$

$${\omega }_0=m{\text{mit}\Omega }_{\theta 0}(r_s)-n{\text{mit}\Omega }_{\varphi 0}(r_s)$$

$${\omega }_0=-n({\omega }_E+{\omega }_{\ast e}),$$

$${\omega }_0=-n{\omega }_E-n\left(1-\frac{{\eta }_i{\lambda }_{\theta i}}{1+{\eta }_i}\right){\omega }_{\ast i},$$

$${\omega }_{\ast i}=-\frac{(d{p}_i/dr{)}_{r_s}}{e{n}_e(r_s)R_0{B}_{\theta }(r_s)}$$

Full text

We find that the response of the plasma to the applied RMP, in the immediate vicinity of a given rational surface, is governed by nonlinear, rather than linear, theory. This is the case because the magnetic island widths associated with driven reconnection exceed the linear layer widths, even in cases where driven reconnection is strongly suppressed by plasma rotation. We find that the natural frequency at a given rational surface (i.e., the helical frequency at which the locally resonant component of the RMP would need to propagate in order to maximize driven reconnection) is offset from the local ${\bf E}\times{\bf B}$ frame in the ion diamagnetic direction. The size of the offset is mostly determined by neoclassical poloidal flow damping. Finally, our analytic nonlinear response model is found to be largely consistent with the simulations of Hu, Nazikian et al, and also correctly predicts the RMP-induced mode penetration threshold in DIII-D discharge #158115.

\NF

1 Introduction

Tokamak discharges operating in high-confinement mode (H-mode) [1] exhibit intermittent bursts of heat and particle transport, emanating from the outer regions of the plasma, that are known as type-I edge localized modes (ELMs) [2]. ELMs are fairly harmless in present-day tokamaks possessing carbon plasma-facing components. However, large ELMs can cause a problematic influx of tungsten ions into the plasma core in tokamaks possessing tungsten plasma-facing components [3]. Moreover, it is estimated that the heat load that ELMs will deliver to the tungsten plasma-facing components in a reactor-scale tokamak will be large enough to cause massive tungsten ion influx into the core, and that the erosion associated with this process will unacceptably limit the lifetimes of these components [4]. Consequently, developing robust and effective methods for ELM control is a high priority for the international magnetic fusion program.

The most promising method for the control of ELMs is via the application of static resonant magnetic perturbations (RMPs). Complete RMP-induced ELM suppression was first demonstrated on the DIII-D tokamak [5]. Subsequently, either mitigation or compete suppression of ELMs has been demonstrated on the JET [6], ASDEX-U [7], KSTAR [8], and EAST [9] tokamaks.

ELMs are thought to be caused by peeling-ballooning instabilities, with intermediate toroidal mode numbers, that are driven by the strong pressure gradients and current density gradients characteristic of the edge region of an H-mode tokamak discharge [10], which is known as the pedestal region. The initial observations of RMP-induced ELM suppression were interpreted as an indication that the magnetic field in the pedestal is rendered stochastic by the applied RMP, leading to greatly enhanced transport via thermal diffusion along magnetic field-lines [5, 11]. This explanation was quickly abandoned because no significant reduction in the electron temperature gradient in the pedestal is observed during RMP-induced ELM suppression experiments, whereas a very significant reduction would be expected in the presence of stochastic fields. It is now generally accepted that response currents generated within the pedestal, as a consequence of plasma rotation, play a crucial role in the perturbed equilibrium in the presence of RMPs, and that these currents act to prevent the formation of RMP-driven magnetic island chains—a process known as shielding—and, thereby, significantly reduce the stochasticity of the magnetic field [12].

This paper concentrates on a particular (but completely typical) DIII-D H-mode discharge (#158115 [13]) in which ELMs were successfully suppressed by an externally applied $n=2$ RMP. In this discharge, the relative phase of magnetic perturbations generated by two sets of external field coils is modulated sinusoidally at a frequency of 1 Hz, causing the amplitudes of the helical harmonics of the applied RMP that resonate in the pedestal region to modulate in a cycloidal manner at the same frequency. The application of the RMP is observed to generate two distinct plasma responses. The first response—which is known as the density pump-out—is characterized by a reduction in the pedestal density, accompanied by a much smaller reduction in the pedestal temperature, whose magnitude varies smoothly with the amplitude of the edge-resonant components of the RMP. The density pump-out is not observed to be associated with ELM suppression. The second response—which is known as mode penetration—occurs when the amplitude of the edge-resonant harmonics of the RMP exceeds a certain critical value, is associated with the sudden formation of a locked magnetic island chain at a rational surface that lies close to the top of the pedestal, and is accompanied by a sudden shift in the edge ion toroidal rotation. Mode penetration is observed to be strongly correlated with ELM suppression. If the amplitude of the edge-resonant harmonics of the RMP falls below a second, somewhat smaller, critical value then the locked magnetic island chain is expelled from the plasma, and the edge ion toroidal rotation returns to its original value—this process is known as mode unlocking. Mode unlocking is invariably accompanied by the resumption of ELMs. Note that we can be sure that mode penetration is associated with the formation of a locked magnetic island chain because when mode unlocking occurs the chain is observed to spin-up and decay away (see Figure 5 in [13]). The fact that mode penetration is a necessary condition for ELM suppression is a reasonable inference from existing DIII-D experimental data [13]. It is the hypothesis of this paper that mode penetration is also a sufficient condition. Incidentally, it is not difficult to understand why the formation of a locked island chain close to the top of the pedestal might give rise to ELM suppression. The flattening of the temperature and density profiles across the island region reduces the pressure gradient at the top of the pedestal, which is likely to move the plasma farther from the peeling-ballooning stability threshold. Moreover, the formation of a locked island is likely to enhance any nonlinear interaction between the applied RMP and the peeling-ballooning mode (essentially by allowing the associated magnetic fields to phase lock to one another).

A complete numerical simulation of DIII-D discharge #158115 would entail running a nonlinear extended (because the code would need to incorporate diamagnetic and neoclassical effects) full-magnetohyrodynamical (MHD) (because ELMs cannot be modeled using reduced-MHD) code in toroidal geometry (because ELMs cannot be modeled in cylindrical geometry) for approximately $10^{\,8}$ Alfvén times (which corresponds to 1 second of experimental time). Unfortunately, this is completely impossible. (The most advanced, current-day, nonlinear, toroidal, extended-MHD codes are typically only capable of running for $10^{\,4}$ Alfvén times [14].) Hence, to make further progress, some sort of reduced model is required. If we accept that mode penetration is a sufficient condition for the achievement of ELM suppression then the problem is greatly simplified, because we do not need to directly model ELMs, and, thus, we can perform a nonlinear reduced-MHD calculation in cylindrical geometry. This is precisely the approach taken in a recent paper by Hu, Nazikian et al. [15]. In this paper, computer simulations—made using the cylindrical, multi-harmonic, five-field, nonlinear, initial-value code, TM1 [16, 17, 18]—of DIII-D discharge #158115 find that the formation of RMP-driven magnetic island chains at the bottom and the top of the pedestal can account, in a quantitative fashion, for the observed density pump-out, as well as the mode-penetration-induced ELM suppression threshold [15]. The simulations also find that driven magnetic reconnection is strongly shielded in the middle of the pedestal. (Incidentally, the five fields in TM1 are the poloidal magnetic flux, the electron density, the perpendicular ion vorticity, the parallel ion velocity, and the electron temperature.) It should be noted that the simulations presented in [15] use the experimental density, temperature, safety-factor, and radial electric field profiles, apply an appropriately modulated $n=2$ RMP with the experimental spectrum of resonant harmonics, and simulate the response of the plasma for the equivalent of 1 second of experimental time.

Unfortunately, even the reduced model of [15] requires many tens of hours of cpu time (corresponding to many days of actual time) to perform a complete simulation. It is the primary aim of this paper to employ analytic theory to further simplify the model in such a manner that a complete simulation can be performed in a matter of minutes in real time. A secondary aim is to gain a more exact understanding of the physical mechanism that underlies RMP-induced ELM suppression in the DIII-D tokamak. Two analytic theories are employed in this study. The first is the cylindrical, single-harmonic, four-field, linear, resonant plasma response model of [19, 20]. The second is the cylindrical, single-harmonic, four-field, nonlinear, resonant plasma response model of [21, 22, 23]. (Note that, unlike the model used in [15], for the sake of simplicity, the models used in this paper do not evolve the electron temperature profile separately from the density profile. Hence, they are four-field, rather than five-field, models.)

2 Preliminary Analysis

2.1 Plasma Equilibrium

Consider a large aspect-ratio, low-$\beta$, tokamak plasma whose equilibrium magnetic flux surfaces map out (almost) concentric circles in the poloidal plane. Such a plasma is well approximated as a periodic cylinder. Suppose that the minor radius of the plasma is $a$. Standard right-handed cylindrical coordinates ($r$, $\theta$, $z$) are adopted. The system is assumed to be periodic in the $z$-direction, with periodicity length $2\pi\,R_{0}$, where $R_{0}\gg a$ is the simulated plasma major radius. It is convenient to define the simulated toroidal angle $\phi=z/R_{0}$.

The equilibrium magnetic field is written ${\bf B}({\bf r})=[0,\,B_{\theta}(r),\,B_{\phi}]$. The associated equilibrium plasma current density takes the form ${\bf j}({\bf r})=[0,\,0,\,j_{\phi}(r)]$, where

[TABLE]

The safety factor,

[TABLE]

parameterizes the helical pitches of equilibrium magnetic field-lines. In a conventional tokamak plasma, $|q(r)|$ is of order unity, and is a monotonically increasing function of $r$.

2.2 Plasma Response

Consider the response of the plasma to a static RMP. Suppose that the RMP has $|m|$ periods in the poloidal direction, and $n>0$ periods in the toroidal direction. (Note that $m$ is positive if $q$ is positive, and vice versa.) It is convenient to express the perturbed magnetic field and the perturbed plasma current density in terms of a magnetic flux-function, $\psi(r,\theta,\phi,t)$. Thus,

[TABLE]

where

[TABLE]

This representation is valid provided that $|m|/n\gg a/R_{0}$ [24].

As is well known, the response of the plasma to the applied RMP is governed by the equations of perturbed, marginally-stable (i.e., $\partial/\partial t\equiv 0$), ideal magnetohydrodynamics (MHD) everywhere in the plasma, apart from a relatively narrow (in $r$) region in the vicinity of the rational surface, minor radius $r_{s}$, where $q(r_{s})=m/n$ [24].

It is convenient to parameterize the RMP in terms of the so-called vacuum flux, ${\mit\Psi}_{v}(t)=|{\mit\Psi}_{v}|\,{\rm e}^{-{\rm i}\,\varphi_{v}}$, which is defined to be the value of $\hat{\psi}(r,t)$ at radius $r_{s}$ in the presence of the RMP, but in the absence of the plasma. Here, $\varphi_{v}$ is the helical phase of the RMP, and is assumed to be constant in time. Likewise, the response of the plasma in the vicinity of the rational surface to the RMP is parameterized in terms of the so-called reconnected flux, ${\mit\Psi}_{s}(t)=|{\mit\Psi}_{s}|\,{\rm e}^{-{\rm i}\,\varphi_{s}}$, which is the actual value of $\hat{\psi}(r,t)$ at radius $r_{s}$. Here, $\varphi_{s}(t)$ is the helical phase of the reconnected flux.

The intrinsic stability of the $m$/$n$ tearing mode is governed by the tearing stability index [25],

[TABLE]

where $\hat{\psi}(r)$ is a solution of the marginally-stable, ideal-MHD equations, for the case of an $m$/$n$ helical perturbation, that satisfies physical boundary conditions at $r=0$ and $r=a$ (in the absence of the RMP). According to resistive-MHD theory [25, 26], if ${\mit\Delta}^{\prime}>0$ then the $m$/$n$ tearing mode spontaneously reconnects magnetic flux at the rational surface to form a helical magnetic island chain. In the following, it is assumed that ${\mit\Delta}^{\prime}<0$, so that the $m$/$n$ tearing mode is intrinsically stable. In this situation, any magnetic reconnection that takes place at the rational surface is due solely to the RMP.

2.3 Linear Response Regime

In this paper, we shall examine two different plasma response regimes at the rational surface. The first of these is the so-called semi-collisonal regime [20, 27, 28]. This is a linear, two-fluid, low-collisionality regime in which the reconnected magnetic flux induced by the RMP is governed by [21]

[TABLE]

Here,

[TABLE]

is the linear layer width, whereas

[TABLE]

are the hydromagnetic timescale, resistive diffusion timescale, electron diamagnetic frequency, and ion sound radius, respectively, at the rational surface. Moreover, $s=(d\ln q/d\ln r)_{r=r_{s}}$ is the local magnetic shear, ${\cal A}$ the amplification factor (i.e., the factor by which the radial magnetic field at the rational surface due to the RMP is enhanced with respect to its vacuum value due to equilibrium plasma currents external to the rational surface), $m_{i}$ the ion mass, and $e$ the magnitude of the electron charge. Furthermore, $\sigma(r)$, $n_{e}(r)$, $\rho(r)\equiv m_{i}\,n_{e}(r)$, $T_{e}(r)$, $p_{e}(r)\equiv n_{e}(r)\,T_{e}(r)$ are the equilibrium plasma electrical conductivity, electron number density, mass density, electron temperature, and electron pressure profiles, respectively. Finally,

[TABLE]

where ${\mit\Omega}_{\theta}(r,t)$ and ${\mit\Omega}_{\phi}(r,t)$ are the plasma poloidal and toroidal angular velocity profiles, respectively. [To be more exact, ${\mit\Omega}_{\theta}(r,t)$ and ${\mit\Omega}_{\phi}(r,t)$ are the poloidal and toroidal angular velocity profiles of an imaginary fluid that convects reconnected magnetic flux at rational surfaces. It is assumed that changes in these velocity profiles are mirrored by changes in the actual plasma velocity profiles. A magnetic island convected by the imaginary fluid propagates at its so-called natural frequency. The relationship between the natural frequency and the ${\bf E}\times{\bf B}$ frequency is specified in Section 2.6.] It is helpful to define the viscous diffusion timescale at the rational surface,

[TABLE]

where $\mu(r)$ is the equilibrium plasma (perpendicular) viscosity profile.

It should be noted that the analysis of [20], combined with the experimental data listed in Table 1, leads to the conclusion that the appropriate cylindrical, four-field, linear, plasma response regime at the rational surface is the so-called SCi (first semi-collisional) regime. [In particular, the dimensionless parameters $c_{\beta}$, $D$, $P$, and $Q$ that control the plasma response, according to the analysis of [20], are calculated at the -8/2 and the -11/2 rational surfaces in the pedestal of DIII-D discharge #158115 [13] in Table 2. In the first case, the fact that $D>1$ implies that we should consult Figure 3 in [20]. According to this figure, the fact that $c_{\beta}\,D<Q<D$ and $P>1$ indicates the appropriate response regime at the -8/2 surface is the SCi regime. In the second case, the fact that $c_{\beta}^{\,1/3}<D<1$ implies that we should consult Figure 2 in [20]. According to this figure, the fact that $Q<D$ and $P>Q^{\,3}/D^{\,6}$ indicates the the appropriate response regime at the -11/2 surface is the SCi regime. The response and layer thickness in the SCi regime are listed in Table 1 of [20].] Equation (7) is a slightly simplified implementation of the response of the plasma in the vicinity of the rational surface in the SCi regime.

Incidentally, the experimental data listed in Table 1 is derived from the safety-factor, electron number density, electron temperature, and ${\bf E}\times{\bf B}$ frequency profiles shown in Figure 2 of [15]. The electron diamagnetic frequency, magnetic shear, and $\eta_{i}$ values listed in the table are calculated directly from these profiles. The electron and ion temperature profiles are assumed to be the same. The values of the perpendicular momentum diffusivity are obtained from the TRANSP code [29]. The values of the effective ion charge number, $Z_{\rm eff}$, come from line emission spectroscopy.

Equation (7) can be conveniently rewritten as an island width evolution equation [21],

[TABLE]

and an island phase evolution equation,

[TABLE]

Here,

[TABLE]

is the full (radial) width of the magnetic island chain that forms at the rational surface. (Incidentally, it is assumed that $W\ll r_{s}$.) Moreover,

[TABLE]

is termed the vacuum island width. Finally,

[TABLE]

is the helical phase of the island chain relative to the RMP.

The semi-collisional response regime holds when $W<\delta_{\rm SC}$. That is, when the magnetic island width falls below the linear layer width.

2.4 Nonlinear Response Regime

The second response regime investigated in this paper is the so-called Rutherford regime [30]. This is a nonlinear regime in which the reconnected magnetic flux induced by the RMP is governed by two equations. The first of these is the Rutherford island width evolution equation [26],

[TABLE]

where ${\cal I}=0.8227$. The second governing equation is the so-called no-slip constraint [24],

[TABLE]

according to which the island chain is forced to co-rotate with the aforementioned (see Section 2.3) magnetic-flux-convecting imaginary fluid at the rational surface.

The Rutherford response regime holds when $W>\delta_{\rm SC}$. That is, when the magnetic island width exceeds the linear layer width.

2.5 Plasma Angular Velocity Evolution

It is easily demonstrated that zero net electromagnetic torque can be exerted on magnetic flux surfaces located in a region of the plasma that is governed by the equations of marginally-stable, ideal-MHD [24]. Thus, any electromagnetic torque exerted on the plasma by the RMP develops in the immediate vicinity of the rational surface, where ideal-MHD breaks down. The net poloidal and toroidal electromagnetic torques exerted in the vicinity of the rational surface by the RMP take the forms [24, 30]

[TABLE]

respectively.

We can write

[TABLE]

where ${\mit\Omega}_{\theta\,0}(r)$ and ${\mit\Omega}_{\phi\,0}(r)$ are the equilibrium poloidal and toroidal plasma angular velocity profiles, respectively, whereas ${\mit\Delta\Omega}_{\theta}(r,t)$ and ${\mit\Delta\Omega}_{\phi}(r,t)$ are the respective modifications to these profiles induced by the aforementioned electromagnetic torques. The modifications to the angular velocity profiles are governed by the poloidal and toroidal angular equations of motion of the plasma, which take the respective forms [24, 31]

[TABLE]

and are subject to the spatial boundary conditions [24]

[TABLE]

Here, $\tau_{\theta}(r)$ is the neoclassical poloidal flow-damping time profile [32], and $\tau_{\phi}(r)$ the neoclassical toroidal flow-damping time profile. The neoclassical toroidal flow-damping is assumed to be generated by non-resonant components of the applied RMP [33, 34]. The factor $(1+2\,q^{\,2})$ in (2.5) derives from the fact that incompressible poloidal flow has a poloidaly-varying toroidal component that effectively increases the plasma mass being accelerated by the poloidal flow-damping force [31]. It turns out that, in the presence of strong poloidal and toroidal flow-damping, the modifications to the plasma poloidal and toroidal angular velocity profiles are localized in the vicinity of the rational surface [34]. Assuming that this is the case, it is a good approximation to replace $q$, $\rho$, $\mu$, $\tau_{\theta}$, and $\tau_{\phi}$ in (2.5) and (27) by their values at the rational surface. We are, nevertheless, assuming that the localization width greatly exceeds the linear layer width (or the island width, in the nonlinear case).

Equations (13), (24), and (25) imply that

[TABLE]

where

[TABLE]

is the so-called natural frequency of the $m$/$n$ tearing mode. In other words, $\omega_{0}$ is the helical phase velocity of a naturally unstable $m$/$n$ tearing mode in the absence of the RMP.

2.6 Natural Frequency

According to the cylindrical, single-helicity, four-field, linear analysis of [20], the appropriate natural frequency for the linear response model is

[TABLE]

where $\omega_{E}=E_{r}(r_{s})/[R_{0}\,B_{\theta}(r_{s})]$ is the ${\bf E}\times{\bf B}$ frequency at the rational surface, and $E_{r}(r)$ is the equilibrium radial electric field profile.

According to the cylindrical, single-helicity, four-field, nonlinear analysis of [23], the appropriate natural frequency for the nonlinear response model is

[TABLE]

where

[TABLE]

is the ion diamagnetic frequency at the rational surface, $\eta_{i}=(d\ln T_{i}/d\ln n_{e})_{r=r_{s}}$, $T_{i}(r)$ the equilibrium ion temperature profile, $p_{i}(r)=n_{e}(r)\,T_{i}(r)$ the equilibrium ion pressure profile, and the dimensionless parameter $\lambda_{\theta\,i}$ is specified in B.

3 Linear Response Model

3.1 Unnormalized Linear Response Model

According to Sections 2.3 and 2.5, the complete linear response model takes the form

[TABLE]

where $q_{s}=m/n$.

3.2 Normalized Linear Response Model

It is helpful to define the typical semi-collisional magnetic reconnection timescale,

[TABLE]

Let $\hat{r}=r/a$, $\hat{t}=t/\tau_{\rm SC}$, $\hat{\omega}_{0}=\omega_{0}\,\tau_{\rm SC}$, $\hat{\omega}_{\theta}(\hat{r},\hat{t})=-m\,{\mit\Delta\Omega}_{\theta}(r,t)\,\tau_{\rm SC}$, $\hat{\omega}_{\phi}(\hat{r},\hat{t})=n\,{\mit\Delta\Omega}_{\phi}(r,t)\,\tau_{\rm SC}$, $\hat{W}=W/\delta_{\rm SC}$, $\hat{W}_{v}=W_{v}/\delta_{\rm SC}$. The normalized linear response model reduces to

[TABLE]

where $\hat{r}_{s}=r_{s}/a$, $\hat{\mit\Delta}^{\prime}=(-{\mit\Delta}^{\prime}\,r_{s})/(2\,|m|)$, $\epsilon_{a}=a/R_{0}$, $\nu_{\theta}=\tau_{\rm SC}/\tau_{\theta}$, $\nu_{\phi}=\tau_{\rm SC}/\tau_{\phi}$, $\nu_{\mu}=(\tau_{\rm SC}/\tau_{V})\,(r_{s}/a)^{2}$, $\hat{W}_{0}=W_{0}/\delta_{\rm SC}$, and

[TABLE]

3.3 Solution of Plasma Angular Equations of Motion

We can solve (46)–(49) by writing [35]

[TABLE]

where

[TABLE]

The solution is exact in the limit $N\rightarrow\infty$. Here, $J_{m}(z)$ is a standard Bessel function, and $j_{m,n}$ denotes the $n$th zero of the $J_{m}(z)$ Bessel function [36]. It is easily demonstrated that

[TABLE]

and [37]

[TABLE]

Hence, we obtain

[TABLE]

where

[TABLE]

and

[TABLE]

3.4 Final Form of Normalized Linear Response Model

The normalized linear model reduces to the following closed set of equations:

[TABLE]

where $b_{f}={\cal A}\,\hat{W}_{v}^{\,2}$, and $L=\hat{W}_{0}^{-4}$.

If we define

[TABLE]

then the previous set of equations yield

[TABLE]

4 Nonlinear Response Model

4.1 Unnormalized Nonlinear Response Model

According to Sections 2.4 and 2.5, the complete nonlinear response model takes the form

[TABLE]

4.2 Normalized Nonlinear Response Model

The normalized form of the nonlinear response model is

[TABLE]

4.3 Final Form of Normalized Nonlinear Response Model

The normalized nonlinear model reduces to the following closed set of equations:

[TABLE]

Here, we have reused the analysis of Section 3.3 to solve the plasma angular equations of motion.

If we express the previous equations in terms of the variables $X$ and $Y$, which are defined in (71) and (72), then we get

[TABLE]

where

[TABLE]

Here, a slight modification has been made to (97) and (98) in the linear regime, $\hat{W}<1$, (in which case they are not valid anyway) in order to render them non-singular at the origin of the $X$-$Y$ plane.

5 Results

5.1 Introduction

According to [15], the density pump-out in DIII-D discharge #158115 is due to mode penetration at the $m=-11$/$n=2$ rational surface, which lies at the bottom of the pedestal. The formation of a magnetic island chain at the -11/2 surface leads to a local flattening of the plasma temperature and density profiles via parallel transport along magnetic field-lines. In order to flatten the profiles, the island width must exceed certain critical values that depend on the ratios of the relevant parallel and perpendicular diffusivities at the rational surface [38]. For the case of the electron temperature, making use of the analysis of [38], as well as the data in Table 1, we estimate the critical width to be

[TABLE]

where $v_{e}=\sqrt{T_{e}/m_{e}}$, $\epsilon_{s}=r_{s}/R_{0}$, $\chi_{\perp}$ is the perpendicular diffusivity, and $m_{e}$ the electron mass. For the case of the density, we estimate the critical width to be

[TABLE]

where $v_{i}=\sqrt{T_{i}/m_{i}}$. If both profiles are flattened across the island chain then, as a consequence of the fact that the equilibrium density gradient at the -11/2 surface greatly exceeds the temperature gradient (i.e., because $n_{e}/(T_{i}\,\eta_{i})\simeq 3.5\,\,[10^{\,19}\,{\rm m}^{-3}/({\rm keV})]$; see Table 1, as well as Figure 2 of [15]), the flattening naturally produces a much larger reduction in the pedestal density (measured in units of $10^{\,19}\,{\rm m}^{-3}$) than in the pedestal temperature (measured in units of keV).

According to [15], ELM suppression in DIII-D discharge #158115 is due to mode penetration at the $m=-8$/$n=2$ rational surface, which lies at the top of the pedestal. The formation of a magnetic island chain at the -8/2 surface leads to a local flattening of the plasma temperature and density profiles. As before, the island width must exceed certain critical values to flatten the profiles. We estimate the critical width required to flatten the electron temperature profile to be [cf. (103) and Table 1]

[TABLE]

whereas the critical width required to flatten the density profile is [cf. (104) and Table 1]

[TABLE]

The flattening of the temperature and density profiles at the -8/2 rational surface is presumed to be sufficient to prevent the plasma in the pedestal from ever exceeding the peeling-ballooning stability threshold, which leads to ELM suppression.

5.2 Linear Simulations

Table 1 shows measured and estimated physics parameters at the $m=-8$/$n=2$ and $m=-11$/$n=2$ rational surfaces in DIII-D discharge #158115. Incidentally, all experimental minor radii quoted in this paper are flux-surfaced-averaged minor radii, rather than minor radii on the outboard mid-plane. Making use of the analysis contained in the Appendices, these parameters can be used to derive the input parameters for the linear model [i.e., (73)–(77)] that are listed in Table 3. Note that the dimensionless toroidal flow-damping rate has been set to zero (mostly because there is insufficient data to calculate its value). Note, further, that our knowledge of the edge current profile in DIII-D discharge #158115 is insufficient to allow us to calculate the parameters $\hat{\mit\Delta}^{\prime}$ and ${\cal A}$ with any degree of accuracy. Fortunately, the penetration threshold exhibits virtually no dependence on the parameter $\hat{\mit\Delta}^{\prime}$, which can, therefore, safely be given its vacuum value of unity. With less justification, the amplification parameter, ${\cal A}$, is also given its vacuum value unity. In DIII-D discharge #158115, the relative phase of RMPs generated by two sets of external field coils is modulated sinusoidally at a frequency of 1 Hz. This causes the amplitudes of the resonant harmonics of the applied RMP to modulated in a cycloidal manner at the same freqeuncy.

Figure 1 shows a linear simulation of the response of the plasma at the -8/2 rational surface in DIII-D discharge #158115 to an RMP whose magnitude is switched on at $t=0$, and then modulated cycloidally at a frequency of 1 Hz. All simulations in this paper are performed with 200 velocity harmonics [i.e., $N=200$ in (77) and (101)]. The simulation data presented in Figure 1 is qualitatively similar to the experimental data shown in Figure 1 of [15], as well as the TM1 simulation data displayed in Figure 2 of the same paper. In particular, if the amplitude of the applied RMP rises above a certain threshold value then there is a bifurcation from a shielded solution characterized by $\hat{W}\ll\sqrt{b_{f}}$ to an unshielded solution characterized by $\hat{W}\sim\sqrt{b_{f}}$. This bifurcation, which is known as mode penetration, is accompanied by a sudden reduction in the natural frequency to zero, as well as sudden shifts in the plasma poloidal and toroidal rotation at the rational surface. Figure 2 shows the mode penetration process in more detail. Furthermore, referring again to Figure 1, if the amplitude of the applied RMP falls below a second smaller threshold value then there is a bifurcation from an unshielded solution to a shielded solution. This bifurcation, which is known as mode unlocking, is accompanied by the recovery of the natural frequency to its unperturbed value, as well as sudden shifts in the plasma poloidal and toroidal rotation at the rational surface. Figure 3 shows the mode unlocking process in more detail. Note that, immediately after unlocking, the magnetic island chain at the rational surface spins-up and decays, before eventually re-locking to the RMP in a fixed helical phase relation.

As is apparent from Figure 1, mode penetration is triggered as soon as the natural frequency has been reduced to approximately one half of its unperturbed value [24]. Moreover, about one third of the reduction in the natural frequency associated with mode penetration is due to a shift in the local plasma poloidal rotation, whereas two thirds is due to a shift in the local toroidal rotation. In the absence of neoclassical poloidal flow-damping, only about 2% (i.e., a fraction $\epsilon$—see Table 3) of the change in the natural frequency would be due to a shift in local plasma toroidal rotation. The large increase in the fraction of the overall frequency change due to toroidal velocity shifts, in the presence of neoclassical poloidal flow-damping, is a consequence of the fact that the majority of charged particles in the DIII-D pedestal are trapped in banana orbits, and cannot, therefore, freely rotate in the poloidal direction. Hence, neoclassical poloidal flow-damping is comparatively strong in the pedestal region of the DIII-D tokamak (compared to that in the plasma core). (The same is likely to be true in the pedestal regions of all relatively large, conventional aspect-ratio tokamaks.) The mode-penetration-induced shift in the plasma toroidal rotation shown in Figure 1 corresponds to a toroidal velocity shift in the co-current direction of about $30\,{\rm km/s}$ (see Table 4), which is in agreement with experimental observations [15]. The mode-penetration-induced shift in the plasma poloidal rotation shown in Figure 1 corresponds to a poloidal velocity shift in the ion diamagnetic direction of only about $2\,{\rm km/s}$. Furthermore, we would expect the latter velocity shift to be strongly localized in the vicinity of the rational surface [34]. Hence, it is not surprising that the poloidal velocity shift is not observed experimentally [39].

According to Figures 2 and 3, the sudden collapse/recovery of the natural frequency associated with mode penetration/mode unlocking takes place on a timescale of about a millisecond (see Table 4), and is accompanied by a sudden shift in the local plasma poloidal rotation that takes place on the same timescale. The mode penetration/mode unlocking -induced shift in the local plasma toroidal rotation takes place on a significantly longer timescale (at least, 10 milliseconds).

There is one major difference between the simulation data shown in Figure 1 and the experimental and TM1 simulation data shown in [15]. According to Figure 1, mode penetration at the -8/2 rational surface in DIII-D discharge #158115 occurs when $b_{f}$ exceeds the critical value $649$, which corresponds to a vacuum radial field at the -8/2 surface of $b_{v}=45$ gauss. (See Table 4.) However, the penetration threshold inferred from experimental data is more like $b_{v}=6$ gauss [15]. This discrepancy is probably related to the fact that, according to Figures 1 and 2, mode penetration at the -8/2 rational surface in DIII-D discharge #158115 occurs when $W\sim 4\,\delta_{\rm SC}$. In other words, when the magnetic island width exceeds the linear layer width. Given that the linear model is only valid when $W<\delta_{\rm SC}$, we conclude that mode penetration at the -8/2 rational surface in DIII-D discharge #158115 is actually governed by nonlinear physics. This particular conclusion is not consistent with the TM1 simulations described in [15], according to which mode penetration at the -8/2 rational surface seems to be governed by linear physics (because of the absence of island pulsations in the TM1 simulations—see Section 5.3). One possible explanation for this disagreement is that plasma perpendicular viscosity is artificially increased by a large factor in TM1 simulations, in order to mimic the effect of strong neoclassical poloidal flow-damping, but such an increase may also artificially increase the linear layer width. Note, from Table 4, that the semi-collisional layer width at the -8/2 rational surface is only $4\,{\rm mm}$, which is similar to the ion sound radius. On the other hand, the critical island width above which mode penetration is triggered is about $1.6\,{\rm cm}$. Incidentally, after mode penetration has occurred, the -8/2 island width rises well above the values required to locally flatten the temperature and density profiles. (See Section 5.1.)

Figure 4 shows a linear simulation of the response of the plasma at the -11/2 rational surface in DIII-D discharge #158115 to an RMP whose magnitude is switched on at $t=0$, and then modulated cycloidally at a frequency of 1 Hz. The simulation data presented in Figure 4 is qualitatively similar to the TM1 simulation data shown in Figure 2 of [15]. In particular, it is clear from Figure 4 that there is insufficient plasma rotation at the -11/2 rational surface to enable the effective shielding of driven magnetic reconnection. In other words, $\hat{W}\sim\sqrt{b_{f}}$ at all times. The peak -11/2 island width is sufficiently large to locally flatten the temperature and density profiles. (See Section 5.1.) Note, however, that $W>\delta_{\rm SC}$, which implies that driven reconnection at the -11/2 rational surface in DIII-D discharge #158115 is actually governed by nonlinear physics.

5.3 Nonlinear Simulations

Making use of the analysis contained in the Appendices, the experimental data given in Table 1 can be used to derive the input parameters for the nonlinear model [i.e., (97)–(102)] that are listed in Table 5. Note that the dimensionless toroidal flow-damping rate has again been set to zero.

Figure 5 shows a nonlinear simulation of the response of the plasma at the -8/2 rational surface in DIII-D discharge #158115 to an RMP whose magnitude is switched on at $t=0$, and then modulated cycloidally at a frequency of 1 Hz. The simulation data presented in Figure 5 is qualitatively similar to the experimental data shown in Figure 1 of [15]. As before, if the amplitude of the applied RMP rises above a certain threshold value then mode penetration occurs. In other words, there is bifurcation from a shielded to an unshielded solution, accompanied by a sudden reduction in the natural frequency to zero, as well as sudden shifts in the plasma toroidal and poloidal rotation at the rational surface. Figure 6 shows the mode penetration process in more detail. Furthermore, again referring to Figure 5, if the amplitude of the applied RMP falls below a second smaller threshold value then mode unlocking occurs. In other words, there is a bifurcation from an unshielded solution to a shielded solution, accompanied by the recovery of the natural frequency to its unperturbed value, as well as sudden shifts in the plasma toroidal and poloidal rotation at the rational surface. Figure 7 shows the mode unlocking process in more detail. As in the linear case, mode penetration is triggered when the natural frequency has been reduced to about half of its original value. Moreover, about two thirds of the change in the natural frequency associated with mode penetration/mode unlocking is due to a shift in the local toroidal plasma rotation, and about one third to a shift in the local poloidal plasma rotation.

The main difference between the nonlinear simulations discussed in this section and the linear simulations discussed in Section 5.2 lies in the nature of the shielded solution. In the linear simulations, the shielded solution consists of a narrow magnetic island chain that has a fixed helical phase shift of about $+\pi/2$ with respect to the locally resonant component of the RMP. On the other hand, in the nonlinear simulations, the shielded solution consists of a narrow island chain whose helical phase continually increases in time, and whose width periodically falls to zero, at which times its helical phase jumps by $-\pi$ radians. This type of pulsating island solution was first predicted in [30], is discussed in detail in [21, 22, 23], and has been observed experimentally (see Figure 29 in [40]). Note, from Figures 5 and 6, that the width of the pulsating island chain exceeds the linear layer width (i.e., $\hat{W}>1$) during most of its cycle, which justifies the nonlinear approach employed in this section.

Figure 4 implies that mode penetration at the -8/2 rational surface in DIII-D discharge #158115 occurs when $b_{f}$ exceeds the critical value $80$, which corresponds to a vacuum radial field at the -8/2 surface of $b_{v}=5.5$ gauss (see Table 4). This estimate for the penetration threshold is close to the experimentally inferred value of 6 gauss [15].

According to Section 2.6, the natural frequency in the linear regime is

[TABLE]

where $\omega_{\perp\,e}=\omega_{E}+\omega_{\ast\,e}$. On the other hand, the natural frequency in the nonlinear regime is

[TABLE]

For the case of the -8/2 rational surface in DIII-D discharge #158115, $\eta_{i}=1.9$ and $\lambda_{\theta\,i}=0.272$, which implies that

[TABLE]

in the nonlinear regime. The fact that the natural rotation of a nonlinear magnetic island chain is offset in the ion diamagnetic direction, rather than the electron diamagnetic direction, relative to the local ${\bf E}\times{\bf B}$ frame, leads to a much smaller prediction for the natural frequency in the nonlinear regime relative to that in the linear regime. This is the main reason why the predicted penetration threshold in the nonlinear regime is so much smaller than that in the linear regime. According to linear physics, we would expect mode penetration at the -8/2 surface to be triggered when $\omega_{\perp\,e}$ passes through zero. Moreover, once mode penetration has occurred, and the natural frequency becomes zero, we would expect $\omega_{\perp\,e}$ to be pinned to zero at the rational surface. On the other hand, according to nonlinear physics, we would expect mode penetration at the -8/2 surface to be triggered when some frequency offset from the ${\bf E}\times{\bf B}$ frequency in the ion diamagnetic direction passes through zero. Moreover, once mode penetration has occurred, and the natural frequency becomes zero, we would expect the offset frequency to be pinned to zero at the rational surface. In fact, experimental RMP-induced ELM suppression data from the DIII-D tokamak is not consistent with $\omega_{\perp\,e}$ being the trigger frequency [41], but rather some frequency offset from $\omega_{\perp\,e}$ in the ion diamagnetic direction, which constitutes strong evidence that mode penetration at the top of the pedestal in DIII-D RMP-induced ELM suppression experiments is not governed by linear physics. Incidentally, there is clear experimental evidence that the natural frequency of a nonlinear magnetic island chain is offset in the ion diamagnetic direction relative to the local ${\bf E}\times{\bf B}$ frame [42, 43].

Figure 8 shows a nonlinear simulation of the response of the plasma at the -11/2 rational surface in DIII-D discharge #158115 to an RMP whose magnitude is switched on at $t=0$, and then modulated cycloidally at a frequency of 1 Hz. The simulation data presented in Figure 8 is qualitatively similar the TM1 simulation data shown in Figure 2 of [15]. In particular, it is clear from Figure 8 that there is not enough plasma rotation at the -11/2 rational surface to enable the effective shielding of driven magnetic reconnection. In other words, $\hat{W}\sim\sqrt{b_{f}}\equiv\hat{W}_{f}$ at all times. Note, however, that $W>\delta_{\rm SC}$, which confirms that driven reconnection at the -11/2 rational surface in DIII-D discharge #158115 is governed by nonlinear physics.

6 Summary and Discussion

This paper investigates the plasma response to an externally generated, static, $n=2$, RMP, whose amplitude is modulated cycloidally at a frequency of 1 Hz, at two rational surfaces located in the pedestal of DIII-D H-mode discharge #158115 [13]. The first rational surface is the $m=-8$/$n=2$ surface, and lies at the top of the pedestal. The second is the $m=-11$/$n=2$ surface, and lies at the bottom of the pedestal. According to the nonlinear cylindrical reduced-MHD simulations of [15], mode penetration at the -11/2 surface is responsible for the so-called density pump-out. Moreover, mode penetration at the -8/2 surface is hypothesized to be responsible for RMP-induced ELM suppression. This paper examines mode penetration at the -8/2 and -11/2 surfaces using two analytic plasma response models. The first model is linear in nature, and the second nonlinear.

Our linear response model is based on the analysis of [20], in which the four-field model [44] is used to find all possible linear, two-fluid, drift-MHD, resonant plasma response regimes when a static RMP is applied to a large aspect-ratio tokamak plasma. We deduce that the particular response regime that is appropriate at both the -8/2 and -11/2 rational surfaces in DIII-D discharge #158115 is the so-called first semi-collisional regime (SCi). Incidentally, it has long been recognized that linear (tearing) layer physics in high-temperature tokamak plasmas is semi-collisional in nature, rather than collisional or collisionless [27, 45, 46]. Our linear response model does not incorporate the screening effect due to magnetic curvature that was discovered in [47]. However, this effect may be negated by parallel thermal transport [48]. More importantly, the curvature screening effect is a prediction of collisional layer physics, and the true layer physics in tokamak plasmas is semi-collisional. Indeed, [45] found that, in a semi-collisional layer, the effect of the perturbed bootstrap current is much greater in magnitude than, and opposed to, the effect of magnetic curvature. A much more serious deficiency in our linear layer model emanates from the fact that the layer width is similar to the ion sound radius. Given that the ion and electron temperatures in the pedestal of DIII-D discharge #158115 are almost equal, this implies that the layer width is also similar to the ion gyroradius. Unfortunately, the finite ion gyroradius width is not taken into account in the layer analysis of [20]. Moreover, it is known that, in situations in which the ion gyroradius is similar to, or exceeds, the linear layer width, finite ion orbit width effects can significantly modify the layer response [45, 49].

The linear response model adopted in this paper is augmented by plasma poloidal and toroidal equations of motion that determine how the quasi-linear electromagnetic locking torque that develops at the rational surface, in response to the applied RMP, modifies the local plasma rotation. The equations of motion take plasma perpendicular viscosity, neoclassical poloidal flow-damping, and neoclassical toroidal flow-damping into account. When the equations of motion are combined with the linear response model, a closed set of equations is obtained; these equations are solved numerically.

The linear response model is only valid when the width of the RMP-induced magnetic island chain at the rational surface falls below the linear layer width. In the opposite situation, in which the driven island width exceeds the linear layer width, the linear response model must be replaced by a nonlinear response model. It turns out that the appropriate nonlinear response model is, in many ways, simpler than the linear response model, given that it essentially consists of the Rutherford island width evolution equation [26] combined with the no-slip constraint [24].

The main conclusion of our linear simulations, which are described in Section 5.2, is that a linear response model is inapplicable at both the -8/2 and -11/2 rational surfaces in DIII-D discharge #158115. The problem is that the semi-collisional layer widths at the rational surfaces are so small (a few mm, in both cases) that rotational shielding is not sufficient to reduce the driven magnetic island widths below the layer widths.

Our nonlinear simulations, which are described in Section 5.3, give results that are similar to the experimental results described in [15]. At the -8/2 rational surface, driven magnetic reconnection is strongly screened by plasma rotation as long as the resonant component of the radial magnetic field remains below a certain threshold value. However, as soon as the radial field exceeds the threshold value, which is about 6 gauss, there is a sudden and irreversible breakdown of screening, accompanied by rapid shifts in the local plasma toroidal and poloidal angular velocities. On the other hand, at the -11/2 rational surface, the plasma rotation is not large enough to screen driven magnetic reconnection.

7 Future Work

There are a number of improvements that could be made to the analytic model described in this paper. Such improvements include; using a more accurate neoclassical model; taking into account the coupling of different rational surfaces via mode-penetration-induced changes in the plasma rotation, density, and temperature, profiles in the pedestal [15]; taking into account the coupling of different poloidal harmonics of the RMP due to toroidicity, the Shafranov shift, and flux surface shaping [50]; employing a more realistic plasma equilibrium; including island saturation terms [51, 52, 53], the perturbed bootstrap current [38], and the perturbed ion polarization current [54], in the Rutherford equation; and taking into account orbit-squeezing effects due to the strong shear in the radial electric field that is typically present in H-mode tokamak pedestals [55].

It is a well-known fact that RMP-induced ELM suppression is only effective when the value of $q_{95}$ (i.e., the safety-factor at the 95% magnetic flux-surface) lies in certain narrow windows [13, 56]. In future work, we intend to use an improved version of the analytic model of RMP-induced ELM suppression presented in this paper to investigate the dependence of the ELM suppression threshold on $q_{95}$. Such an investigation will inevitably entail thousands of simulations, and is only feasible with the type of highly-reduced analytic model described in this paper.

Acknowledgements

This research was funded by the U.S. Department of Energy under contract DE-FG02-04ER-54742. The author would like to thank R. Nazikian, Q.M. Hu, and C. Paz-Soldan for helpful discussions.