Non-linear modeling of the plasma response to RMPs in ASDEX Upgrade

The edge localized modes (ELMs), occurring at the plasma edge in tokamaks, induce large transient heat loads on plasma facing components. In ITER, the heat loads on the divertor due to ELMs are expected to be intolerable for materials if unmitigated, hence motivating a broad effort on developing and understanding reliable ways to mitigate the ELMs. One of the promising methods to control ELMs is the application of non-axisymmetric resonant or non-resonant magnetic perturbations (RMPs or NRMPs) by dedicated coils. The original idea is that a small applied magnetic perturbation ($\delta B/{{B}_{\text{toroidal}}}\sim {{10}^{-4}}$) is expected to induce magnetic islands on the resonant surfaces characterized by the safety factor q  =  m/n (where m and n are respectively the poloidal and toroidal mode numbers, with the dominant n fixed by the RMP configuration). At the plasma edge where rational surfaces are close to each other, consecutive island chains are likely to overlap and induce an ergodic layer [1]. The radial particle and heat transport being enhanced in the ergodic zone, RMPs should be able to slightly deconfine the plasma edge and reduce the pedestal pressure gradient under the ELM-triggering threshold (so far understood in the framework of the peeling–ballooning or P–B theory) [2, 3]. In addition, RMPs are also believed to lower the P–B stability boundary limit, resulting in more frequent, smaller ELMs [4, 5].

When applying RMPs, the mitigation or suppression of type-I ELMs was successfully obtained in DIII-D, ASDEX Upgrade (AUG), JET, MAST, KSTAR, NSTX and EAST tokamaks [6–12]. However the conditions to obtain a strong mitigation or complete ELM suppression have been proven to be more complicated than this simple picture, due to the strong screening of RMPs by plasma flows. Even though the understanding of the plasma response to RMPs has been much improved in recent years [13–24], the interaction between plasma and RMPs has to be further studied to better interpret the current experiments and to be capable to predict the effect of RMPs in ITER. In particular, in addition to the resonant response, the role of the so-called 'edge kink response' or 'peeling response' has to be assessed, where the amplification of the external field perturbation is observed, probably resulting from the interaction between the applied perturbation and a marginally stable kink mode [25]. At low collisionality, the ELM suppression in DIII-D and the strongest ELM mitigation in AUG were recently shown to be related to the excitation of the edge kink (or peeling) response [26–28].

This paper aims at better understanding the role of both the resonant and the kink responses on the ELM mitigation by RMPs in AUG at low collisionality. In this respect, the plasma response to RMPs was modeled with the non-linear resistive MHD code JOREK [29, 30], using the data extracted from different RMP configurations tested in ASDEX Upgrade experiments. In section 2, the experimental discharges considered and the comparison between the experimental data and the kinetic and rotation profiles in simulation are described. In section 3, the generic features of the plasma response to RMPs are discussed. In section 4, the impact of the differential phase between upper and lower RMP coils (varied in experiments and simulations accordingly) on the plasma response is assessed. The evolution of the pedestal profiles is then described in section 5, and discussions and conclusions are provided in section 6.

2.1. Experimental discharges

The data used in modeling are extracted from AUG shots #30826 and #31128 (figure 1), where the effect of magnetic perturbations was studied at low collisionality. In these experiments, the currents in the two rows of 8 RMP-coils [31, 32] were set such that magnetic perturbations with a dominant toroidal mode number n  =  2 were applied. Small n  =  6, 10, 14... sideband modes were also induced due to the geometry of the coils, however vacuum calculations showed that their amplitude is very small as compared to the n  =  2 amplitude [33], thus they were neglected in the modeling described in this paper.

Zoom In Zoom Out Reset image size

