Gyrokinetic investigation of the nonlinear interaction of Alfvén instabilities and energetic particle-driven geodesic acoustic modes

This paper presents a study of the interaction between Alfv (cid:2) en modes and zonal structures, considering a realistic ASDEX Upgrade equilibrium. The results of gyrokinetic simulations with the global, electromagnetic, particle-in-cell code ORB5 are presented, where the modes are driven unstable by energetic particles with a bump-on-tail equilibrium distribution function, with radial density gradient. Two regimes have been observed. At low energetic particle concentration, the zonal structure (identiﬁed as an energetic particle-driven geodesic acoustic mode) is more unstable than the Alfv (cid:2) en mode. In the regime at high energetic particle concentration, the Alfv (cid:2) en mode is more unstable than the zonal structure. The interplay between the modes leads to a modiﬁcation of their growth rates as well as to a modiﬁcation of their saturation levels. The theoretical explanation of the mode interaction is given in terms of three-wave coupling of the energetic particle-driven geodesic acoustic mode and Alfv (cid:2) en


Introduction
The next generation of fusion relevant machines (ITER [1,2], DEMO [3]) will be characterized by a large population of energetic particles (EPs).These are fusion products and charged particles generated by external power sources (like neutral beam injection, NBI).EPs are characterized by velocities much higher than the thermal velocity of the particles of the background plasma, v EP v th,i,e , where the subscripts i and e refer respectively to the main ion population and to the electrons present in the plasma.The transport and confinement of the EPs in Tokamaks and the physics related to EPs is of primary importance in order to achieve self heating plasma.In fact EPs can drive unstable, via wave-particle interactions, symmetry breaking electromagnetic perturbations whose presence can redistribute the EPs population, expelling them out of the plasma before they can thermalize [4].This consequently can lead to a less effective heating.In addition, the violent migration of the EPs toward the walls, caused by generated short but intense perturbations dubbed abrupt large-amplitude events [5], can possibly damage the machine.The understanding of the intensity of the fields of the saturated induced instabilities together with the transport of the EPs represent a key topic to reach fusion.
In fusion devices the EPs have a characteristic periodic motion (bounce frequency, transit frequency,...) of the same order of magnitude as that typical of the shear Alfvén waves (SAWs [6]).When the drive provided by the EPs exceeds the damping of these plasma fluctuations, a broad spectrum of Alfvén waves can be excited [7,8].These driven instabilities are classified into two types: Alfvén eigenmodes (AMs) and energetic particle modes (EPMs).The AMs are characterized by frequencies located inside the frequency gaps of the SAWs continuum spectrum, created by the field geometry and by plasma non-uniformities [9,10].The EPMs are non-normal modes of the SAWs continuum spectrum, emerging as discrete fluctuations at the frequency where wave-EP power is maximized [11].The linear and nonlinear dynamics of SAWs driven instabilities has been reviewed in Ref. [4].
Additionally, EPs can drive unstable modes with frequency comparable with that of the geodesic acoustic modes (GAM, [12]), characterized by (m, n) = (0, 0) scalar potential and (m, n) = (1, 0) up-down anti-symmetric density perturbation (being m and n respectively the poloidal and toroidal mode numbers).These driven modes are the energetic particle geodesic acoustic mode (EGAM), excited via free energy associated with velocity space gradients in the EPs distribution, as it has been shown analytically in Ref. [13].Studies on the nonlinear dynamics of EGAMs have been recently presented in [14,15].
The comprehension of the dynamics of the great zoology of modes present in plasma, their interaction and the redistribution of EPs, is crucial to understand the properties of burning plasma.The study of the saturation levels of these instabilities is fundamental to be able to be predictive regarding the intensity of the fields that will be present in future reactors.This motivate the interest in the study of the nonlinear evolution of these instabilities.Traditionally, the two main saturation mechanisms of such instabilities, have been recognized as A) reduction of the drive due to the redistribution of the EP population [16]; B) the transfer of energy to other modes via mode-mode coupling mediated by the thermal plasma nonlinearities [17].More recently, a new mechanism has also been studied, consisting in mode-mode coupling mediated by the EP nonlinearity.This mechanism has been found to be responsible, for example, of the excitation of zonal structures by Alfvén modes, and it has been called forced-driven excitation [18,19].
In the present paper, the interaction between AMs and EGAMs is studied showing the results of numerical simulations obtained with the nonlinear, gyrokinetic, electromagnetic, PIC code ORB5 [20].The dynamics that EGAM (m, n) = (0, 0) and dominant AM (m, n) = (2, 1) exhibit in numerical simulations where only one toroidal mode number is retained (n = {i} , i = 0, 1) is compared with the dynamics observed in simulations where both toroidal mode numbers are kept.The nonlinearities have been maintained only in the EPs dynamics, while the background particle species (Deuterium and electrons) follow their unperturbed trajectories.The EPs have a double-bump-on-tail distribution function in velocity space with radial density gradient.The background plasma species have Maxwellian distribution functions.The paper is structured as follows.In Sec.2, the model of ORB5 is presented.In Sec.3 the equilibrium in use in the simulations is shown, together with the simulations details.In Sec.4.1 we investigate the coupling of AMs and EGAMs in a simplified configuration (circular flux surfaces) to simplify the physics and more easily compare with analytical theory.In Sec.4.2, the analytical interpretation is provided.In Sec.4.3, the application to a more realistic configuration (ASDEX Upgrade experimental magnetic equilibrium) is described.Here the results of numerical simulations obtained with a realistic scenario, the so-called NLED-AUG case [21,22], are presented.This section represents an extension of the studies detailed in Ref. [23,24] and in Ref. [25].There, the dynamics, respectively, of EGAMs and AMs has been individually investigated, retaining only one toroidal mode number.In this work the main novelty is that, in order to study the interaction of the EGAM, with the dominant Alfvén mode (AM), (m, n) = (2, 1) (which is identified as an EPM), both the toroidal modes n = {0, 1} are retained in the performed simulations.

