Effect of magnetic island geometry on ECRH / ECCD and consequences to the NTM stabilization dynamics

In the majority of codes that model ECCD-based NTM stabilization, the analysis of the EC propagation and absorption is performed in terms of the axisymmetric magnetic field, ignoring effects due to the island topology. In this paper, we analyze the wave propagation, absorption and current drive in the presence of NTMs, as well as the ECCD-driven island growth, focusing on the effect of the island geometry on the wave deposition. A primary evaluation of the consequences of these effects on the NTM evolution is also made in terms of the modified Rutherford equation.


Introduction
NTMs limit tokamak performance leading to loss of confinement or even disruption [1].Various techniques are under development for the control of these modes, with most important the application of localized ECCD.In the majority of models for the simulation of ECCD-based NTM stabilization, the analysis is performed on the basis of the unperturbed magnetic topology, taking into account only the flux surface of interest and ignoring any effects from the islands.However, the island topology introduces significant changes in the magnetic field and the plasma profiles as compared to the axisymmetric case, which may play a crucial role in the wave deposition and current drive.
In this work, we study the EC propagation, absorption and current drive in the presence of NTMs, using a ray-tracing code and a magnetic configuration that includes islands, and the resulting NTM dynamics are computed with the modified Rutherford equation.We focus on the effect of the island topology on the ECCD efficiency via the magnetic field, the flattening of the radial plasma profiles within the island region and the volumes of the flux surfaces onto which the wave is deposited.In Sec. 2 we introduce the non-axisymmetric magnetic configuration, in Sec. 3 we describe the ray-tracing scheme, which includes the computation of the volumes in the island topology, whereas in Sec. 4 the main results are presented and the implications of these results for the NTM dynamics are discussed.

Analytic description of the perturbed equilibrium
In order to derive an analytic model for the tokamak magnetic equilibrium including islands, we utilize the Clebsch form for the magnetic field, B = ∇ψ t × ∇θ + ∇ϕ × ∇ψ p , which is expressed in the toroidal coordinate system (r, θ, ϕ), with r the minor radius (R the major radius), θ the poloidal angle, ϕ the toroidal angle, and ψ p , ψ t the poloidal and toroidal fluxes.For a more explicit expression, we consider ψ t , ψ p as functions of r, θ, ϕ and transform to orthonormal toroidal coordinates a e-mail: ioanna@astro.auth.grRegarding the axisymmetric part of the magnetic field, a very convenient and widely-used representation is what is known as the "vacuum magnetic field" [2] where B t , B p are the toroidal and poloidal components, B 0 the toroidal field on the magnetic axis, ϵ A (r) = r/R 0 the inverse aspect ratio and q(r) = dϕ/dθ the safety factor.Regarding the perturbation, since magnetic islands appear at flux surfaces where q assumes rational values q(r s ) = m/n, with r s the minor radius of the O-point and m, n the toroidal and poloidal windings after which the field lines close, there should be two components of the magnetic field, a helical B h and a radial B r .For a perturbation amplitude ε mn , in a first approach, these fields at r = r s are For an implementation of the NTM topology we have to specify ε mn (r), and in order to preserve the divergence-free of the magnetic field we formulate the topology as a perturbation ψ p 1 to the background poloidal flux, i.e. ψ p = ψ p 0 + ψ p 1 .For the perturbation term we set the following form From Eqs. ( 1) and ( 4), the radial magnetic field is where ξ = mθ − nϕ and the exponential form of the radial dependence was chosen with the aim to make the perturbation local in radial direction.The magnetic lines b ρ (s) = ( ρ r (s), ρ θ (s), ρ ϕ (s) ) = (x (s) , y (s) , z (s)) are defined by and we set ds = √ dr 2 + r 2 dθ 2 + (R 0 + r cosθ) 2 dϕ 2 , so that |dρ/ds| = 1 and the parameter s is the arc-length along the field line.The three differential equations (6) are solved numerically in Cartesian coordinates by using a 4 th -order Runge-Kutta, adaptive step-size method.Figure 1 shows the poloidal cross section of the perturbed magnetic topology for the mode (m, n) = (3, 2).The cross-section is calculated by determining the surface of section in the poloidal plane at ϕ = 0.The plasma parameters are chosen as relevant to ITER, major axis R 0 = 6.2 m, minor axis α = 1.9 m, toroidal field on the magnetic axis B 0 = 5.51 T, and a monotonously increasing safety factor is assumed [2] with χ = r 2 /α 2 , so that q(0) = 1 at the magnetic axis and q(α) = 4 at the edge.The applied perturbation strength is ϵ 32 = 0.005, corresponding to an island width of approximately 16cm.
In order to calculate the volumes into which the wave energy is deposited, it is important to have an analytical labeling of the flux surfaces in the island region.For an adequate description of the island's morphology we derive a novel flux-surface label Ω that starts from a more detailed description of the helical and the radial field than just using their values at r = r s , and also contains less approximations.For this reason, when computing the helical and radial field perturbation terms introduced above, we do not assume that r = r s everywhere.The expression for Ω is then given by Ω can be normalized with a linear transformation in order to assume values in [−1, 1] within the separatrix.The width w of the magnetic island now being given by Apart from the changes in magnetic topology, an NTM causes the plasma pressure to assume a constant value inside the the island chain.This effect in the pressure profile leads in turn to a flattening in the electron density and temperature profiles.In order to model this situation, we assume the density and temperature profiles to be parabolic functions outside the island, as in the unperturbed case, constant in the island region and equal to the values at the outer island boundary.