In the discharge #30826 presented in figure 1(a), the differential phase between the magnetic field perturbation generated by upper and lower RMP-coils is slowly rotated between $\Delta \Phi =+{{90}^{\circ}}$ and $\Delta \Phi =-{{90}^{\circ}}$ by varying the current applied in the coils. In this shot, the strongest ELM mitigation (characterized by the largest ELM frequency, see figure 1(a)) is obtained for $\Delta \Phi$ between $+{{90}^{\circ}}$ and $+{{60}^{\circ}}$ and it is correlated with the strongest density pumpout observed. The increase in ELM frequency is probably induced by the reduction of the pedestal stability limit due to the application of RMPs [4, 5]. Note that a transient ELM-free phase (in blue) is observed for $\Delta \Phi =-{{30}^{\circ}}$ to $-{{90}^{\circ}}$: this phenomenon is not fully understood, but a possible explanation is that the stability limit may be enhanced in these RMP configurations, allowing the pedestal to grow until it reaches the new stability limit [32]. This ELM-free state was only obtained in transient phases so far. After this shot, several discharges with similar plasma configuration and parameters have been performed, applying steady magnetic perturbations for different $\Delta \Phi$ values. The strongest ELM mitigation was obtained for $\Delta \Phi =+{{90}^{\circ}}$, in the discharge #31128 given in figure 1(b). Since #31128 was better diagnosed, the equilibrium profiles before RMP application of this discharge were used as input for our simulations. Starting from this equilibrium plasma, steady magnetic perturbations were applied for $\Delta \Phi =+{{90}^{\circ}}$ (phase in shot #31128), $+{{60}^{\circ}}$, $+{{30}^{\circ}}$, ${{0}^{\circ}}$, $-{{30}^{\circ}}$, $-{{60}^{\circ}}$ and $-{{90}^{\circ}}$, to scan the effect of the magnetic perturbation spectrum, comparably to the discharge #30826 presented in figure 1(a).

Initial average temperature $T=\left({{T}_{\text{e}}}+{{T}_{\text{i}}}\right)/2$ and electron density n_e profiles (plotted in figure 2(a) with error bars of the measurements as a function of the normalized poloidal flux ${{\psi}_{n}}$) are fitted from experimental data. In our model, it is assumed that ion and electron temperature are equal: $T={{T}_{\text{e}}}={{T}_{\text{i}}}$. Note that the density profile was chosen to be close to the fitted experimental profile, while the temperature gradient (and thus the pressure gradient) was enhanced in the pedestal within the error bars in order to make the plasma peeling–ballooning unstable. This steepening thus induced a smaller temperature pedestal width as compared to the density pedestal width in our simulations. The equilibrium reconstruction and the $F{{F}^{\prime}}$ profile calculation are made with the equilibrium code CLISTE [34]. The resulting q-profile is plotted in figure 2(b), with the position of the resonant surfaces q  =  m/2 (for an n  =  2 perturbation) marked as black diamonds.

Zoom In Zoom Out Reset image size

The JOREK extended reduced MHD model described in [20] is used for modeling. The equations are solved for the following variables: magnetic flux ψ, mass density $\rho ={{m}_{\text{i}}}{{n}_{\text{e}}}$ (with m_i the ion mass and n_e the electron density), temperature $T={{T}_{\text{e}}}={{T}_{\text{i}}}$, parallel velocity ${{V}_{\parallel}}$, electric potential u, toroidal current $j={{ \Delta }^{\ast}}\psi$ and vorticity $W=\nabla _{\bot}^{2}u$. Perpendicular transport is modeled with diffusive heat and particle terms (with heat and particle diffusivity coefficients reduced in the pedestal to reproduce the Edge Transport Barrier) balanced with heat and particle sources, such that the temperature and density profiles do not evolve in time in axisymmetric simulations (n  =  0 only). Since it is challenging to properly reproduce realistic sources in modeling, for simplicity the particle source was chosen to be constant over the plasma and the heat source profile was taken to follow the initial temperature profile. These sources, as well as the heat and particle diffusive coefficients, are kept constant along the simulation. For numerical reasons, the central resistivity is taken  ∼10 times larger than the Spitzer value and the resistivity profile follows a T^−3/2 dependence.

A simulation is run as follows: in a first step, the axisymmetric equilibrium (n  =  0 only) is calculated, and the (n  =  2) perturbations are added in a second step. Note that larger n modes (leading to ELM eruptions) were not included in this modeling, so as to consider only the plasma response to RMPs. The interaction between ELMs and RMPs will be the subject of future works. Firstly, in the simulations, it is important to consistently model the plasma flows in order to properly reproduce the radial force balance. As described in [20], a source of parallel rotation mimicking the experimental rotation profile, the two-fluid diamagnetic flows as well as the neoclassical friction constraining the poloidal velocity are included in the model to self-consistently describe these plasma flows.