The Model
ORB5 [20] is a nonlinear, global, electromagnetic, particle in cell (PIC) code which solves the gyrokinetic Vlasov-Maxwell system of equations [26,27], accounting for the presence of collisions and sources.The code uses a system of straight field-line coordinates, (s, θ * , ϕ).The poloidal flux ψ, normalized at its value at the edge ψ 0 , plays the role of radial coordinate (s = ψ/ψ 0 , 0 ≤ s ≤ 1).ϕ is the toroidal angle while the poloidal magnetic θ * angle is defined as: being q(s) the safety factor profile, θ the geometric poloidal angle and B the background magnetic field, linked to the magnetic potential A 0 through the equation B = ∇ × A 0 In ORB5 all the physical quantities are normalized to four reference parameters: the mass and charge of the main ion species (m i and q i = eZ i , being e the elementary charge and Z i the atomic number of the i-th ion species), the values of the electron temperature at a radial location s 0 and of the magnetic field on-axis (T e (s 0 ) and B 0 respectively).The derived units are obtained from these four parameters.For example, the time is given in units of the inverse of the ion cyclotron frequency ω ci = q i B 0 /(m i c) (being c the speed of light in vacuum), the velocities are multiple of the ion sound velocity c s = T e (s 0 )/m i , and so on.
In ORB5 is possible to consider both analytical equilibrium, comprising circular magnetic surfaces, and ideal-MHD equilibria.This is a solution of the Grad-Shafranov equation, calculated through the CHEASE code [28].
The species distribution function f s is divided into a prescribed time-independent background distribution function F 0,s and into a time dependent part δf s .The subscript s refers to the particle species (that is s = i, e, f , respectively background ions, electrons and fast particles).The time dependent part of the distribution function is sampled with numerical particles, called markers, representing a portion of the phase space.The gyrokinetic Vlasov equation for the perturbed distribution function is: for the bulk ions, F 0 (R, E) only, while for a bump-on-tail F 0 (R, E, v ).Details can be found in Ref. [29].In Eq.2 R is the gyrocenter trajectories, v and v ⊥ are the parallel and perpendicular velocities respect to the background magnetic field.The equations of motions in mixed-variable formulation of the gyrocenter characteristics (R, E, v , µ) are [30]: where the terms proportional to δ are the nonlinear terms, corresponding to the perturbed equations of motion and ... is the gyro-average.Equations ( 3) to ( 5) are formulated in mixedvariables formulation [31].Here, the perturbed magnetic potential A has been split into its symplectic and hamiltonian parts: A = A s + A h .In Equations ( 3) to (5) the modified vector potential is present A * , the modified parallel magnetic filed B * and the unit vector in the direction of the magnetic filed, are present: The equations of motion are coupled with the following field equations: the quasineutrality condition, the parallel Ampère's law and the ideal Ohm's law: where: The integrals in eqs.( 8) and ( 9) represent respectively the mixed-variable gyrocenter density and the mixed-variable gyrocenter current.These integrals are calculated in the phase-space volume dW = B * dv dµdα (being α the gyro-phase).ρ is the particle gyroradius, while ρ s = √ m s T s /(q s B) is the thermal gyroradius.The equations are solved through the mixed-variable pullback algorithm, presented in Ref. [31], which has been able to mitigate the so-called cancellation problem, [32].The typical modes of interest are mainly aligned with the magnetic field line m nq(s).So a filter is applied to the Fourier coefficients of the perturbed density and current [33].In this way all the nonphysical modes introduced by charge and current deposition are filtered out.For each toroidal mode n ∈ [n min , n max ] only the poloidal modes m ∈ [−nq(s) ± ∆m], where ∆m is the width of the retained poloidal modes.
In the present paper we will study the nonlinear interaction of modes, analyzing simulations where only the EPs full dynamic is retained.This means that only in the EPs equations of motion, the nonlinear terms will be present, that is the terms proportional to δ in eqs.( 3) and (4).Electrons finite Larmor radius effects are neglected, given ρ e −→ 0.