Ray tracing, wave absorption and current drive
The Hamiltonian equations for the ray propagation in the GO approach are [3] where r is the position of the ray, k is the wave vector, τ is a parameter denoting the length along the ray path, and the Hamiltonian function H is given by the dispersion relation of the wave, where ϵ h is the Hermitian part of the dielectric tensor and I the unit tensor.Assuming cold plasma propagation, one can adopt cold plasma response and get the expression for the Hamiltonian in [4].The calculation of the wave absorption along the ray path can be done in terms of the imaginary part of the wave vector, as determined from the mode-specific dispersion relation [5] dP where v g is the group velocity.Knowing the power of the ray along its path, the power dP w deposited in a small radial interval can be calculated from Eq. ( 12), and division by the corresponding volume dV s , which is contained between the two flux surfaces enclosing the radial interval, yields the absorbed wave power density.The total driven current per wave power defines the current drive efficiency η CD EC-17 Workshop Following that, the computation of the ECCD can be done in terms of the linear adjoint method [6].
For the accurate estimation of the wave power density, the calculation of the plasma volume between two nearby flux surfaces is necessary, especially in case the local topology is expected to affect power deposition.We first need to calculate the total volume V s contained inside one flux surface where ξ = mθ − nϕ is the angle coordinate perpendicular to the helical line through the O-point and r 1 , r 2 , ξ 1 and ξ 2 are the integration limits for r and ξ, which are determined as follows: A given fluxsurface has a unique value of Ω = Ω i .The equation Ω(r s , ξ) = Ω i determines the limits ξ 1 , ξ 2 in the helical (ξ) direction and at the radius r = r s , and, for ξ in its integration interval [ξ 1 , ξ 2 ], the equation Ω(r, ξ) = Ω i determines the limits r 1 (ξ), r 2 (ξ) in the radial direction as a function of ξ.This yields ] .
The absorbed power per unit volume can then be evaluated as dP w /dV s = (dP w /dτ)/(dV s /dτ).

Numerical Results
In this section, we present the numerical results on the effect of the magnetic island geometry on the EC wave propagation, power deposition and current drive.The wave is injected from the outermost flux surface (r = α), with the toroidal and poloidal launching angles such that the ray is targeted to a region very close to the O-point, and the EC resonance is located around the center of the island.We present results for the plasma and wave parameters of ITER: R 0 = 6.2 m, α = 1.9 m, B 0 = 5.51 T, n e 0 = 10 20 m −3 , n e α = 10 19 m −3 , T e (0) = 10 KeV, T e (α) = 1 KeV, ω/2π = 160 GHz (the fundamental O-mode), P w0 = 10 MW.Regarding the parameters of the NTM generating the magnetic islands, we throughout consider the mode (3, 2), with a magnetic perturbation strength ε 32 = 5 • 10 −3 , which corresponds to an island width of approximately 16cm.In Fig. 1(b), we show ray paths in the presence of the NTM for two different poloidal launching angles.The figure also shows the corresponding ray paths in the absence of the NTM, i.e. in the unperturbed equilibrium, for a comparison of the two cases.The waves are injected at the poloidal angle θ l = 30 o , with the toroidal launching angle being ϕ l = 0 o .For the assumed profiles of B, n e , T e , and the value of ω, the region of the EC resonance is located around R = 6 m and marked by two vertical lines.The ray propagation is seen not to be much affected by the magnetic field of the island, since the relative magnetic perturbation strength (∼ 10 −3 ) is small.
As explained above, with a magnetic island present, the wave power is absorbed in smaller fluxsurface volumes, with the direct consequence that the ECCD is strongly increased.This is shown in Fig. 2 for a typical case, where the wave is injected at the poloidal angle θ l = 30 o and the toroidal launching angle is 0 o .In the presence of an NTM, the absorbed power density is found to be nearly 2 times larger than the one in the unperturbed case.The driven current density has a similar behavior, attaining a value almost 2.5 times larger than in the absence of the NTM perturbation.
Based on the previous results, the question arises whether the effects of the island geometry on the ECCD efficiency also play a significant role in NTM stabilization.Here, only a preliminary investigation of this issue is made.The NTM dynamics can sufficiently be described by the modified Rutherford equation, which is the established model for the dynamic evolution of NTMs where W is the half-width of the magnetic island, τ r = 936.4r 2 s is the resistive time-scale of the plasma, ∆ ′ = −m is the classical TM stability index, ∆ ′ BS is the stability index corresponding to the bootstrap  current, which becomes significant for high values of the plasma pressure, and ∆ ′ CD is the term which represents the stabilizing effect of the ECCD on the island.The ∆ ′ -indices are defined as ] .
In the above, ϵ A is the aspect ratio of the tokamak, β P = 2µ 0 P/B 2 is the ratio of the plasma pressure over the magnetic pressure, L q and L p are the shear lengths of the safety factor and the plasma pressure, d CD is the wave power deposition width ( j CD0 and g wd are defined below) and W d ,W b , W pol are characteristic threshold values for the island width [7].
Apart from the total driven current, two more characteristics of the ECCD are important for NTM stabilization, the spatial distribution of the driven current density and its efficiency in replacing the missing bootstrap current.The current density is assumed here to exhibit a Gaussian distribution in space, j CD = j CD0 exp , where r mis is the small distance of misalignment of the ECCD deposition center from the island's O-point.The efficiency of the ECCD in stabilizing the NTM is described by the parameter g wd = 0. [7].The phase diagram of the modified Rutherford equation is shown in Fig. 3 for toroidal launching angle 0 o and the three different cases of (i) no ECCD applied, (ii) ECCD applied but calculated in the absence of the island, and (iii) ECCD applied and calculated by taking into account the island topology.