For the source of parallel rotation (considered to be close to the toroidal rotation), the (main) ion toroidal rotation frequency ${{ \Omega }_{\text{tor}}}$ was extracted from experimental measurements at the outboard midplane, assuming that main ions and impurity ions rotate similarly in toroidal direction. The toroidal rotation frequency was taken as a flux function (${{ \Omega }_{\text{tor}}}=f(\psi )$), consistent with the experimental observations to lowest order [35]. The source of toroidal rotation was then implemented as a constraint in the viscous term of the parallel momentum equation: $\mu \Delta \left({{V}_{\parallel}}-{{ \Omega }_{\text{tor}}}\centerdot 2\pi R\right)$, where μ denotes the viscosity coefficient. This way, the equilibrium parallel rotation profile (plotted in figure 3(a)) perfectly matches the experimental profile at the outboard midplane. The parallel rotation profile at the inboard midplane, corresponding to ${{V}_{\parallel}}={{ \Omega }_{\text{tor}}}(\psi )\centerdot 2\pi R$, is also given. Note that a more complex implementation taking into account the neoclassical toroidal viscous torque [36–38] would be necessary to model the toroidal rotation more accurately. Such implementation will be the subject of future works.

Zoom In Zoom Out Reset image size

As for the main ion poloidal rotation, it cannot be directly extracted from experiments, but previous studies showed that the poloidal rotation is very close to the neoclassical prediction in the pedestal [39], with a good match between experimental data and calculations with the NEOART code [39, 40]. The neoclassical poloidal velocity calculated with JOREK (figure 3(b)) shows good agreement with the NEOART calculation using the same steepened temperature and density profiles as input. The JOREK value (black dashed line) is slightly larger than the NEOART value in the pedestal due to the effective charge Z_eff taken as 1 in JOREK. Note that both JOREK and NEOART poloidal neoclassical velocities would have been smaller if the temperature and density profiles had not been steepened, thus the neoclassical poloidal rotation is possibly overestimated in JOREK. The neoclassical friction is implemented in JOREK as a neoclassical tensor $\nabla \centerdot {{ \Pi }_{\text{i},\text{neo}}}\propto {{\mu}_{\text{neo}}}\left({{V}_{\theta}}-{{V}_{\theta,\text{neo}}}\right)$ [20]. The resulting neoclassical velocity in JOREK ${{V}_{\theta}}$ (plotted in blue) is close to ${{V}_{\theta,\text{neo}}}$ but slightly differs due to other terms (among others the stress tensor) constraining the velocity in the poloidal momentum equation.

The resulting equilibrium radial electric field E_r, calculated from the equilibrium force balance, is given in figure 3(c). The E_r radial profile in JOREK matches rather well the E_r estimated from experimental profiles using the NEOART calculation of the neoclassical poloidal velocity. Note that the E_r value estimated from experimental profiles would be smaller by a factor of two in the pedestal if the pedestal temperature profile was not steepened to make the plasma P–B unstable (figure 3(c), black dashed line). In both cases, this represents a large radial electric field well in the pedestal (and thus a large electron perpendicular rotation), capable to strongly screen RMPs at the pedestal top, as discussed below.

Once equilibrium (n  =  0) flows are established, n  =  2 perturbations are added in the simulation. The n  =  2 magnetic flux perturbation induced by RMP coils is beforehand calculated in the vacuum with the VACFIELD code [41], and applied in JOREK as boundary condition for the n  =  2 magnetic flux perturbation. This constraint does not allow for the modification of the magnetic perturbation at the boundary and thus may slightly damp the field amplification inside the plasma. In future works, RMP coils will be directly implemented in the JOREK-STARWALL version of the code [42], allowing for field amplification at the boundary. The vacuum field and the JOREK boundary are presented in figure 4(a). The magnetic perturbation at the boundary is progressively increased to its nominal value in a thousand Alfvén times, thus RMPs progressively penetrate into the plasma taking into account the plasma response. The magnetic energy of the n  =  2 mode, plotted in figure 4(b), is growing due to the external RMP application, until it saturates. The 2D-plot of the n  =  2 magnetic flux perturbation penetrating in the plasma (in the $\Delta \Phi =-{{90}^{\circ}}$ case) as well as the n  =  2 response toroidal current perturbation on resonant surfaces and in the SOL, are given in figures 4(c) and (d).