Equilibrium
The discharge #31213@0.84shas been chosen as a base case for linear and nonlinear EPs simulations.The uniqueness of this scenario, the so called NLED-AUG case [21], is due to the fact that it exhibits a neutral beam (NB) induced fast-ion β comparable to that of the background plasma.In addition, the fast ions have energy 100 times larger than the thermal background.These unexplored corner of plasma parameters has been chosen to match the realistic ratios of plasma parameters that are going to be met in future fusion machines and to obtain a scenario where the transport of fast particles and the induced mode dynamics can be mainly attributed to the presence of the EPs being the effects of the background plasma minimized [22].This scenario is rich of nonlinear physics.A TAE burst is observed to trigger EGAMs, suggesting the coupling of these modes via the velocity space (EPs avalanches) and via mode-mode coupling processes.The great variety of nonlinear physics present here, together with the fact that the mode dynamics is mainly mediated by the EPs, makes of this an important scenario for the validation of theoretical tools and codes.
The safety factor profile has a reversed shear, with a minimum located at the radial position s 0.5, in the amount of 2.2.In Fig. 1 the q-profile of the NLED-AUG case is shown (blue curve).The background plasma temperature profiles and the electron density profile of the NLED-AUG case are shown in Fig. 2.   density profile and temperature profiles of the NLED-AUG case.The main ion species and the EPs are deuterium plasma.The equilibrium quasi neutrality is fulfilled by keeping constant the radial electron density profile and varying the EPs concentration, together with the deuterium concentration, satisfying: n e = i Z i n i (being n the density profile of the i-th species).The bulk plasma species (electrons and deuterium) have Maxwellian distribution functions, while the EPs have a double bump-on-tail distribution function (as in [23], [24], [29], [34]), because an anisotropy in velocity is needed to drive unstable EGAMs, [15]: The local maximum of the distribution function, in the velocity space, is located at v = v ,0 , while tval is the width in velocity space of the two shifted Maxwellians.In Sec.4.1, a circular magnetic equilibrium will be considered and will be characterized by a safety factor profile having a radial dependence close to that of the NLED-AUG case (see Fig. 1, orange dotted line).The EPs will have an on-axis radial density profile (see Fig. 2 on the right, the orange dotted curve).Through this approximation, we have that the EPs provide a contribution with constant sign to the linear growth rate of the driven modes, through the spatial derivative of their distribution function [8,25,35], being: In Eq.13, the poloidal flux ψ plays the role of radial coordinate, f i and q i are, respectively, the distribution function and charge of the i−th species of the plasma.In Sec.4.3 we will adopt the realistic ASDEX Upgrade magnetic equilibrium.The EPs will have an off-axis radial density profile (see Fig. 2 on the right, the green curve), as modelled by TRANSP [36].
4 Modification of the AM saturation in the presence of ZS