01013-p.5
For the plasma and wave parameters we use the same ITER-relevant values as in the previous section, and for the parameters newly introduced in this section we set β P = 0.5, L q /L p = 1, w d = 0.01r s , w b = 0.02r s , w pol = 0.015r s , r mis = 1 cm.Fig. 3(a) shows that (i) without ECCD applied the island is unstable to grow in a wide range of widths, (ii) with ECCD applied but calculated without taking the island topology into account, there is still a region of unstable widths, much smaller now in extent, and (iii) with ECCD applied and calculated by taking the island geometry into account, there is an effective NTM stabilization, since dW/dt has becomes negative for all widths.The reason for the effectiveness is the peaking of the ECCD profile when the island topology is taken into account.

Summary and Discussion
In this paper we studied EC wave propagation, absorption and current drive in the presence of magnetic islands, using ray tracing and a novel approach to the modeling of a tokamak magnetic field that includes islands, and investigated the effects of the island geometry on the NTM dynamic evolution.Up to now, this issue has been analyzed by ignoring any effects from the island topology on the ECCD efficiency.The analytic model for the magnetic islands superimposed onto the tokamak equilibrium field is more realistic than the standard models adopted so far in the literature, using a novel fluxsurface labeling in the island region which improves the precision in the determination of the fluxsurface volumes and leads to increased accuracy in the calculation of the power deposition.
Our main results imply that: (a) In most cases, the propagation until the O-point is little affected by the magnetic island topology, due to the very small magnitude of the perturbation, with any effects owed to the flattening of the plasma pressure within the island region.(b) The absorbed wave power and the driven current are affected significantly by the island geometry, namely enhanced, the energy is absorbed and redistributed by the plasma particles in flux-surface volumes much smaller than the ones in an axi-symmetric topology.(c) The use of an island topology allows a more accurate estimate of the NTM stabilization effect of ECCD, and it turns out that the increased ECCD peaking in the islands topology results in a much better efficiency of NTM stabilization, i.e. stabilization can be achieved with less wave power than the one estimated in axisymmetric geometry.
Concerning the basic assumptions made in our modeling approach, the routines for the computation of the wave propagation, absorption and current drive are based on linear methods that do not account for (a) wave diffraction, (b) trapping of particles in the island chain, or (c) nonlinear effects due to the very large values that the absorbed power density may attain.Also, the effect of the island geometry on the NTM dynamics has been studied in terms of a quasi-self-consistent coupling of the wave propagation solver to the equation governing the island dynamics.The development of a fully self-consistent model, either in terms of changing the magnetic topology at certain time slices during the wave propagation, or with ray tracing in time-dependent fields, is under way.Other studies have been made [8], simulating the ECCD process in more detail than with the linear models used here, yet without accounting for the consequences of the island topology on ECCD, which find that the transport properties enhance the current drive efficiency compared to the axisymmetric case.

Fig. 1 .
Fig. 1.Poloidal cross section in Cartesian coordinates (a) for the entire poloidal plane and (b) zoomed into the magnetic island area, for the (3, 2) mode and with ITER-like parameters.

Fig. 2 .
Fig. 2. Deposited wave power per unit volume dP w /dV s (a) and the corresponding driven current density (b) as a function of R, for the cases of propagation with a (3,2) NTM mode present and without perturbation.

ig. 3 .
Phase diagram of the modified Rutherford equation for m/n = 3/2 and ϕ l = 0 o .