Zoom In Zoom Out Reset image size

In all $\Delta \Phi$ configurations, the plasma response shows similar features: magnetic island chains are formed on the resonant surfaces q  =  4/2, 5/2 and 6/2, as shown in the Poincaré plot in (${{\psi}_{\text{norm}}}$, θ) coordinates given in figure 5(b). This is due to the fact that the perpendicular electron velocity (figure 5(c)) is sufficiently small in this region of the plasma to prevent the formation of large screening currents on these resonant surfaces [14, 17, 20]. However in the pedestal, the electron perpendicular flow is very large (large radial electric field well), thus screening currents are large on the resonant surfaces q  =  7/2 and 8/2 and only very small magnetic islands can be seeded on these surfaces. At the very edge (for $q\geqslant 9/2$ and ${{\psi}_{\text{norm}}}\geqslant 0.97$), the resistivity $\eta \propto {{T}^{-3/2}}$ is large enough to prevent the screening current formation; subsequently, an ergodic layer is induced by RMP penetration at the very edge. This ergodic layer comes with the formation of lobe structures near the X-point, observable in the Poincaré plot in (R, Z) coordinates (figure 5(a)) and also observed in experiments [43].

Zoom In Zoom Out Reset image size

Even though these general observations are valid for all $\Delta \Phi$ configurations, some important characteristics such as the size of the ergodic layer and the lobe structures as well as the kinking of the field lines, differ sensitively depending on $\Delta \Phi$, i.e. on the applied spectrum, as discussed in section 4.

Particularly at the plasma edge, the plasma response substantially depends on the differential phase between the coils applying the RMPs. First, as shown in figure 6 for three representative simulations, a different kinking of the field line (due to the ideal plasma response—so called kink response) is observed in the different configurations. The kinking at the edge around the midplane ($\theta =0$) is maximal for $\Delta \Phi =-{{90}^{\circ}}$ and minimal for $\Delta \Phi ={{0}^{\circ}}$, while the kinking near the X-point ($\theta \approx -\pi /2$) is maximal for $\Delta \Phi =+{{90}^{\circ}}$ and minimal for $\Delta \Phi =-{{90}^{\circ}}$.

Zoom In Zoom Out Reset image size

Consistently, a similar trend can be observed in the 3D-profiles of density and temperature. At the edge of the outboard midplane, the distortion of the density profile (figure 7) when following the toroidal direction is maximal for $\Delta \Phi =-{{90}^{\circ}}$ (the radial variation between the maximal and minimal toroidal excursion reaches $\Delta R\approx 8$ mm on q  =  9/2) while this displacement is smaller for $\Delta \Phi ={{0}^{\circ}}$ ($\Delta R\approx 4$ mm) and $\Delta \Phi =+{{90}^{\circ}}$ ($\Delta R\approx 3.5$ mm).

Zoom In Zoom Out Reset image size

In the same way, the distortion of the density profile near the X-point shows evidence of the field line deformation at the vicinity of the X-point. Figure 8, which follows a vertical line passing through the X-point along the toroidal direction, shows that the vertical distortion of the density is largest for $\Delta \Phi =+{{90}^{\circ}}$ ($\Delta Z\approx 7$ mm on q  =  9/2, while $\Delta Z\approx 4$ mm for $\Delta \Phi ={{0}^{\circ}}$ and $\Delta Z\approx 2$ mm for $\Delta \Phi =-{{90}^{\circ}}$).

Zoom In Zoom Out Reset image size