Basic physics of the AM/ZS interaction
In this section, we focus on the basic physics of the interaction between AMs and Zonal Structures (ZS) indicating, with this term, axysimmetric perturbations in general, that can be: Zero Frequency Zonal Flows (ZFZF), geodesic acoustic modes (GAM) or energetic particle driven geodesic acoustic modes (EGAM).To this aim, we consider here a simplified configuration, where the flux surfaces are circular and concentric.This allows us to neglect the secondary correction due to the geometry.
For computational reasons, in Sec.4.1 the electron mass has been then taken 500 times lighter than the ions mass m e = m i /500.In Tab.2, the values of some important parameters used in the simulations are shown.In Fig. 3 we show the simulation results obtained in two regimes: at low ( n EP / n e = 0.0379) and high ( n EP / n e = 0.114) EPs concentration (being .. the volume average).For each regime, we present the mode dynamics observed in simulations where only a single toroidal mode was retained, that is only n = {0} or n = {1}.These correspond respectively to the green and red curves in Fig. 3.We compare this evolution with the dynamics observed in simulations where both toroidal modes are present n = {0, 1} (blue and orange curves in Fig. 3).
The measured growth rates and frequencies will be provided in units of Alfvén frequencies measured on axis: ω A0 = 1 R 0 B 0 4πρ m,0 , being ρ m,0 the value of the background plasma density on axis: ρ m,0 = m i n i + m e n e .
In the regime at low EPs concentration, we observe the ZS (green curve in Fig. 3) to be more unstable than the dominant AM (red curve in Fig. 3).Their mode structure is shown in Fig. 4 and the measured growth rates and frequencies are reported in Tab.3.The ZS is identified as an EGAM while the AM (mainly peaked around s 0.5) has a frequency sitting slightly below the SAW continuum branch (m, n) = (2, 1).
Table 3: Growth rates and frequencies of the dominant modes in the simulations where only one toroidal mode number is retained.Low EPs concentration.
In the simulation where both the toroidal modes are present, we observe the EGAM dynamics (blue curve in Fig. 3) to be practically unaffected by the presence of the AM in the linear phase.On the contrary, the AM dynamics (orange curve in Fig. 3) appears to be driven by the EGAM.This results in an increase of the AM drive and of its frequency, lying now on the Alfvén continuum.

These values have been measured in the temporal domain t[ω −1
ci ] ∈ [28000; 39000]: These values of growth rate and frequency have been measured in the linear phase of the mode dynamics.They have been labelled with the superscript N L (standing for nonlinear ), to underline that these values arise from the nonlinear interaction between the AM and the EGAM, as opposed to what happened in the simulations where a single toroidal mode was present.Another qualitative proof of the mode-mode interaction is provided by Fig. 5. There, we show the observed mode structure for each toroidal mode number in simulations where both are present n = {0, 1}.We show the mode structure at a time before the saturation level is reached (around t ∼ 40000 ω −1 ci ).In Fig. 5 we observe the peaks of the electrostatic potential of (m, n) = (2, 1), located near the nodes of the electrostatic potential (m, n) = (0, 0).This suggests that the energy is going from the EGAM to the AM and the acting interaction mechanism is mode-mode coupling, as will be explained below.
In the regime at higher EP concentrations (Fig. 3) we observe the AM to be more unstable than the ZS: When the two modes are together, we observe the AM dynamics to be practically unchanged, instead the mode (m, n) = (0, 0), that here we identify as a ZFZF, has a growth rate: The AM pumps the ZS, with a forced-driven mechanism, that was analytically derived in Ref. [18].In Ref. [18] a pumping TAE with its complex conjugate was found to be responsible for the drive of the ZFZF, growing with: γ N L ZF ZF = 2γ AM .In Sec.4.2, we extend the analytical calculation of Ref. [18], to show that a general AM can force-drive a ZS and to explain the inverse mechanism, namely the excitation of an AM by a ZS, which is what is observed in the regime at low EPs concentration.

