Calculating RF current condensation with self-consistent ray-tracing
Richard Nies, Allan H. Reiman, Eduardo Rodriguez, Nicola Bertelli, Nathaniel J. Fisch
CCalculating RF current condensation with self-consistent ray-tracing
R. Nies,
1, 2, a) A. H. Reiman,
1, 2, b) E. Rodriguez,
1, 2, c) N. Bertelli, d) and N. J. Fisch
1, 2, e) Department of Astrophysical Sciences, Princeton University, Princeton, NJ, 08543 Princeton Plasma Physics Laboratory, Princeton, NJ, 08540 (Dated: May 14, 2020)
By exploiting the nonlinear amplification of the power deposition of RF waves, current condensation promises newpathways to the stabilisation of magnetic islands. We present a numerical analysis of current condensation, couplinga geometrical optics treatment of wave propagation and damping to a thermal diffusion equation solver in the island.Taking into account the island geometry and relativistic damping, previous analytical theory can be made more preciseand specific scenarios can be realistically predicted. With this more precise description, bifurcations and associatedhysteresis effects could be obtained in an ITER-like scenario at realistic parameter values. Moreover, it is shown thatdynamically varying the RF wave launching angles can lead to hysteresis and help to avoid the nonlinear shadowingeffect.
I. INTRODUCTION
Reliable mitigation and avoidance of disruptions is criticalto the success of ITER and potential future tokamak powerplants. Sudden loss of plasma confinement poses a seri-ous threat to machine components through high heat loads,EM forces and runaway electrons . In the JET tokamakequipped with an ITER-like wall, 95% of natural disruptionsare preceded by magnetic islands , making their stabilisationan essential task.Magnetic islands are suppressed by generating a stabilisingresonant component of the magnetic field, as is typically doneby driving current at the island O-point. Current is generallydriven directly by RF waves , such as electron-cyclotron(EC) and lower-hybrid (LH) waves. Additionally, by deposit-ing power and thereby heating the plasma, RF waves canalso modify the ohmic current profile by decreasing the lo-cal resistivity. Thus, both RF heating and current drive canbe used to stabilise magnetic islands, as has been investi-gated theoretically and experimentally . For instance,stabilisation of magnetic islands with ECCD is planned inITER .Typically, both EC and LH waves are damped on fast su-perthermal electrons. Therefore, the damping rate stronglydepends on temperature through the electron population inthe tail of the distribution function. This high sensitivity ofpower deposition to temperature can result in a positive feed-back loop, where the magnetic island is heated by the RFwave, the elevated temperature leads to an increased powerdeposition, and so forth. This nonlinear effect, called currentcondensation , can lend further help in stabilising islands,as narrower power deposition and current profiles centred onthe island’s O-point can be achieved . Furthermore, cur-rent condensation can lead to bifurcations, where the islandtemperature would increase without bound, if not for other a) [email protected] b) [email protected] c) [email protected] d) [email protected] e) fi[email protected] limiting effects such as depletion of the wave power or stifftemperature gradients .The present work extends previous analytical studies with a numerical approach of current condensation. A geo-metrical optics treatment of wave propagation and damping iscoupled with a solver of the thermal diffusion equation in theisland geometry, as presented in Sec. II. The calculation iter-ates between solution of the thermal diffusion equation andcalculation of the power deposition along the ray trajectoriesin the presence of the perturbed temperature to obtain a self-consistent solution of the nonlinear coupled system.The geometrical optics approach allows for the inclusion ofrelativistic effects in the damping. We show in Sec. III thatthe previously developed theory of current condensation can be generalised to account for these relativistic effects, aswell as the island geometry. In Sec. IV, we obtain values andtrends in the bifurcation threshold consistent with previouswork. Furthermore, the same calculations show that a bifurca-tion and associated hysteresis in the island temperature couldbe obtained in ITER-like H-mode and L-mode scenarios, atrealistic values of input power, diffusion coefficient and islandtemperature perturbation. A constant diffusion coefficient washowever used, an approximation which holds for low temper-ature perturbations, below the ITG instability threshold. Anestimate suggests that this is indeed justified for the lowesttemperature perturbations we observed at a bifurcation, al-though more detailed calculations will need to be performedin the future.Finally, we show in Sec. V that a bifurcation and hystere-sis can be obtained by varying the RF wave launching angles,which could be a pathway for future experimental verifica-tion of current condensation. Additionally, the launching an-gles can be adjusted to circumvent the nonlinear shadowingeffect . II. COUPLING OF RAY-TRACING AND MAGNETICISLAND MODEL
The numerical approach presented below aims to simulatethe nonlinear dynamics of current condensation, yielding self-consistent temperature and power depositions. The newlydeveloped code OCCAMI (Of Current Condensation Amid a r X i v : . [ phy s i c s . p l a s m - ph ] M a y Magnetic Islands) couples the ray-tracing code GENRAY with a driven heat diffusion equation solver for the magneticisland. The ray-tracing computes the RF wave propagationand damping. The ensuing power deposition is used to solvethe steady-state diffusion equation to obtain the temperatureprofile in the island. The temperature is then given back as aninput to the ray-tracing code. This process is repeated untilconvergence in the island temperature is attained.The temperature and power deposition thus obtained areself-consistent and allow us to investigate current condensa-tion. This numerical treatment expands on previous analyticalwork where the initial power deposition was assumed to beconstant or exponentially decreasing around a peak , al-though the latter also incorporated self-consistent depletion ofthe wave energy. Furthermore, geometric effects inherent tomagnetic islands are now included in the heat diffusion equa-tion solver, whereas a slab model had previously been used tokeep the problem analytically tractable .We now present in more detail the coupling of the ray-tracing for wave propagation and absorption with the heat dif-fusion equation solver. A. Ray-tracing for wave propagation and power deposition
The code GENRAY simulates the propagation and ab-sorbtion of electromagnetic waves in the geometrical opticsapproximation. The coupling of the island model with theray-tracing calculations occurs through the temperature pro-file. We are therefore assuming that the perturbation in themagnetic field B r associated with a magnetic island is small,with a negligible impact on ray propagation and absorption.Axisymmetry of the plasma is lost in the presence of mag-netic islands, whence the island temperature profile becomesa function of not only a radial coordinate, but also the heli-cal angle ζ = θ − N / M ϕ , with the poloidal (toroidal) angles θ ( ϕ ) and poloidal (toroidal) mode numbers M ( N ). How-ever, the code GENRAY assumes axisymmetry, by requiring aone-dimensional temperature profile as input. We incorporatethe island-linked asymmetry through an effective temperatureprofile for each ray. This effective profile corresponds to thetemperature profile that the ray experiences as it propagatesthrough the plasma, as illustrated in Fig. 1. This assumes thatthe ray trajectory is not significantly altered by the change intemperature between iterations. In the case of multiple passes,i.e. when the ray trajectory traverses a given radius at multi-ple points in its trajectory, the effective profile is taken to bethat experienced by the ray in its last traversal of the island.This is a reasonable approximation when most of the poweris absorbed in a single pass, i.e. we assume the damping onall passes but the last must have been negligible. We will con-sider EC waves in this study, for which this is generally thecase.The solution for the temperature shown in Fig. 1 corre-sponds to a case where there is a single ray propagatingthrough a locked island. For a rotating island, it is necessaryto calculate the total power deposited in the island throughone rotation. This has been implemented for a fast rotating is- Figure 1. Effective Temperature profile input to GENRAY, as afunction of normalised radius ρ N = √ ψ t , with the toroidal flux ψ t .The upper plot shows the 2-dimensional temperature profile with a ( M = , N = ) magnetic island. The dotted lines are temperaturecontours, while the solid line represents the ray trajectory. The lowerplot shows different cuts in the upper plot: at ζ = , π (O-point), ζ = π / , π / land, where the diffusion time is much longer than the islandrotation time. Then, the total power deposition can be approx-imated by averaging the power deposition along multiple raytrajectories sampling the island at different phases. B. Heat diffusion equation solver in island geometry
The power deposition obtained from the ray-tracing calcu-lation is used to update the island temperature. Integrating thesteady-state diffusion equation once (a detailed derivation canbe found in appendix A), we obtain ∂ u ∂ σ = − P dep ( σ ) n χ ⊥ T s σ E ( σ ) − ( − σ ) K ( σ ) W M π r r R . (1)Here, u = ( T − T s ) / T s is the normalised island temperature,with the temperature at the separatrix T s , P dep ( σ ) is the powerdeposited within the island flux surface σ ( σ = , E ( σ ) and K ( σ ) are re-spectively the complete elliptical integrals of the first and sec-ond kind, W is the island width, r r is the minor radius at theresonant surface, R the tokamak major radius and M the is-land’s poloidal mode number. The temperature at the sepa-ratrix T s is assumed to be constant in our simulations. Thiscorresponds to the case where, for example, the EC poweris initially deposited radially inward from the island, at radii r < r r − W /
2, and then redirected outward to r ∼ r r . As the to-tal power deposited within the flux surface at r = r r + W / r ≥ r r + W / .More generally, the temperature perturbation at the separatrixmay be negligibly small, but this is not always the case.Furthermore, the perpendicular heat diffusion coefficient χ ⊥ is taken to be constant in this study. The interplay betweencurrent condensation effects and a variable heat diffusion co-efficient in the form of stiff gradients was investigated analyt-ically in Rodríguez, Reiman, and Fisch . A correspondingnumerical treatment with OCCAMI is left for future work.As it stands, the diffusion equation has been reduced to a1st order ODE (Eq. 1). It can thus be readily solved for theisland temperature profile, given the power deposited fromGENRAY and the boundary condition u ( σ = ) =
0. Thisis done numerically with a 4th order Runge-Kutta integrator.Note that the second boundary condition for the original 2ndorder diffusion equation, ∂ u / ∂ σ ( σ = ) =
0, was used whenintegrating the originally 2nd order diffusion equation, and isrequired for regularity.The obtained temperature is then fed back to GENRAYthrough an updated effective temperature profile and the raypropagation and power deposition are calculated anew. Thiscycle is repeated until convergence in the island temperatureis reached, i.e. when the relative change in the normalisedtemperature perturbation between iterations is below a giventhreshold ε (in this study, ε = · − was chosen). III. EFFECTS OF RELATIVISTIC DAMPING
In this section, we show how the sensitivity of damping totemperature in the classical case, w , can be generalised toaccount for relativistic damping effects, leading to the defi-nition of an effective w eff . Approximate formulas for the O1and X2 modes are presented. Spatial variation of the dampingwithin the island must also be taken into account, leading tothe introduction of an average w eff . To keep this study self-consistent and motivate our analysis of w eff , we first presenta brief summary of the linear theory of resonant wave damp-ing, and of the basic theory of relativistic effects in electroncyclotron wave damping. A. Linear theory of resonant wave damping
In addition to calculating the wave propagation, ray-tracingcodes such as GENRAY also generally provide the damp-ing coefficient of the wave, obtained from the anti-hermitianpart of the dielectric tensor. Formulas for the dielectric ten-sor in various regimes are obtained in the linear regime ofwave damping, i.e. assuming the distribution function tobe Maxwellian. However, this is only an approximation asthe distribution function is modified by the wave interactionthrough quasi-linear diffusion of particles in velocity space, making the absorption of RF wave energy a nonlinear process.Nevertheless, the linear theory is a suitable approximation inmultiple scenarios, as shown below, and we thus make use ofit in the analysis and simulations presented in this study.In certain cases, e.g. at low wave power, the electron dis-tribution remaining Maxwellian is a valid assumption. More-over, for EC waves, the linear theory was found to be valid andnearly independent of RF wave power, in the non-relativisticlimit and for diffusion in the perpendicular velocity v ⊥ only.The classical (non-relativistic) resonance condition betweenEC wave and electrons is k (cid:107) v (cid:107) = ω − n Ω , with n the harmonicnumber of the resonance, k (cid:107) the wavenumber parallel to themagnetic field, ω the wave frequency, Ω = eB / m e the cy-clotron frequency and v (cid:107) the parallel velocity of resonant elec-trons. Then, if the particles diffuse in v ⊥ only, they remain inresonance with the wave as the classical resonance conditionis independent of v ⊥ . The analysis will be more complicatedin the relativistic case, where the resonance condition dependson v ⊥ (see below, Eq. 2), or when diffusion of particles alsooccurs in the v (cid:107) direction.Evidence for the linear theory’s more general validity canbe found in a study comparing the current drive obtainedfrom linear theory with results of Fokker-Planck calculations,which solve the nonlinear damping problem, and with exper-iment: good agreement was shown for a large range of pa-rameters in DIII-D . The linear theory was also shown toyield accurate results for an ITER benchmark scenario in astudy comparing multiple ray-tracing codes (including GEN-RAY), which compute the damping from linear theory, withtwo Fokker-Planck solvers . Nevertheless, additional FokkerPlanck calculations need to be performed for a broader rangeof parameters spanning those investigated in our study. B. Relativistic resonance and damping
As shown in Sec. III A, the classical resonance conditionbetween EC wave and electrons is independent of v ⊥ . Thespatial damping rate is obtained by integrating over the dis-tribution function of resonant electrons in velocity space, andwill thus be proportional to the population of electrons withparallel velocity satisfying the resonance condition. In thelinear regime, i.e. when the distribution function in the par-allel velocity is Maxwellian, the spatial damping rate of ECwaves thus obeys α ∝ e − w , with the thermal velocity v T , and w = ( ω − n Ω ) / ( k (cid:107) v T ) (e.g. Swanson ).However, relativistic effects on the damping cannot be ne-glected in realistic scenarios. Indeed, the classical resonancecondition needs to be modified to take into account the rela-tivistic mass increase (e.g. Fidone, Granata, and Meyer ): ω − n Ω γ = k (cid:107) v (cid:107) . (2)In the relativistic case, the resonance follows an ellipse in ( v (cid:107) , v ⊥ ) space, due to the factor γ = ( − ( v (cid:107) + v ⊥ ) / c ) − / ,according to Eq. 2. In particular, this sets an important con-straint for resonance on the low-field side ( n Ω < ω ), n Ω ω ≥ (cid:113) − N (cid:107) , (3)where N (cid:107) = k (cid:107) c / ω is the parallel refractive index. As Ω ∝ B ∼ / R , Eq. 3 leads to a relativistic boundary, rendering part ofthe tokamak’s low-field side inaccessible to heating and cur-rent drive with EC waves. This is most apparent for low val-ues of the parallel refractive index, for which the relativisticboundary is close to the resonance Ω / ω =
1, and the wavetypically damps very strongly in a narrow spatial region.In the following, we present an approximate form of thedamping coefficient α for the O1 mode, as is appropriate forITER. The corresponding formulas for the X2 mode can befound in Appendix D. For the O1 mode, we assume absorptionof the wave’s L-polarisation to be negligible, and the wave’sR-polarisation to be generally small. Then the damping is duemainly to the electric field component along the backgroundmagnetic field. In that case, we can approximate the dampingcoefficient as α ≈ | E z / E | ω cN ε (cid:48)(cid:48) , (4)with the longitudinal polarisation | E z / E | , the refractive in-dex N ≈ (cid:113) − ω p / ω , the electron plasma frequency ω p = (cid:112) n e e / ( ε m e ) and electron density n e . The anti-hermitiancomponent ε (cid:48)(cid:48) of the dielectric tensor, derived by Fidone,Granata, and Meyer by integrating over the resonant ellipsein velocity space, is reproduced here: ε (cid:48)(cid:48) ≈ πω p Ω R / N / (cid:107) N ⊥ µ S (cid:16) − N (cid:107) (cid:17) / (cid:34) I / ( ξ ) (cid:32) + N (cid:107) Ω R ω (cid:33) − I / ( ξ ) (cid:18) ξ + N (cid:107) Ω R ω (cid:19) (cid:35) e µ (cid:32) − Ω / ω − N (cid:107) (cid:33) , (5)where I ν ( ξ ) are the ν -th order modified Bessel functions ofthe first kind, µ = m e c T , (6) R = (cid:115)(cid:18) Ω ω (cid:19) − + N (cid:107) , (7) ξ = N (cid:107) R µ − N (cid:107) , and (8) S = H (cid:32)(cid:18) Ω ω (cid:19) − + N (cid:107) (cid:33) (9)is a Heaviside function enforcing the relativistic constraint ofEq. 3. C. Sensitivity of relativistic damping to temperature
As seen above, the linear damping rate of EC waves sat-isfies α ∝ e − w in the classical limit, with w ∝ / T . This strong sensitivity of damping to temperature is essential forthe current condensation effect, with nonlinear effects becom-ing potentially relevant for w ∆ T / T (cid:38) . . The quantity w therefore provides a direct indicator of the sensitivity ofdamping to temperature in the classical limit.We are interested in obtaining the sensitivity of dampingto temperature taking into account relativistic effects in thedamping. We define an effective w eff as w = T ∂ T ( ln α ) , (10)such that w eff measures the sensitivity of damping to tem-perature, in analogy to the classical case. Indeed, Eq. 10can be viewed as a first order correction term in a Taylorexpansion of ln ( α ) . Therefore, for small temperature per-turbations, w ∆ T / T indicates the strength of nonlinear ef-fects, as the damping is approximately amplified by a factorexp ( w ∆ T / T ) . Finite temperature perturbations are treatedin Sec. III D in the regime w ∝ / T .Note that the damping satisfies α ∝ e − w only in the casewhere w ∝ / T , i.e. when w eff possesses the same tempera-ture dependency as w , as can be shown by integrating Eq. 10.However, w ∝ / T does not necessarily imply w eff = w .Using the definition of w eff in Eq. 10, the approximateform of the damping in Eq. 4 and assuming the longitudi-nal polarisation | E z / E | to be independent of temperature (coldplasma approximation), we obtain w ≈ T ∂ T ( ln ( ε (cid:48)(cid:48) )) = − ξ ∂ ξ ( ln ( ε (cid:48)(cid:48) )) . Combined with Eq. 5, this yields − w = + µ (cid:32) − Ω / ω − N (cid:107) (cid:33) + F ( ξ , a ) , (11)where a = N (cid:107) Ω R ω , (12) F ( ξ , a ) = ξ I / ( ξ ) (cid:16) − ξ + a ξ − a (cid:17) + I / ( ξ ) (cid:16) + a + ξ + a ξ (cid:17) I / ( ξ ) ( + a ) − I / ( ξ ) (cid:16) ξ + a (cid:17) . (13)The w eff of Eq. 11 reduces to w = ( ω − Ω ) / ( k (cid:107) v T ) in theclassical limit, consisting of N (cid:107) (cid:29) T / ( m e c ) and N (cid:107) (cid:29) | − ( Ω / ω ) | , as well as µ | − ( Ω / ω ) | (cid:29)
1, as the dampingcoefficient α also reduces to the classical formula in the samelimit (see Appendix C).Our ray-tracing calculations employ the more general ap-proximation of the dielectric tensor for a relativistic elec-tron plasma from Mazzucato, Fidone, and Granata , whichadds higher harmonic corrections to the treatment of Fidone,Granata, and Meyer . However, it is found that Eq. 11 agreeswell with a numerical evaluation of Eq. 10 in situations of in-terest for ITER, as shown in Appendix B.The temperature dependence of w eff defined in Eq. 11 isnontrivial, in contrast to the classical case where w ∝ / T . Itis however shown in Appendix C that w ∝ / T is a suitableapproximation in the limit ξ (cid:29)
1, and N (cid:107) ∼ | − ( Ω / ω ) | . If N (cid:107) (cid:29) | − ( Ω / ω ) | however, w ∝ / T can still hold when µ ( − ( Ω / ω )) (cid:29)
1. The first set of conditions are typicallysatisfied for reasonably high N (cid:107) and not too high temperatures(see Fig. 8), as the condition N (cid:107) (cid:29) | − ( Ω / ω ) | is very re-strictive. Therefore, multiple insights from the theory devel-oped in the case of classical damping remain valid in the rel-ativistic case, e.g. the temperature perturbations necessary toobtain a bifurcation, as shown in Sec. IV. D. Relativistic damping and current condensation
The w eff derived above is now connected to the current con-densation effect. We define the nonlinear amplification pa-rameter Θ NL as the logarithmic change in the damping rate α due to a change in temperature ∆ T from an unperturbed tem-perature T , Θ NL ≡ ln (cid:18) α ( T = T + ∆ T ) α ( T = T ) (cid:19) . (14)For small temperature perturbations, using Eq. 10, we obtain Θ NL = w ∆ T / T , a result valid for arbitrary forms of w eff .The case of finite temperature perturbations can be treated byassuming w ∝ / T is valid (see Appendix C). Then, Eq. 10can be integrated to obtain α = α e − w . Further, defining u ≡ ∆ T / T , Eq. 14 reduces to Θ NL = − w ( T = T ( + u )) + w ( T = T )= w ( T = T ) u + u . (15)As expected, this quantity reduces to Θ NL = w u for smalltemperature perturbations u (cid:28)
1, as was assumed in previouswork .The nonlinear amplification parameter Θ NL proves useful toascertain whether nonlinear effects like current condensationcan become relevant in a given scenario. Indeed, Θ NL (cid:38) . Θ NL → w eff ( T = T ) can be obtained fromthe wave damping in the unperturbed temperature profile, anapproximate lower bound on the temperature perturbation u necessary to observe nonlinear effects can be obtained frominverting Eq. 15 without performing the full nonlinear calcu-lation. This motivates the use of the approximated form inEq. 15 instead of inserting the full damping coefficient intoEq. 14, which could not be readily solved for the temperatureperturbation. Henceforth in this study, we will consider w eff to be evaluated at the initial unperturbed temperature T .Although w eff and Θ NL are useful local quantities, thedamping rate can vary significantly within a given magneticisland, especially for large island widths. Indeed, while wegenerally assume the initial island temperature to be flat, thequantities Ω / ω , N (cid:107) and thus also w eff will in general be func-tions of position within the island. Thus, a suitable island average of w eff must be found. We define w ≡ − ln (cid:0) (cid:104) exp (cid:0) − w (cid:1) (cid:105) σ ≤ (cid:1) , (16)where (cid:104) f (cid:105) σ ≤ is the mean value of f within the magnetic is-land ( σ ≤ α = α exp (cid:0) − w (cid:1) when w ∝ / T (see Appendix C for region of validity). IV. BIFURCATIONS AND HYSTERESIS IN ITER
Current condensation can lead to bifurcations, where thenonlinear amplification of temperature leads to a runaway ef-fect for the temperature and power deposited in the island.The temperature continues increasing , until another limit-ing physical mechanism comes into play, leading to satura-tion. For example, the wave may have deposited all of itspower , or the temperature increase might be limited by stifftemperature gradients .Depending on the case, the limiting effects might lead thetemperature to either smoothly transition to higher values, orexperience a discontinuous jump at the bifurcation . In thelatter case, hysteresis phenomena can be observed, as the jumpfrom low to high temperature will not occur at the same pa-rameter values as that from high to low temperature. Hys-teresis curves can be traced out e.g. by varying the RF inputpower (as considered in this section), the island width or the RF wave launching angles, as shown in Sec. V.It is shown in this section that bifurcations can be obtainedin ITER-like scenarios at realistic parameter values of tem-perature perturbation, diffusion coefficient and input power.Furthermore, values of the bifurcation threshold are shown tobe consistent with previous work. A. Simulations setup
The following simulations are based on ITER-like H- andL-mode scenarios, with temperature profiles shown in Fig. 2.A large island of width W N = . ρ N = √ ψ t , with the toroidal flux ψ t ) is introduced at the q = ρ N = . χ ⊥ = . s − . Such small values are justified as tur-bulent transport is reduced due to the flattened temperaturein the island region . This approximation will thus breakdown when the Ion Temperature Gradient (ITG) threshold is Figure 2. Temperature profiles for H and L-mode scenarios, beforeand after local flattening due to the magnetic island ( W = . , r r = . T s = . . et al. ). exceeded, i.e. when κ c ≤ − RT ∂ T ∂ r ≈ Ra u W N / , (17)with the major radius R = R + r r , major radius at the magneticaxis R = . a = . u = ∆ T ( σ = ) / T s ,and ITG threshold κ c ≈ . Then, effects of turbu-lent transport can be neglected when the temperature pertur-bations remain below u / W N (cid:46)
1, in which case the use of alow constant diffusion coefficient is justified. The maximallyallowed temperature perturbation can also be estimated fromthe temperature profile without island flattening. Assumingthe temperature gradients to be limited by ITG in this case,we can estimate the ITG threshold to be exceeded when theisland temperature reaches the temperature of the profile with-out island. For the H-mode profiles in Fig. 2, this would allowtemperature perturbations up to u ∼ .
25, which is consistentwith the previous estimate u ∼ W N = . . Furthermore, turbulence wasfound to enhance the transport of fast electrons accelerated byEC waves in regimes where it would not greatly impact thatof bulk electrons . Thus, one aim of this section is to ob-tain bifurcations and hysteresis behaviour at low temperatureperturbations, below the ITG threshold.The simulations shown below were obtained by a coarsescan in the parameter space of launcher position (from anupper launcher case to halfway between upper and equato- rial launcher) and launching angles (poloidal launching an-gle 112 ◦ ≤ α ≤ ◦ , measured from positive ˆ Z , and toroidallaunching angle 34 ◦ ≤ β ≤ ◦ , measured from negative ˆ R through launcher, in steps of size 1 ◦ ). Those simulations notdisplaying a bifurcation were discarded. To trace the hystere-sis curve, the power was gradually increased up to 20 MW,the maximal EC power available in ITER, and decreased backto low powers. Close to the bifurcation, the island tempera-ture is very sensitive to the power deposited, whence a smallstep size of ∼
20 W was used. The relative error in the islandtemperature at the bifurcation can be estimated as being of thesame order as the relative change in u in the last step beforethe bifurcation, which is maximally 0 .
4% in the simulationsshown here.Furthermore, the cases where the relativistic boundary(Eq. 3) is located within the island were excluded. In suchcases, the region where power can be deposited in the islandis shrunk, such that higher temperature perturbations are nec-essary to observe a bifurcation for a given value of w eff . Inprevious work, the cases of deposition starting at the islandcentre and at the island edge were considered .The island’s phase is locked such that the ray goes throughits O-point. Adjustment of a locked island’s phase to depositpower at the island’s O-point has been achieved in DIII-D withexternal magnetic perturbations .The EC wave propagation and damping were obtainedfrom GENRAY , using a cold dispersion relation for thewave propagation and Mazzucato, Fidone, and Granata ’sapproximation of the dielectric tensor for a relativistic elec-tron plasma for the wave damping. B. Bifurcation threshold in ITER
The observed values of w (as defined in Eq. 16) and ofthe normalised temperature perturbation in the island centreat the bifurcation u B , are shown in Fig. 3. In the limit of verysmall w , the damping is already strong in the linear regime,such that nonlinear effects cannot help to focus or draw in thepower deposition, hence no bifurcation is observed. Whengoing to very large w , too little power is deposited in theisland for nonlinear effects to be relevant. Thus, an interme-diate region where bifurcations can be observed is found, e.g.5 . (cid:46) w (cid:46) . w can be obtained in the L-mode scenario(up to w ≈ T s is lower (Fig. 2), yielding higher effective powers in the RHSof the diffusion equation (Eq. 1). Analogously to Rodríguez,Reiman, and Fisch , we define a normalised power densityas P ≡ PW w / ( V island nT s χ ⊥ ) , with the input wave power P , average density n , separatrix temperature T s and island vol-ume V island = π W r r R / M . Higher effective powers can thusbe achieved e.g. by reducing the island temperature and den-sity.As can be readily seen in Fig. 3, there is a strong correlationbetween the temperature perturbation at the bifurcation u B andthe w value. The nonlinear amplification parameter at the Figure 3. Observed values of w and of the normalised temperatureperturbation in the island centre at the bifurcation u B for ITER-likeH-mode and L-mode scenarios. Increasing w helps to obtain bifur-cations at lower temperature perturbations. bifurcation Θ BNL ≡ w u B / ( + u B ) (18)is shown in Fig. 4 as a function of normalised power P . Mostof the data points are in the region Θ BNL ≈ . − .
0, with atrend of a small decrease with increasing input power P atthe bifurcation. Due to the near constancy of Θ BNL , it can beeffectively used as a threshold parameter, below which no bi-furcation can be obtained. It is especially useful as w isobtained from a simple ray tracing calculation in the unper-turbed temperature profile. Therefore, a lower bound on thetemperature perturbation necessary to observe a potential bi-furcation can be obtained from solving Eq. 18 for u B , withoutsolving the full nonlinear problem.A small amount of spread in the value of Θ BNL at a given P can be observed. Some deviation is not surprising, as the be-haviour for strong variations of w eff within the island may notbe fully captured by the averaged w eff (Eq. 16), and the tem-perature dependence of w ∝ / T assumed in the definitionof Θ NL (Eq. 14) is only approximate.The values of Θ BNL agree well with previous results. Using aconstant power deposition and no wave depletion, Reiman andFisch find that a bifurcation occurs when Θ BNL = w u B ≈ . Θ BNL = w u B ≈ . find a bifurcationwhen Θ BNL = w u B ≈ . − .
5, for deposition starting at theisland edge. In comparison, the values of Θ BNL ≈ . − . . More-over, the tendency of Θ BNL to slowly decrease with increasinginput power P observed in Fig. 4 is also consistent with pre-vious analysis .Even though no detailed optimisation was performed in the Figure 4. Nonlinear amplification parameter at the bifurcation Θ BNL = w u B / ( + u B ) as a function of the normalised wave input powerdensity P at the bifurcation. The values of Θ BNL ≈ . − . P . simulations shown in this study, it can be seen from Fig. 3that bifurcations were obtained at temperature perturbationsdown to u = .
24 and 0 .
19 for the H-mode and L-mode tem-perature profiles, respectively. This suggests that bifurcationscould be obtained before the ITG threshold is exceeded at u ∼ . − .
25 (Eq. 17) for realistic parameter values. TheITG threshold criterion was only roughly estimated, however;and more realistic calculations which include stiff-gradient ef-fects self-consistently will thus need to be undertaken in thefuture.The lowest u values at bifurcation were obtained for up-per launcher cases, close to ITER’s planned upper launcherposition. In these cases, w eff tends to be approximately flatinside the island as the ray propagation in the poloidal planeoccurs mostly in the ˆ Z -direction, and the resonance is thusapproached slowly. To put the necessary temperature pertur-bations into perspective, values of u ≥ . W / a ∼ . .Large toroidal launching angles ( β ≥ ◦ ) were chosen inthis study to obtain higher w eff , yielding stronger nonlinear ef-fects, although ITER’s upper launcher is planned to operate ata smaller toroidal launching angle of 20 ◦ . Moreover, singlerays were used in this study, instead of gaussian beam pro-files. A detailed investigation of current condensation effectsfor the planned ITER upper launcher steering mirrors, also in-corporating gaussian beams represented by multiple rays, aswell as stiff gradient effects, is left for future work. V. DYNAMIC VARIATION OF RF WAVE LAUNCHINGANGLES
Current condensation can lead to hysteresis as elevated is-land temperatures can draw in and maintain the power deposi-tion close to the island centre, instead of having it be depositedcloser to the island edge or even outside of the island. A hys-
Figure 5. Hysteresis in normalised temperature at island centre u from variation of the poloidal launching angle α . Circle data pointsshow the ascending part of the hysteresis curve, going from the lowerto the upper branch (decreasing α ), while crosses show the descend-ing part, going from the upper branch to the lower branch (increasing α ). teresis curve can be traced by e.g. varying the wave poweror island width . Another way to obtain hysteresis is todynamically vary the toroidal or poloidal launching angle, thesecond of which we will demonstrate in this section.Furthermore, we will show that the dynamic variation ofthe poloidal launching angle can help to circumvent the shad-owing effect, a nonlinear inhibition effect . At high islandtemperatures, a significant fraction of the wave power may bedamped at the island edge before the wave reaches the islandcentre. This not only leads to reduced island temperaturesbut also to possibly destabilising currents driven close to theisland separatrix. One way to bring the power deposition to-ward the island O-point despite the shadowing effect involvespulsing the input wave power . We will show in this sec-tion that the shadowing effect can also be avoided by varyingthe poloidal launching angle, such that damping at the islandedge is reduced. Note that similar results can be achieved byvarying the toroidal launching angle.We again investigate the ITER-like H-mode case of Sec. IV,with a launcher situated close to ITER’s planned upperlauncher position , R = Z = . / W N = . ρ N ), aconstant diffusion coefficient χ ⊥ = . s − and 20 MW ofEC wave power are considered. The toroidal launching angleis held fixed at β = ◦ , while the poloidal launching angle isvaried from α = . ◦ down to 138 ◦ , and back up, in smallsteps of 0 . ◦ to accurately trace the hysteresis curve. Theresulting hysteresis curve for the island temperature is shownin Fig. 5.The power deposition at several points along the hysteresiscurve is shown in Fig. 6. Initially, at high poloidal launch-ing angles, little to no power is deposited inside the island(e.g. α = . ◦ in Fig. 6). The poloidal launching angleis decreased, moving the power deposition into the island, un-til a bifurcation is reached at α ≈ . ◦ . Increasing the Figure 6. Evolution of power deposition profile during hysteresisfrom variation of the poloidal launching angle. Vertical gray linesindicate the magnetic island edges. At first, almost all power is de-posited outside of the island (solid curve, α = . ◦ ). The powerdeposition is then moved into the island by decreasing the poloidallaunching angle, until enough power is deposited for a bifurcation tooccur at α ≈ . ◦ (dashed to dashdotted curves). Increasing thepoloidal launching angle back to α = . ◦ (dotted curve), hys-teresis behaviour is displayed, as most of the power is still depositedwithin the island. poloidal launching angle back to α = . ◦ , the power isstill deposited in the island, i.e. the system displays hysteresisbehaviour.The power deposition after the transition to the upperbranch at α = . ◦ is not centred on the island’s O-point ( ρ N ≈ . α = . ◦ ). This results in a higher central island temper-ature, as can be seen in Fig. 5.The hysteresis in the poloidal launching angle shown abovecould be preferable to that in the EC wave power for an ex-perimental investigation of current condensation. Comparedto the hysteresis in power, it allows operation at the maximumEC power, which was shown in Fig. 4 to lead to bifurcationsat lower values of Θ BNL ≡ w u B / ( + u B ) , i.e. bifurcationscould be obtained at lower temperatures. In the case shownhere, the temperature at the bifurcation is small for the H-mode case, u B = .
26, although no significant optimisation orlarge parameter scans were performed.Note that the upper branch solution of Fig. 5 will undoubt-edly be modified by stiff gradient effects, due to the largetemperature perturbations ( u >
4) that will trigger ITG in-stabilities. However, the mechanism to circumvent shadow-ing presented here is more generally valid. Current conden-sation including stiff gradient effects has been investigatedanalytically , with a corresponding numerical treatment leftfor future work. VI. CONCLUSION
In this study, we presented a numerical treatment of cur-rent condensation effects, coupling a ray-tracing code for thewave propagation and damping, with a heat diffusion equationsolver to obtain the temperature in the magnetic island. Thisallows to investigate current condensation in realistic scenar-ios, in particular including the island geometry and relativisticeffects in the damping. These were shown to lead to a general-isation of the bifurcation threshold identified in previous ana-lytical work . Furthermore, bifurcations were obtained forrealistic parameter values ( P =
20 MW, χ ⊥ = . s − ) inITER-like H- and L-mode scenarios, at low temperature per-turbations u ∼ .
2. Stiff gradient effects could be negligibleat such low temperature perturbations; more detailed calcula-tions which self-consistently include stiff gradient effects arehowever needed to demonstrate this. Finally, we showed thatdynamically varying the poloidal launching angle can lead toa current condensation induced hysteresis and can remedy tothe nonlinear inhibition brought about by the shadowing ef-fect. Current condensation could enable the stabilisation oflarge magnetic islands, leading to improved disruption avoid-ance. Therefore, an optimisation study for ITER, includingstiff gradient effects, the planned ITER launcher position andthe use of multiple rays to represent gaussian beams, is inpreparation.
ACKNOWLEDGMENTS
The simulations shown in this study were run on the PPPLresearch cluster. This work was supported by U.S. DOE DE-AC02-09CH11466 and de-sc0016072. The data that supportthe findings of this study are available from the correspondingauthor upon reasonable request.
Appendix A: Derivation of diffusion equation in islandgeometry
We assume the magnetic island to be symmetric in the ra-dial coordinate r , an approximation for narrow islands in alarge aspect ratio tokamak of circular cross-section. We de-fine the island flux coordinate σ , ranging from 0 in the centreto 1 at the separatrix, as r − r r = ± W (cid:113) σ − sin ( M ζ / ) , (A1)with the normalised radius at the resonant surface r r , the he-lical angle ζ = θ − N ϕ / M ∈ [ − π / M , π / M ) , the poloidal(toroidal) angles θ ( ϕ ) and poloidal (toroidal) mode numbers M ( N ). A new angular coordinate η ∈ [ − π /, π ) is defined assin ( M ζ / ) = σ sin ( η ) . (A2) Using the definition of η , Eq. A1 can be rewritten as r − r r = W σ cos η . (A3)We now investigate the steady-state diffusion equation in theisland, ∇ · ( n χ ∇ T ) = − p , (A4)with χ the heat diffusion coefficient, n the plasma density, and p the power density. Integrating over the island volume up tothe flux surface σ , (cid:90) π − π d η (cid:90) π d ϕ (cid:90) σ d σ (cid:48) J ∇ · ( n χ ∇ T ) (A5) = − (cid:90) π − π d η (cid:90) π d ϕ (cid:90) σ d σ (cid:48) J p , where, in a large aspect ratio approximation, the Jacobian J isgiven by J − = ∇ σ · ∇ η × ∇ ϕ ≈ MW r r R (cid:112) − σ sin ησ . (A6)Then, assuming parallel diffusion to be significantly strongerthan perpendicular diffusion, χ ⊥ (cid:28) χ (cid:107) , the temperature canbe assumed to be equilibrated on flux surfaces and thus be-comes a function of σ only, T = T ( σ ) . Defining P dep ( σ ) asthe triple integral on the right hand side, which represents thetotal power deposited inside the flux surface σ , and taking thedensity n and cross-field thermal diffusivity χ ⊥ to be constant,we obtain ∂ T ∂ σ · π (cid:90) π − π d η J | ∇ σ | = − P dep ( σ ) n χ ⊥ . (A7)Evaluating the integral on the left hand side for large aspectratios, we obtain (cid:90) π − π d η J | ∇ σ | ≈ R r r WM ( W / ) E ( σ ) − (cid:0) − σ (cid:1) K ( σ ) σ . (A8)Combining Eqs. A7 and A8 results in the diffusion equationin the island geometry, ∂ u ∂ σ = − P dep ( σ ) n χ ⊥ T s σ E ( σ ) − ( − σ ) K ( σ ) ( W / ) V island , (A9)where V island = π R r r W / M is the island volume. The ex-pression in Eq. A9 is equivalent to Eq. 1, presented in themain text. The term containing elliptic integrals in Eq. A9is an island geometric term. Previous studies of currentcondensation used a slab model of the island, for whichthe diffusion equation reduces to ∂ u ∂ σ = − P dep ( σ ) n χ ⊥ T s ( W / ) V island . (A10)The equations Eq. A9 and A10 are very similar, but for theadded island geometric term in the former.0 Figure 7. Comparison of w eff obtained from Eq. 11 (Analytical) andfrom numerical evaluation of GENRAY damping coefficient (Nu-merical). Two cases are considered, at low toroidal launching an-gles ( β = ◦ ), with N (cid:107) ∼ .
3, and higher toroidal launching angles( β = ◦ ), with N (cid:107) ∼ .
5. The w values are normalised by factors15 .
62 and 2 .
94 for β = ◦ and β = ◦ , respectively. Appendix B: Validity of w eff formula We compare the formula for w eff in Eq. 11 with a numeri-cal finite differences evaluation of Eq. 10 using the dampingcoefficient from the ray-tracing code GENRAY. The dampingformula by Mazzucato, Fidone, and Granata is used in theray-tracing calculations.The ITER-like H-mode profile of Sec. IV is used, and twocases with differing launching angles are considered, to ob-tain the case of low N (cid:107) ∼ . N (cid:107) ∼ . T ≈ − w values have been normalised by factors 15 .
62 and 2 .
94 for β = ◦ and β = ◦ , respectively. The agreement is excel-lent, with Eq. 11 only slightly overestimating the value of w eff .Note that the relativistic constraint (Eq. 3) is very prominentfor the low N (cid:107) ∼ . w eff . Appendix C: Temperature dependence and classical limit of w eff The temperature dependence of w eff in Eq. 11 is compli-cated due to the modified Bessel functions in the F ( ξ , a ) term(Eq. 13). In the remainder of this Appendix, we will show that w ∝ / T , as was assumed e.g. in Eq. 15, is a valid approx-imation in the limit of ξ (cid:29) | − ( Ω / ω ) | ∼ N (cid:107) .If | − ( Ω / ω ) | (cid:28) N (cid:107) , the 1 / T proportionality can still holdprovided that µ | − ( Ω / ω ) | (cid:29) ξ (cid:29)
1, which is sensible as µ (cid:29)
1, even for thermonucleartemperatures. However, reasonably large N (cid:107) and R also haveto be assumed for ξ (cid:29) N (cid:107) . Indeed, low N (cid:107) values willgenerally display significant damping close to the relativisticboundary (Eq. 3), where R → ξ (cid:29)
1, the modified Bessel functions can beapproximated as I ν ( ξ (cid:29) ) ≈ e ξ (cid:112) πξ ∞ ∑ k = ( − ) k b k ( ν ) ξ k , (C1) b k ( ν ) = ( ν − )( ν − ) ... ( ν − ( k − ) ) k !8 k . Then, Eq. 13 reduces to F ( ξ (cid:29) , a ) ≈ (C2) ξ ( a − ) − ( a − )( a − ) ξ + a − a + ξ + a − ξ ( a − ) − ξ ( a − )( a − ) + ξ ( − a ) . Consider the parameter ε = ( a − ) . The limit ε → Y →
1, where Y ≡ Ω / ω .From Eq. 12, ε = Y (cid:114) − Y − N (cid:107) − , (C3)such that ε (cid:28) | Y − | (cid:28) N (cid:107) . This condition ispart of the classical limit. Note that although the relativisticboundary (Eq. 3) constrains the damping to occur for Y ≥ − N (cid:107) on the tokamak low field side, the condition | Y − | (cid:28) N (cid:107) is a much stronger constraint.First, consider the case where | Y − | ∼ N (cid:107) . Then, F ( ξ , a ) ≈ ξ and, from Eq. 11, − w ≈ + µ (cid:32) − Y − N (cid:107) + N (cid:107) R − N (cid:107) (cid:33) . (C4)In the regime where | Y − | ∼ N (cid:107) and µ (cid:29)
1, the first term isnegligible, such that w ∝ µ ∝ / T , as desired.We now treat separately the regime | Y − | (cid:28) N (cid:107) . First,rewrite Eq. 11 as a function of ξ and a , − w = + a ξ N (cid:107) (cid:32)(cid:114)(cid:16) − N (cid:107) (cid:17) (cid:16) − N (cid:107) / a (cid:17) − (cid:33) + F ( ξ , a ) . (C5)Furthermore, Eq. C2 can be rewritten as F ( ξ , a = + ε ) ≈ ξ − ( ξ ε ) + ( ξ ε ) + ( ξ ε ) + ( ξ ε ) + . (C6)Inserting this expression into Eq. C5, we obtain − w ≈ ξ ε (cid:16) N (cid:107) − (cid:17) + ( ξ ε ) − ( ξ ε ) − ( ξ ε ) + ( ξ ε ) + . (C7)1Using | − Y | (cid:28) N (cid:107) and defining δ = − Y , Eq. C3 yields ξ ε ≈ µδ and ε ≈ δ ( − N (cid:107) ) / N (cid:107) . Then, Eq. C7 reduces to w ≈ µδ N (cid:107) − ( µδ ) − ( µδ ) − ( µδ ) + ( µδ ) + . (C8)The classical limit can then be obtained in the ordering µδ (cid:29)
1, such that w ≈ µδ N (cid:107) − = w . (C9)This indeed corresponds to the exact classical limit, obtainedby inserting α ∝ µ − / exp ( µδ / ( N (cid:107) )) into Eq. 10. Notethat we previously omitted the weak µ − / dependence of theclassical damping, assuming that the exponential term woulddominate, an invalid assumption when the resonance is ap-proached, δ →
0. In this limit, ( ω − Ω ) / ( k (cid:107) v T ) →
0, and thestrong sensitivity of damping to temperature is lost as damp-ing occurs on the bulk thermal electrons. From Eq. C9, the w ∝ / T proportionality holds strictly in the classical limitwhen the first term dominates, i.e. µδ / N (cid:107) (cid:29) µδ (cid:29) . There, only the ε (cid:48)(cid:48) component of the dielectric tensor was considered, which re-duces to the classical damping with only the limits µ (cid:29) N (cid:107) (cid:29) | − Y | . We are however interested in the ε (cid:48)(cid:48) com-ponent, for which µ | − Y | (cid:29) w T , as shown inFig. 8 as a function of temperature, for different values of N (cid:107) and Y = Ω / ω = − δ . Constancy of w T indicates that the w ∝ / T proportionality holds.First, consider the diamonds in Fig. 8, representing the tem-perature values for which ξ =
5. It can be seen that forlarger temperatures, the ξ (cid:29) w ∝ / T will break down, with w even going to nega-tive values. Away from the resonance, more specifically when | − Y | ∼ N (cid:107) , we should expect ξ (cid:29) N (cid:107) = . Y = .
96, for which | − Y | / N (cid:107) = .
87, or N (cid:107) = . Y = .
85, for which | − Y | / N (cid:107) = . w T is roughly constant to the left of the diamonds, where ξ (cid:29) ξ (cid:29) N (cid:107) values.Secondly, consider the circles in Fig. 8, representing thepoint on each curve for which µδ =
5. For those cases where ξ (cid:29) | − Y | (cid:28) N (cid:107) , a further requirement µδ (cid:29) w eff to the classical limit. If µδ (cid:29) w eff of Eq. C8 will have a more complicatedtemperature dependence. Indeed, for small δ = − Y (e.g. N (cid:107) = . Y = .
99, for which | − Y | / N (cid:107) = . Figure 8. Change of w T as a function of temperature, for twovalues N (cid:107) = . , . Y = Ω / ω . Circles indicatetemperatures for which µ ( − Y ) =
5, while diamonds correspond totemperatures for which ξ = a regime where w T is constant, to one where the behaviouris more complicated.Summarising, the assumption w ∝ / T thus strictly holdsfor ξ (cid:29)
1, and when simultaneously | − Y | ∼ N (cid:107) . If | − Y | (cid:28) N (cid:107) however, w ∝ / T can still hold for µ | − Y | (cid:29) µ ( − Y ) (cid:29) N (cid:107) ). In the simulationsshown in Sec. IV, typical values are µ ∼ , N (cid:107) ∼ . δ ∼ .
1. Then ξ ∼
80 and | − Y | / N (cid:107) ∼ .
75, and we can assume w ∝ / T to be valid, according to Eq. C4. Appendix D: Formulas for w eff in the X2 mode Whereas the O1 mode considered in Sec. III is most rele-vant for ITER, the X2-mode is more relevant for several exist-ing tokamaks, like DIII-D or AUG. Therefore, we repeat herethe procedure in Sec. III and derive a formula for w eff for theX2-mode.Just like for the O1 mode, the components of the dielectrictensor for a relativistic electron plasma can be obtained fromFidone, Granata, and Meyer . The first diagonal componentis reproduced here: ε (cid:48)(cid:48) = πω p Ω (cid:18) R N (cid:107) (cid:19) / N ⊥ S (cid:113) − N (cid:107) I / ( ξ ) e µ (cid:32) − Ω / ω − N (cid:107) (cid:33) , (D1)2with R = (cid:115)(cid:18) Ω ω (cid:19) − + N (cid:107) , (D2) ξ = N (cid:107) R µ − N (cid:107) , (D3) S = H (cid:32)(cid:18) Ω ω (cid:19) − + N (cid:107) , (cid:33) (D4)and the other symbols previously defined for Eq. 5. The damp-ing coefficient can be approximated as α ≈ ω cN ε (cid:48)(cid:48) . (D5)Then, using Eq. 10, w ≈ T ∂ T ( ln ( ε (cid:48)(cid:48) )) in conjunction withEq. D1, we obtain − w = µ (cid:32) − Ω / ω − N (cid:107) (cid:33) + ξ I / ( ξ ) I / ( ξ ) − . (D6) REFERENCES A. H. Boozer, “Theory of tokamak disruptions,” Physics of Plasmas ,058101 (2012). A. H. Boozer, “Magnetic surface loss and electron runaway,” PlasmaPhysics and Controlled Fusion , 024002 (2019). M. Lehnen, K. Aleynikova, P. B. Aleynikov, D. J. Campbell, P. Drewelow,N. W. Eidietis, Y. Gasparyan, R. S. Granetz, Y. Gribov, N. Hartmann, E. M.Hollmann, V. A. Izzo, S. Jachmich, S. H. Kim, M. Koˇcan, H. R. Koslowski,D. Kovalenko, U. Kruezi, A. Loarte, S. Maruyama, G. F. Matthews, P. B.Parks, G. Pautasso, R. A. Pitts, C. Reux, V. Riccardo, R. Roccella, J. A.Snipes, A. J. Thornton, and P. C. de Vries, “Disruptions in ITER andstrategies for their control and mitigation,” Journal of Nuclear MaterialsPLASMA-SURFACE INTERACTIONS 21, , 39–48 (2015). S. N. Gerasimov, P. Abreu, I. Carvalho, R. Henriques, G. Artaserse, P. Bu-ratti, M. Baruzzo, I. Coffey, E. De la Luna, T. Hender, R. Felton, P. Lomas,P. McCullen, M. Maslov, F. Rimini, W. Schippers, G. Szepesi, U. Kruezi,M. Tsalas, S. Moradi, L. Piron, and L. Zakharov, “Overview of disrup-tions with JET-ILW,” in (2018). N. J. Fisch, “Confining a Tokamak Plasma with rf-Driven Currents,” Phys-ical Review Letters , 873–876 (1978). N. J. Fisch and A. H. Boozer, “Creating an Asymmetric Plasma Resistivitywith Waves,” Physical Review Letters , 720–722 (1980). N. J. Fisch, “Current generation with low-frequency waves,” Physics of Flu-ids , 27 (1981). C. Karney and N. Fisch, “Currents driven by electron cyclotron waves,”Nuclear Fusion , 1549–1557 (1981). C. F. F. Karney and N. J. Fisch, “Efficiency of current drive by fast waves,”The Physics of Fluids , 116–126 (1985). C. F. F. Karney, N. J. Fisch, and F. C. Jobes, “Comparison of the theory andthe practice of lower-hybrid current drive,” Physical Review A , 2554–2556 (1985). P. T. Bonoli and R. C. Englade, “Simulation model for lower hybrid currentdrive,” Physics of Fluids , 2937–2950 (1986). N. J. Fisch, “Theory of current drive in plasmas,” Reviews of ModernPhysics , 175–234 (1987). C. F. F. Karney, N. J. Fisch, and A. H. Reiman, “Green’s function for RF-driven current in a toroidal plasma,” AIP Conference Proceedings , 430(1989). J. Decker, “ECCD for Advanced Tokamak Operations Fisch-Boozer versusOhkawa Methods,” AIP Conference Proceedings , 447–454 (2003). Y. R. Lin-Liu, V. S. Chan, and R. Prater, “Electron cyclotron current driveefficiency in general tokamak geometry,” Physics of Plasmas , 4064–4071 (2003). R. Prater, “Heating and current drive by electron cyclotron waves,” Physicsof Plasmas , 2349–2376 (2004). V. Chan and G. Guest, “Stabilization of tearing modes in tokamaks usingelectron cyclotron heating,” Nuclear Fusion , 272–276 (1982). A. H. Reiman, “Suppression of magnetic islands by rf driven currents,”Physics of Fluids , 1338 (1983). C. C. Hegna and J. D. Callen, “On the stabilization of neoclassical mag-netohydrodynamic tearing modes using localized current drive or heating,”Physics of Plasmas , 2940–2946 (1997). Q. Yu, S. Günter, G. Giruzzi, K. Lackner, and M. Zabiego, “Modeling ofthe stabilization of neoclassical tearing modes by localized radio frequencycurrent drive,” Physics of Plasmas , 312–322 (1999). D. De Lazzari and E. Westerhof, “On the merits of heating and current drivefor tearing mode stabilization,” Nuclear Fusion , 075002 (2009). N. Bertelli, D. De Lazzari, and E. Westerhof, “Requirements on localizedcurrent drive for the suppression of neoclassical tearing modes,” NuclearFusion , 103007 (2011). J. Pratt, G. T. A. Huijsmans, and E. Westerhof, “Early evolution of electroncyclotron driven current during suppression of tearing modes in a circulartokamak,” Physics of Plasmas , 102507 (2016). K. Hoshino, M. Mori, T. Yamamoto, H. Tamai, T. Shoji, Y. Miura,H. Aikawa, S. Kasai, T. Kawakami, H. Kawashima, M. Maeno, T. Mat-suda, K. Oasa, K. Odajima, H. Ogawa, T. Ogawa, T. Seike, T. Shi-ina, K. Uehara, T. Yamauchi, N. Suzuki, and H. Maeda, “Avoidance of${\mathit{q}}_{\mathit{a}}$=3 disruption by electron cyclotron heatingin the JFT-2M tokamak,” Physical Review Letters , 2208–2211 (1992). D. Kislov, V. Alikaev, Y. Esipchuk, A. Kakurtn, A. Kislov, D. Martynov,G. Notkin, K. Razumova, A. Sushkov, and V. Volkov, “The m=2, n=1 modesuppression by ECRH on the T-10 tokamak,” Nuclear Fusion , 339–350(1997). G. Gantenbein, H. Zohm, G. Giruzzi, S. Günter, F. Leuterer, M. Maraschek,J. Meskat, Q. Yu, ASDEX Upgrade Team, and ECRH-Group (AUG),“Complete Suppression of Neoclassical Tearing Modes with Current Driveat the Electron-Cyclotron-Resonance Frequency in ASDEX Upgrade Toka-mak,” Physical Review Letters , 1242–1245 (2000). A. Isayama, Y. Kamada, S. Ide, K. Hamamatsu, T. Oikawa, T. Suzuki,Y. Neyatani, T. Ozeki, Y. Ikeda, K. Kajiwara, and the JT-60 team, “Com-plete stabilization of a tearing mode in steady state high- β p H-mode dis-charges by the first harmonic electron cyclotron heating/current drive onJT-60U,” Plasma Physics and Controlled Fusion , L37–L45 (2000). H. Zohm, G. Gantenbein, A. Gude, S. Günter, F. Leuterer, M. Maraschek,J. P. Meskat, W. Suttrop, Q. Yu, A. U. Team, and E. Group (AUG), “Thephysics of neoclassical tearing modes and their stabilization by ECCD inASDEX Upgrade,” Nuclear Fusion , 197–202 (2001). H. Zohm, G. Gantenbein, A. Gude, S. Günter, F. Leuterer, M. Maraschek,J. Meskat, W. Suttrop, and Q. Yu, “Neoclassical tearing modes and theirstabilization by electron cyclotron current drive in ASDEX Upgrade,”Physics of Plasmas , 2009–2016 (2001). R. J. L. Haye, S. Günter, D. A. Humphreys, J. Lohr, T. C. Luce, M. E.Maraschek, C. C. Petty, R. Prater, J. T. Scoville, and E. J. Strait, “Control ofneoclassical tearing modes in DIII–D,” Physics of Plasmas , 2051 (2002). The TEXTOR team, I. G. J. Classen, E. Westerhof, C. W. Domier, A. J. H.Donné, R. J. E. Jaspers, N. C. Luhmann, H. K. Park, M. J. van de Pol, G. W.Spakman, and M. W. Jakubowski, “Effect of Heating on the Suppression ofTearing Modes in Tokamaks,” Physical Review Letters , 035001 (2007). E. Westerhof, A. Lazaros, E. Farshi, M. d. Baar, M. d. Bock, I. Classen,R. Jaspers, G. Hogeweij, H. Koslowski, A. Krämer-Flecken, Y. Liang, N. L.Cardozo, and O. Zimmermann, “Tearing mode stabilization by electroncyclotron resonance heating demonstrated in the TEXTOR tokamak andthe implication for ITER,” Nuclear Fusion , 85–90 (2007). H. Zohm, G. Gantenbein, F. Leuterer, M. Maraschek, E. Poli, and L. Urso,“Control of NTMs by ECCD on ASDEX Upgrade in view of ITER appli-cation,” Plasma Physics and Controlled Fusion , B341–B347 (2007). M. Maraschek, “Control of neoclassical tearing modes,” Nuclear Fusion ,074007 (2012). E. Kolemen, A. Welander, R. La Haye, N. Eidietis, D. Humphreys, J. Lohr,V. Noraky, B. Penaflor, R. Prater, and F. Turco, “State-of-the-art neoclas- sical tearing mode control in DIII-D using real-time steerable electron cy-clotron current drive launchers,” Nuclear Fusion , 073020 (2014). G. Ramponi, D. Farina, M. A. Henderson, E. Poli, G. Saibene, andH. Zohm, “ITER ECRH-ECCD System Capabilities for Extended PhysicsApplications,” Fusion Science and Technology , 193–201 (2007). L. Figini, D. Farina, M. Henderson, A. Mariani, E. Poli, and G. Saibene,“Assessment of the ITER electron cyclotron upper launcher capabilities inview of an optimized design,” Plasma Physics and Controlled Fusion ,054015 (2015). F. Poli, E. Fredrickson, M. Henderson, S.-H. Kim, N. Bertelli, E. Poli,D. Farina, and L. Figini, “Electron cyclotron power management for con-trol of neoclassical tearing modes in the ITER baseline scenario,” NuclearFusion , 016007 (2018). A. Reiman and N. Fisch, “Suppression of Tearing Modes by RadioFrequency Current Condensation,” Physical Review Letters (2018),10.1103/PhysRevLett.121.225001. E. Rodríguez, A. H. Reiman, and N. J. Fisch, “RF current condensation inmagnetic islands and associated hysteresis phenomena,” Physics of Plasmas , 092511 (2019). E. Rodríguez, A. H. Reiman, and N. J. Fisch, “RF current condensationin the presence of turbulent enhanced transport,” Physics of Plasmas ,042306 (2020), publisher: American Institute of Physics. S. Jin, N. J. Fisch, and A. H. Reiman, “Pulsed RF Schemes forTearing Mode Stabilization,” arXiv:2003.10896 [physics] (2020), arXiv:2003.10896. A. Smirnov and R. Harvey, “Calculations of the current drive in DIII-D withthe GENRAY ray tracing code,” Bull. Am. Phys. Soc. , p. 1837 (1995). C. Petty, R. Prater, J. Lohr, T. Luce, W. Fox, R. Harvey, J. Kinsey, L. Lao,and M. Makowski, “Detailed measurements of the electron cyclotron cur-rent drive efficiency on DIII-D,” Nuclear Fusion , 1366–1375 (2002). R. Prater, D. Farina, Y. Gribov, R. Harvey, A. Ram, Y.-R. Lin-Liu, E. Poli,A. Smirnov, F. Volpe, E. Westerhof, and A. Zvonkov, “Benchmarking ofcodes for electron cyclotron heating and electron cyclotron current driveunder ITER conditions,” Nuclear Fusion , 035006 (2008). D. G. Swanson,
Plasma Waves (Elsevier, 2012) google-Books-ID:0GMwFeDDZ9UC. I. Fidone, G. Granata, and R. L. Meyer, “Role of the relativistic mass vari-ation in electron cyclotron resonance wave absorption for oblique propaga- tion,” The Physics of Fluids , 2249–2263 (1982). E. Mazzucato, I. Fidone, and G. Granata, “Damping of electron cyclotronwaves in dense plasmas of a compact ignition tokamak,” The Physics ofFluids , 3745 (1998). A. Snicker, E. Poli, O. Maj, L. Guidi, A. Köhn, H. Weber, G. Conway,M. Henderson, and G. Saibene, “The effect of density fluctuations on elec-tron cyclotron beam broadening and implications for ITER,” Nuclear Fu-sion , 016002 (2018). L. Bardóczi, T. L. Rhodes, T. A. Carter, N. A. Crocker, W. A. Peebles,and B. A. Grierson, “Non-perturbative measurement of cross-field thermaldiffusivity reduction at the O-point of 2/1 neoclassical tearing mode islandsin the DIII-D tokamak,” Physics of Plasmas , 052507 (2016). L. Bardóczi, T. L. Rhodes, A. Bañón Navarro, C. Sung, T. A. Carter, R. J.La Haye, G. R. McKee, C. C. Petty, C. Chrystal, and F. Jenko, “Multi-field/-scale interactions of turbulence with neoclassical tearing mode mag-netic islands in the DIII-D tokamak,” Physics of Plasmas , 056106(2017). X. Garbet, P. Mantica, F. Ryter, G. Cordey, F. Imbeaux, C. Sozzi, A. Manini,E. Asp, V. Parail, R. Wolf, and t. J. E. Contributors, “Profile stiffness andglobal confinement,” Plasma Physics and Controlled Fusion , 1351–1373(2004). F. J. Casson, E. Poli, C. Angioni, R. Buchholz, and A. G. Peeters, “Ef-fect of turbulence on electron cyclotron current drive and heating in ITER,”Nuclear Fusion , 012002 (2015), arXiv: 1407.7767. F. Volpe, A. Hyatt, R. La Haye, M. Lanctot, J. Lohr, R. Prater, E. Strait, andA. Welander, “Avoiding Tokamak Disruptions by Applying Static MagneticFields That Align Locked Modes with Stabilizing Wave-Driven Currents,”Physical Review Letters (2015), 10.1103/PhysRevLett.115.175002. W. Choi, R. La Haye, M. Lanctot, K. Olofsson, E. Strait, R. Sweeney,F. Volpe, and The DIII-D Team, “Feedforward and feedback control oflocked mode phase and rotation in DIII-D with application to modulatedECCD experiments,” Nuclear Fusion , 036022 (2018). W. Choi and F. Volpe, “Simultaneous iterative learning control of modeentrainment and error field,” Nuclear Fusion , 056011 (2019).57