It is thus important to notice that the case for which the strongest ELM mitigation is obtained in experiments (largest ELM frequency for $\Delta \Phi =+{{90}^{\circ}}$) corresponds to the largest (so-called peeling-kink) displacement near the X-point found in this modeling. A similar correlation between the observation of a strong plasma response and a large X-point displacement (calculated with MARS-F) was also found in MAST L-mode plasmas [44].

In addition, this best-ELM-mitigation case in experiments also corresponds to the largest ergodic layer obtained in modeling. Indeed, in figure 6, the ergodic layer induced by RMPs at the edge is slightly thicker in the $\Delta \Phi =+{{90}^{\circ}}$ case (the magnetic field is ergodic for ${{\psi}_{\text{norm}}}\geqslant 0.97$, compared to ${{\psi}_{\text{norm}}}\geqslant 0.98$–0.99 in the other cases), and this is correlated with the presence of slightly longer lobe structures near the X-point. The longer lobes also lead to longer footprint (n  =  2) patterns on the divertor, as plotted in figure 9: for $\Delta \Phi =+{{90}^{\circ}}$, the footprints induced by RMPs on the outer divertor target are slightly larger than for $\Delta \Phi ={{0}^{\circ}}$, and significantly larger than for $\Delta \Phi =-{{90}^{\circ}}$. Note that in the $\Delta \Phi =+{{90}^{\circ}}$ case, the (n  =  2) magnetic perturbation has non-linearly induced a global (n  =  0) displacement of the X-point, and hence a 2 mm-displacement of the radial position of the footprint on the divertor, observable in figure 9. Unfortunately no infra-red thermography data were obtained in this experimental discharge to be compared with the footprints found in modeling.

Zoom In Zoom Out Reset image size

Furthermore, the size of the ergodic layer (and the size of the lobe structures and footprints as well) is actually correlated with the amplitude of the resonant component of the magnetic perturbation. This can be explained by the fact that the resonant components determine the width of the magnetic islands and therefore the width of the overlap region, while being at the same time related to a quantity which determines the lobe length (Poincaré [45] or Melnikov [46] integrals). Indeed, the amplitudes of the (n  =  2) Fourier harmonic of the magnetic flux $|{{\psi}_{mn,\text{res}}}|$ on the resonant surfaces q  =  m/n (plotted in figure 10 as a function of the differential angle between upper and lower RMP coils for m  =  7 to 11), are maximum for the $\Delta \Phi =+90/+{{60}^{\circ}}$ cases. The resonant flux perturbation on q  =  7/2 and 8/2 (black crosses and blue circles in figure 10) is small in all configurations due to the large screening induced by the large perpendicular electron flow at the pedestal top, as previously discussed in sections 2 and 3. However, on the resonant surfaces located close to the separatrix (q  =  m/2 for $m\geqslant 9$), the resonant component is largest for $\Delta \Phi =+90/+{{60}^{\circ}}$ and is reduced as $\Delta \Phi$ is reduced, reaching a 10 times lower value for $\Delta \Phi =-60/-{{90}^{\circ}}$: this means that the width of the stochastic layer at the edge is progressively reduced as $\Delta \Phi$ is reduced during this discharge (figure 1(a)).

Zoom In Zoom Out Reset image size

The 2D-plot of the (n  =  2) component of the radial magnetic field perturbation $\delta {{B}_{r}}=\vec{B}\centerdot \nabla {{\psi}_{0}}/|\nabla {{\psi}_{0}}|\,=$$\left({{\partial}_{R}}{{\psi}_{0}}{{\partial}_{Z}}\psi -{{\partial}_{Z}}{{\psi}_{0}}{{\partial}_{R}}\psi \right)/\left(R|\nabla {{\psi}_{0}}|\right)$ as a function of the poloidal mode number m and the radial direction ${{\psi}_{\text{norm}}}$ (figure 11) highlights particularly well the amplitude of the resonant components (marked by white diamonds) and the kink components, for the different cases $\Delta \Phi =+{{90}^{\circ}}$, ${{0}^{\circ}}$ and $-{{90}^{\circ}}$. In all cases, the resonant component is close to zero at the center due to the strong screening. At the edge (for ${{\psi}_{\text{norm}}}>0.95$), the resonant component (on q  =  m/2 with $m\geqslant 8$) has a large amplitude for $\Delta \Phi =+{{90}^{\circ}}$, a smaller amplitude for $\Delta \Phi ={{0}^{\circ}}$ and a very small amplitude for $\Delta \Phi =-{{90}^{\circ}}$. As for the kink component, composed of modes m  >  n q [47], its amplitude is relatively large in the core in the $\Delta \Phi =-{{90}^{\circ}}$ and $+{{90}^{\circ}}$ cases, but the kink amplitude is especially large at the edge in the $\Delta \Phi =+{{90}^{\circ}}$ case. This edge amplification corresponds to the peeling-kink displacement near the X-point previously discussed in this section.