Wave-wave interaction mediated by the EPs
The theoretical explanation is given in terms of a three-wave coupling of the ZS and AM, mediated by the curvature-pressure coupling term of the EPs.We refer to Ref. [18], where the study of ZFZF driven unstable by TAE was carried out.There, the most unstable mode (a TAE) was shown to trigger the ZFZF, through wave-wave coupling and the ZFZF grows with twice the growth rate of the TAE δφ ZF ZF ∼ e 2γ T AE t .Following the derivation there presented, we want to generalize it to a general AM driving a ZS and to investigate the inverse problem, that is the triggering of an AM driven unstable through wave-wave coupling by an EGAM (being in this case the most unstable mode: γ L EGAM > γ L AM ), beating with the linear AM.To do so, we follow the theoretical derivation exposed in Ref. [18].However, here we start with the Vorticity equation for the AM (and not for the ZS), still neglecting the Reynolds and Maxwell stress, which express the nonlinearity coming from the background plasma (in the present paper only the non-linearity in the EPs are retained): In Eq.17, the first term represent the field line bending (FLB), the second term represents the Inertia Term (IT) and the third the Curvature Coupling Term (CCT).Here e s is the charge of the s-th species, ∂/∂ l is derivative parallel to the background magnetic field, δH m,n is the non-adiabatic response of the EPs associated to the scalar potential with polarization (m, n).
... denotes velocity-space integration and J k = J 0 (γ k ) is the Bessel function of 0-order with argument γ k = kρ, being k the mode wave-vector and ρ the gyro-radius.In Eq.17, the frequency drift ω d appears.It satisfies: where λ dZS and λ dAM define the coordinate transformation from the drift center to the particle gyro-center through the following relations: ω tr = v /(q R 0 ) is the transit frequency.The non-adiabatic response of the EPs is calculated through the nonlinear gyrokinetic equation: where e is the charge of the EPs, the subscript k is the wave-vector corresponding to the (m, n) helicity, δL = δφ − v c δA and Λ k = k=k +k b • k × k is the coupling term.The first term on the right hand side of Eq.20 represents the linear response, while the second gives the nonlinear response.Characterizing the AM component of the CCT with its toroidal and polodial mode numbers (m, n), we have that: This can be rewritten as: Using Eq.18, we can express Eq.22 as: The nonlinear non-adiabatic response of the EPs to the mode (m, n) is obtained from: On the right hand side of Eq.24, the term δH k is the linear response of the EPs to the mode with wave-vector k .This last couples two contributions of modes, expressed in δL k and δH k and the following relations must be satisfied: For the case investigated in the present paper, we want to study how an AM (2, 1) is driven unstable by an EGAM, characterized by a scalar potential dominated by its (0, 0) component, but having δH dominated by the mode numbers (±1, 0).Note that the mechanism described here is valid in general for an AM with any helicity (m AM , n AM ).Finally, with a similar calculation as that provided in Ref. [18], we obtain the following growth rate for the forced-driven AM: Where γ AM and γ EGAM are the linear growth rates of the AM and of the EGAM, and the scalar potential on the left-hand-side of Eq.26, refers to the nonlinearly generated AM, so: This explains how a marginally unstable Alfvén mode can be nonlinearly driven with a growth rate similar to that of the EGAM.

ASDEX Upgrade Equilibrium
In the present section, we show that the results obtained in the previous sections apply also to a realistic ASDEX Upgrade magnetic equilibrium, AUG shot #31213@0.84s.The plasma parameters and profiles are the same used in the previous sections, described in Sec.3.This shot belongs to a series of discharges performed in 2017 in ASDEX Upgrade [22].The plasma scenario has been obtained through plasma parameters previously unexplored in ASDEX Upgrade: β EP /β th 0.2 − 1 and E EP /T th 150 (the subscripts th refers to the thermal plasma).Through these, a scenario where the mode dynamics is dominated by the EPs has been obtained and the turbulence contribution minimized.In these discharges the interaction between EGAMs and AMs has been observed.In order to choose the proper discharge, the beam injection angle dependence of the EGAM has been investigated.The discharge #31213 has been selected to maximize the nonlinear EPs dynamics.This case, named the NLED-AUG case has become the base case for linear and non-linear EP simulations within three European theory projects [21].The NLED-AUG case is obtained with early off-axis NB heating, injected with an angle between the horizontal axis and the beam-line of 7.13 • .As it is shown in the portion of spectrogram in Fig. 6, an EPM-TAE burst is observed to trigger the EGAM.The experimental magnetic equilibrium measured at the time t = 0.84 s is considered in the simulations of the present section.EPs have an off-axis radial density profile (see Fig. 2).A scan against the EPs concentration is presented, 5% ≤ n EP / n D ≤ 10%... indicates the volume average.In this section the electrons have a realistic mass: m e m i /3676.In Fig. 7 the time evolution of the scalar potential is presented, for simulations with low and high EPs concentration.The dominant AM is an EPM, having the scalar potential dominant in its components (m, n) = (2, 1) and being mainly located around the radial position s 0.2 (as in Ref. [25]).
In the simulations where only a single toroidal mode is retained, the frequencies of the dominant modes in the linear phase are: corresponding respectively to ν ZS 45kHz and ν AM 96kHz.The ZS is identified with an EGAM and it results to be in good agreement with what was detailed in Ref. [24].The frequency of the AM lies below the branch of the continuum spectrum (m, n) = (2, 1) and its frequency is in agreement with the experimental results, lying in the EPM-TAE burst.
When both toroidal mode numbers are considered in the simulations, we observe the same tendency that was detailed in Sec.4.1.For high EPs concentration, we still observe the AM to drive the ZS, with a forced driven mechanism, before the saturation is reached.After that the first saturation of the AM is reached, the mode with the frequency of the EGAM is observed to develop.Below n EP / n D = 0.09, instead, we observe the inverse process.An EGAM drives the AM, as it was discussed in Sec.4.2.
In the regime below n EP / n D = 0.09, the difference of the saturation level of the AM with and without ZS is more pronounced (see Fig. 8) with respect to higher concentrations of EPs.We note that in the former regime, the ZS are linearly more unstable than the AMs, whereas in the latter the AMs are more unstable than the ZSs.The ZS in the former regime are identified as EGAMs.