Zoom In Zoom Out Reset image size

A zoom into figure 11 for the $\Delta \Phi =+{{90}^{\circ}}$ case (given figure 12) allows to exhibit the strong coupling between the resonant m component (circled in blue) and the m  +  2 peeling-kink component (circled in green). This suggests that the peeling-kink amplification near the X-point induces the amplification (or prevent the complete screening) of the resonant magnetic perturbation at the edge. The amplified resonant component therefore induces a larger ergodicity at the edge and an increased radial transport. However, considering the 2D-spectrum of the magnetic perturbation ${{\psi}_{mn}}(n=2)$ (figure 13) rather than the spectrum of the radial magnetic field $\delta {{B}_{r}}$ evidences different mode couplings. Indeed, even though the screening is only partial on the resonant surfaces (with the largest resonant component in the $+{{90}^{\circ}}/+{{60}^{\circ}}$ cases, as shown in figure 10), a strong screening of the resonant component ${{\psi}_{mn}}(n=2)$ is observed in all cases (white diamonds in figure 13). As for the peeling-kink component m  >  nq, it is also mostly excited in the $+{{90}^{\circ}}$ case, consistently with the $\delta {{B}_{r}}$ picture, excepted that the m  +  1 component is excited as well as the m  +  2 component.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

To understand what is theoretically the important quantity to consider between $\delta {{B}_{r}}$ and ψ in order to quantify the shielding of the perturbation, we need to calculate the radial displacement of the field lines generated by a magnetic island structure. The following calculation quantifies this radial displacement, assuming a near-circular geometry: the ratio of the radial to the poloidal extent of the island is given by: $\text{d}r/\left(r\text{d}\theta \right)={{B}_{r}}/{{B}_{{{\theta}_{0}}}}$, and the radial derivative of the poloidal flux is: $\text{d}{{\psi}_{0}}/\text{d}r=R\centerdot {{B}_{{{\theta}_{0}}}}$. Thus the displacement in the radial direction (along ${{\psi}_{0}}$) is: $\text{d}{{\psi}_{0}}/r\text{d}\theta =R\centerdot {{B}_{r}}$. So the relevant quantity to consider is rather ψ or the product $R\centerdot {{B}_{r}}$. Subsequently the examination of the B_r components introduce artefact couplings due to the $\cos (\theta )$-dependance (at the first order) of $R\sim {{R}_{0}}\left(1+\epsilon \cos (\theta )\right)$.

Finally when the peeling-kink component is most excited (in the case when ELMs were most strongly mitigated in experiments—for $\Delta \Phi =+60/+{{90}^{\circ}}$), it is also correlated with the largest resonant component at the edge resonant surfaces (due to the poloidal coupling between the m  >  nq kink components and the resonant m component) and thus is correlated to the largest stochastic layer observed in modeling. However the resonant components are globaly shielded (even though the screening is only partial), and the strong coupling between the m and the m  +  2 component observed in $\delta {{B}_{r}}$ spectrum is actually an artefact of the calculation of $\delta {{B}_{r}}$.