Conclusion
The need to understand the nonlinear dynamics of EP-driven instabilities is motivated by the aim to be predictive about future scenarios that will be met in next generation fusion relevant machines.In this scenarios, an EP population is present due to external heating mechanisms or as a product of fusion reactions.Most of the previous theoretical works on the nonlinear dynamics of EPs driven instabilities has been focused on the independent evolution of single modes.The main saturation mechanisms considered in most of the previous studies was the EPs redistribution in real space or in phase space.
In Ref. [19], fully nonlinear magnetohydrodynamic simulations showed that the saturation amplitude of a toroidal Alfvén eigenmode (TAE) is reduced by the non-linear generation of zonal modes.This was analytically explained in Ref. [18], where the formation of a zero frequency zonal structure (ZFZS) was found to be forced-driven by a pumping TAE which couples with its complex conjugate.The induced EPs nonlinearities was found to be dominant over those of the background plasma, which were expressed in the Maxwell and Reynolds stress of the vorticity equation.In both these studies the ZS were not excited directly by the EPs.
This paper has been dedicated to the study of the interaction between two different type of modes driven unstable by the EPs: EGAMs and AMs.The results of numerical simulations obtained with the nonlinear, gyrokinetic, electromagnetic, PIC code ORB5 have been presented.In the present work only the EPs have been allowed to follow their full trajectories, redistributing in the phase space.The nonlinearities arising from the background plasma have been here neglected.The EPs are, in this case, the main responsible of the saturation and mode interactions.Since the thermal plasma nonlinearities are not considered, the same arguments exposed in Ref. [18] have been used to describe the interplay between the modes.The analytical derivation there exposed has been extended to the case here examined.
We observed that, in numerical simulations explained with analytical theory, through three wave coupling mediated by the EPs via the curvature pressure term in the vorticity equation, the nonlinear growth rate of the pumped mode is given by the sum of the linear growth rates of the pumping modes: γ N L n,m = γ L n ,m + γ L n ,m .EGAMs are found to channel energy into the AMs in the regime where the EGAMs linear growth rate is higher than the growth rate of the AM.Realistic ASDEX Upgrade plasma profiles have been considered.These results have been found both in analytical and experimental magnetic equilibrium.Despite the approximation of the EPs distribution function, these theoretical findings are qualitatively consistent with the experimental observations.This work opens the path to a new field of research, where the individual modes can saturate not only due to the EPs redistribution but also due to the mutual interaction.
In future works, simulations with more realistic distribution functions will be presented.Also all plasma nonlinearities will be taken into account, allowing to make comparisons with experimental data of ASDEX Upgrade, where an increased AM activity in the presence of strong EGAMs has also been observed experimentally.

Acknowledgments
Simulations presented in this work were performed on the CINECA Marconi supercomputer within the ORBFAST and OrbZONE projects.
One of the authors, F. Vannini, would like to thank Zhyiong Qiu for useful, interesting discussions and for great help provided in understanding the interaction between AM and ZS and

Figure 1 :
Figure 1: Radial dependence of the safety factor profile of the NLED-AUG case and of circular equilibrium (dashed line) discussed in the Sec.3.

Figure 2 :
Figure 2: Radial dependence of the temperature and density profiles.On the right, the radial density profiles of the EPs are shown, for the NLED-AUG case (green line) and for the present case (orange dotted line).

Figure 4 :
Figure 4: Radial dependence of the absolute value of n = {0} (left) and n = {1} (right) scalar potential components of the simulation at low EPs concentration (see Fig.3 on the left).Here only one toroidal mode number is retained.Color label correspond to different poloidal harmonics m.

Figure 5 :
Figure 5: Radial dependence of the absolute value of n = {0} (left) and n = {1} (right) scalar potential components of the simulation with n = {0, 1} at low EPs concentration (see Fig.3 on the left).Note that the peaks of the n = {1} ES potential correspond to the nodes of the n = {1} ES potential, consistently with a model of the nonlinear forcing of the n = {1} by the n = {0}.Color label correspond to different poloidal harmonics m.

Figure 6 :
Figure 6: Detail of the experimental spectrogram obtained with Mirnov Coil.The magnetic equilibrium at t = 0.84 s has been selected.

Figure 7 :
Figure 7: Modification of the AE dynamics in presence of the ZS, for low (left) and high (right) EPs concentration.

Figure 8 :
Figure 8: Modification of the saturation level of the AE, in presence of an EGAM.Note that strong modification of the AM saturation level in the regime of low EP concentrations, where the EGAM is stronger than the AM.On the right, the ratio of the saturation level of the AM in presence of the ZS against the corresponding value measured when the ZS is absent, is shown.

Table 1 :
482 1.666 2.202 1.055 • 10 8 2.7 • 10 −4 551.6 Constants in use: averaged minor radius, major radius, magnetic field on axis, ion cyclotron frequency, normalized plasma pressure.L x is the normalized size of the plasma system, defined as L x = 2 a 0 cs , being c s the ion sound velocity.In Sec.4.1 and Sec.4.3, we will show the results of simulations obtained with the radial electron

Table 2 :
Main simulations parameters (number of markers, time step, grid points).