The 3D density profile is not only corrugated following the n  =  2 displacement, but the axisymmetric n  =  0 profile is also non-linearly affected by n  =  2 RMPs. Thus a significant reduction of the overall density (so-called density pumpout) is observed in modeling, as shown in figure 14. The outboard midplane density profile is plotted without RMPs and with RMPs in the three configurations $\Delta \Phi =+{{90}^{\circ}}$, ${{0}^{\circ}}$ and $-{{90}^{\circ}}$, 60 ms after switching on the RMPs. The comparison of the $\Delta \Phi =+{{90}^{\circ}}$ case with the experimental profiles (shot #31128, figure 14 (right)) shows that the density pumpout obtained in modeling is actually smaller than the one observed in experiments, but is still significant enough to be underlined. The time scale of pumpout to be fully developed in experiments is around  ∼200 ms. In the modeling, a 'long-time' run (∼150 ms) made for a single case (in order to limit the computational time) shows that the density reduces exponentially when RMPs are applied, until a 3D-steady state equilibrium is reached after 100–150 ms. Already after 60 ms, the density reduction in the pedestal has reached 90–95% of the total density reduction obtained at steady-state. Thus, even though the steady-state density profile would be reduced by 5–10% as compared to the profiles shown in figure 14 (left), other additional mechanisms are necessary to explain the whole pumpout observed in experiments (14 (right)). Physical mechanisms such as the density transport by magnetic flutter [48] and by turbulence [49], not included in this model, are good candidates to explain this discrepancy and remain to be studied in future works. Corrections in the radial electric field, as suggested in [50], might also be considered.

Zoom In Zoom Out Reset image size

In addition, the density pumpout observed in modeling is slightly larger in the $\Delta \Phi =+{{90}^{\circ}}$ case than for other RMP configurations, but the difference is not significant enough as compared to the difference observed in the experiments. So the additional effects quoted above may be necessary to consistently model the actual density pumpout.

As for the temperature evolution due to RMPs in modeling (see the modeled temperature profile at the bottom-left of figure 14), the pedestal temperature is reduced by the RMP application similarly to the pedestal density ($12 \%$ of reduction for both pedestal values), while the central temperature remains almost unaffected by the RMP application (1–2% of reduction at the center, as compared to the  ∼7% of density reduction at the center). As compared to the experimental observations (bottom-right of figure 14), the reduction of the temperature in the pedestal was not observed in experiment, but the unaffected temperature at the center found in modeling is consistent with the experimental measurements.

The modeling of the plasma response to RMPs in ASDEX Upgrade, depending on the applied RMP spectrum (changed by varying the differential phase ▵Φ between upper and lower coil currents), was performed with the resistive non-linear MHD code JOREK. Input profiles and RMP spectrum were chosen to closely match the experimental data. Two-fluid diamagnetic rotation, neoclassical friction and a source of toroidal rotation implemented in the model allowed to reproduce equilibrium rotation profiles similar to the experimental ones. The enhancement of the pedestal density and temperature profiles (in order to make the plasma peeling–ballooning unstable) may nevertheless lead to an overestimated poloidal velocity and thus a slightly overestimated electron perpendicular rotation in the pedestal. Subsequently the shielding of RMPs in the pedestal may be accentuated in the modeling by the large electron perpendicular flow. On the other hand, the resistivity taken in modeling (10 times larger than the Spitzer value) should on the contrary increase the size of magnetic islands and the width of the ergodic layer. In spite of these limitations, the qualitative behaviour of the plasma response to RMPs is modeled consistently with experimental observations.

For all RMP configurations, a similar behaviour was observed in modeling: the screening of magnetic perturbations prevent the formation of significant islands at the pedestal top due to the large electron perpendicular flow. At the very edge (last 1 to $3 \%$ of the normalized poloidal flux inside the separatrix), an ergodic layer is induced by the RMP penetration, due to the local high resistivity (proportional to T^−3/2). RMPs also induce the 3D deformation of the separatrix and the formation of lobe structures near the X-point, which impose an n  =  2 footprint pattern on the divertor target plates.

However, significant differences appear depending on the applied RMP spectrum (varied via $\Delta \Phi$). In the $\Delta \Phi =+{{90}^{\circ}}$ case, the largest kinking of the field lines is observed near the X-point (as compared to other configurations), also inducing the largest 3D displacement of density and temperature profiles around the X-point. This large peeling-kink amplification near the X-point also generates the amplification of the resonant (tearing) response, via the coupling between the (m  >  nq) peeling-kink component with the m resonant component. Due to the amplified resonant component, a larger ergodic layer is observed in this configuration (also going with longer lobe structures) and a larger density pumpout is non-linearly induced. Even though the density pumpout observed in modeling is not as large as in experiments, it is interesting to note that this configuration (with largest peeling-kink response) corresponds to the strongest ELM mitigation observed in experiments, characterized by the largest increase in ELM frequency and the largest density pumpout. Thus, the amplification of the resonant perturbation, induced by the coupling with peeling-kink modes excited by RMPs, probably generates an increased particle transport near the X-point responsible for ELM mitigation. These findings are in agreement with the results obtained in L-mode experiments and MARS-F calculations on MAST. As reported in [44], the density pump-out in the MAST L-mode experiments was found to be related to the displacement near the X-point calculated by MARS-F, which was also correlated with the amplitude of the resonant components, probably due to a similar coupling mechanism between the peeling-kink mode and the resonant perturbation as described here. The magnetic footprints calculated in the MARS-F field were also longer for the configuration with a large X-point displacement (in that case, the even parity one), and only in this configuration the footprints were actually experimentally observed [51].

It is not clear why the ergodicity affects the density profile rather than the temperature profile (which is affected by RMPs only in the pedestal) but it is rather consistent with the experimental observations where strong density pumpout is observed while temperature profile remains unchanged [28, 32]. Note that the width of the ergodic layer is probably too small at the outboard midplane to be observed in experiments, and the displacement near the X-point cannot be measured by diagnostics at this range of parameters. In other experiments at larger magnetic field, the displacement at the midplane could be precisely estimated from ECE diagnostic while rigidly rotating the applied magnetic perturbations [33], but the displacement near the X-point cannot be determined accurately at the moment.

Other modeling has been performed using data from the same shots: ideal MHD modeling with VMEC [52], resistive linear MHD modeling with MARS-F [24] and linear two-fluid resistive MHD modeling with M3D-C¹ [53]. Good qualitative agreement is found between these simulations and the simulations presented here, on the fact that the strongest ELM mitigation obtained in experiments (for $\Delta \Phi$ between $+{{90}^{\circ}}$ and $+{{60}^{\circ}}$) is correlated to the largest excitation of the peeling-kink response near the X-point and thus, to the largest displacement around the X-point. The coupling between the peeling-kink component with the resonant component was also found in MARS-F calculations [24], consistent with the results presented here.

The density pumpout induced by RMPs is however a feature that can only be reproduced in simulations including the non-linear coupling between the (n  =  2) perturbation and the (n  =  0) background profiles. This work shows that this pumpout can be partially explained by the ergodization and a large 3D-displacement near the X-point, yet some physical 'ingredients' such as the direct coupling in the density equation between the current flow along the perturbed field lines and the radial particle transport (as included e.g. in [16, 48]) will be introduced in the model in the near future in order to improve the modeling of the pumpout. Turbulence may also play a role on the increased particle transport induced by RMPs, which cannot be taken into account in this model. More elaborate heat and particle sources, as well as time-dependent heat and particle diffusive coefficients, will also be considered in future studies.

Present and future work will focus on the non-linear interaction between ELMs and RMPs, aiming to study the different couplings of n  =  2 RMPs with unstable modes, depending on the plasma response. In addition, comparisons with DIII-D cases at low collisionality will be made, since the ELM suppression at low collisionality observed in DIII-D [26, 27] is possibly correlated to the excitation of the peeling-kink response near the X-point, similarly to the observations in ASDEX Upgrade.

This work has been carried out within the framework of the EUROfusion Consortium and has received funding from the Euratom research and training programme 2014–2018 under grant agreement No 633053. Part of this work was carried out using the HELIOS supercomputer system at Computational Situational Centre of International Fusion Energy Research Centre (IFERC-CSC), Aomori, Japan, under the Broader Approach collaboration between Euratom and Japan, implemented by Fusion for Energy and JAEA. The work of Pavel Cahyna was supported by Czech Science Foundation grant 16-24724S and by MSMT CR, grant #8D15001. The views and opinions expressed herein do not necessarily reflect those of the European Commission.