Nonlinear dynamics of energetic-particle driven geodesic acoustic modes in ASDEX Upgrade
I. Novikau, A. Biancalani, A. Bottino, Ph. Lauber, E. Poli, P. Manz, G. D. Conway, A. Di Siena, N. Ohana, E. Lanti, L. Villard, ASDEX Upgrade Team
aa r X i v : . [ phy s i c s . p l a s m - ph ] M a r AIP/123-QED
Nonlinear dynamics of energetic-particle driven geodesic acoustic modes in ASDEXUpgrade.
I. Novikau, a) A. Biancalani, A. Bottino, Ph. Lauber, E. Poli, P. Manz, G. D.Conway, A. Di Siena, N. Ohana, E. Lanti, L. Villard, and ASDEX Upgrade Team b) Max-Planck-Institut für Plasmaphysik, 85748 Garching, Germany École Polytechnique Fédérale de Lausanne, Swiss Plasma Center, CH-1015 Lausanne,Switzerland (Dated: 30 March 2020)
Turbulence in tokamaks generates radially sheared zonal flows. Their oscillatory coun-terparts, geodesic acoustic modes (GAMs), appear due to the action of the magnetic fieldcurvature. The GAMs can be driven unstable by an anisotropic energetic particle (EP) pop-ulation leading to the formation of global radial structures, called EGAMs. The EGAMscan redistribute EP energy to the bulk plasma through collisionless wave-particle interac-tion. In such a way, the EGAMs might contribute to the plasma heating. Thus, investigationof EGAM properties, especially in the velocity space, is necessary for precise understand-ing of the transport phenomena in tokamak plasmas.In this work, the nonlinear dynamics of EGAMs without considering the mode interac-tion with the turbulence is investigated with the help of a Mode-Particle-Resonance (MPR)diagnostic implemented in the global gyrokinetic particle-in-cell code ORB5. An ASDEXUpgrade discharge is chosen as a reference case for this investigation due to its rich EP non-linear dynamics. An experimentally relevant magnetic field configuration, thermal speciesprofiles and an EP density profile are taken for EGAM chirping modelling and its com-parison with available empirical data. The same magnetic configuration is used to exploreenergy transfer by the mode from the energetic particles to the thermal plasma includingkinetic electron effects. For a given EGAM level the plasma heating by the mode can besignificantly enhanced by varying the EP parameters. Electron dynamics decreases theEGAM saturation amplitude and consequently reduces the plasma heating, even thoughthe mode transfers its energy to thermal ions much more than to electrons. b) See H. Meyer et al 2009 Nucl. Fusion a) Electronic mail: [email protected] . INTRODUCTION An energetic particle (EP) beam, injected in tokamak plasma, can excite a variety of modes.One of these modes, called energetic-particle-driven geodesic acoustic mode (EGAM) , can havea significant influence on the EP dynamics and on the plasma confinement. First of all, being ex-cited, this mode provides an additional mechanism of the energy exchange between the energeticparticles and the thermal plasma due to the wave-particle interaction. Moreover, because ofits interaction with the turbulence , this mode might be used as an additional knob in the turbu-lence regulation. Significant progress has been made in the last decade in building a theoreticalmodel which can explain the main nonlinear physics of the EGAM . Yet, when one wantsto quantitatively compare with experimental measurements, some of the approximations in theprevious models can be limiting. For example, due to the importance of the wave-particle reso-nances with the thermal and energetic species, a kinetic treatment should be used as in Ref. 12–14.Apart from that, effect of the magnetic surface shape can be considered by employing an exper-imental magnetic equilibrium. In this work, the importance of these effects on an experimentalASDEX-Upgrade (AUG) case by means of numerical simulations with the nonlinear gyrokineticcode ORB5 is investigated. ORB5 has been previously used for nonlinear studies of EGAMsin simplified configurations , and it is used here to compare the nonlinear EGAM dynamicswith experimental measurements.In this work, only the EGAM dynamics, leaving aside its interaction with the turbulence isgoing to be considered. The nonlinear behaviour of this mode, such as the mode chirping ,plasma heating by the EGAMs, and its dependence on the plasma parameters are of interest here.It should be noted here that in contrast to finite- n modes (where n is a toroidal number), whichalso can transfer energy from the EPs to the thermal plasma due to the wave-particle interaction,the EGAM dynamics is practically not associated with additional particle loss, at least if one doesnot consider the mode propagation and topological orbital changes .The EGAM dynamics is going to be investigated in the plasma configuration of the ASDEXUpgrade (AUG) discharge , where a richEP nonlinear dynamics is present . In this discharge, various types of EP-driven instabilities wereidentified, among which there are Alfvén instabilities (see, for example, Ref. 27) and EGAMs.In particular, a clear EGAM chirping was observed. More precisely, at the experimental EGAMfrequency spectrogram (Fig. 1a) one can see at least two different branches of the EGAM evolution2n time. A "short" one is with a ∼ . ∼
45 kHz till ∼ ∼
45 kHz till ∼
58 kHzand frequency saturation with a period of change around 13 ms. (a) (b)
FIG. 1: Experimental EGAM spectrum in the NLED AUG discharge , implemented in ORB5 to study field-particleinteraction in velocity space, is briefly described in Section II A. The mode chirping and corre-sponding comparison with the AUG experiment is presented in Section III. Using two shiftedMaxwellians as an EP distribution function, a nonlinear electrostatic (ES) gyrokinetic simulationwith the code ORB5 produces a relative EGAM up-chirping close to the experimental one. Inthis simulation, EP parameters have been taken consistent with previous works (Refs. 14 and 15),where in the same AUG discharge the ordering of the EGAM frequency was reproduced in linearelectrostatic and electromagnetic (EM) simulations, performed with the codes GENE and ORB5.In Section IV the EGAM linear behaviour is investigated by varying EP parameters. This studyis used later for the analysis of the plasma heating by the EPs through the EGAMs in Section V.Existence of an energy exchange between EPs and thermal particles are emphasized, and dominantresonant thermal species are indicated using the MPR diagnostic. Such kind of plasma heating hasbeen recently demonstrated in Ref. 6 in a realistic 3D equilibrium of the LHD stellarator, usingthe hybrid code MEGA with kinetic treatment of both thermal and energetic ions.By comparing ES simulations with adiabatic electrons and electromagnetic (EM) simulationswith drift-kinetic electrons, it is shown how the inclusion of the electron dynamics affects theEGAM formation and the plasma heating by this mode in Section VI. Although even in the pres-3nce of drift-kinetic electrons the EGAM transfers most of its energy to the thermal ions, electrondynamics still clearly decreases the mode level in the considered AUG discharge, and in such away reduces the energy flow from the EPs to the thermal ions. II. NUMERICAL AUG CONFIGURATION
The AUG plasma is simulated in the gyrokinetic (GK) particle-in-cell (PIC) global codeORB5 , where a realistic magnetic configuration is reconstructed from experimental datausing the Grad-Shafranov solver code CHEASE , without including the magnetic separatrix aslast closed magnetic flux surface (Fig. 1b). As it is described farther in this section, the radialdomain of the plasma configuration is slightly reduced for numerical reasons. The magnetic fieldat the magnetic axis is B = . R = .
67 m, the minor radiusis a = .
482 m. Realistic temperature and density profiles of the thermal species (deuteriumand electrons) are used as well (Figs. 2a-2b). Under "energetic particles" we understand here adeuterium beam, which actually drives EGAMs in a plasma system. This species is described by atwo-bumps-on-tail distribution function (also called as two shifted Maxwellians) in velocity space(Fig. 2c) F , EP ( p z , µ ) = A EP ( ψ ) exp " − m EP T EP (cid:18) p z m EP (cid:19) + µ Bm EP ! − v k , EP T EP cosh (cid:18) p z m EP v k , EP T EP (cid:19) , (1) A EP ( ψ ) = n EP ( ψ )( π ) / T / EP , (2)which leads to the EGAM formation through the wave-particle interaction (due to the mechanismof the inverse Landau damping). Here, p z is the canonical parallel momentum, µ is the magneticmoment, B is the background magnetic field, while v k , EP , T EP being constant input parameters,which specify a shift and width of the bumps respectively. It means, that the EPs are described bya flat temperature profile and a constant in space parallel velocity shift. From Fig. 2b one can seethat the EPs have relatively low temperature, which just indicates the fact that the width of the EPbumps is smaller in comparison to that of the thermal ion Maxwellian. A typical shift of the EPbumps in the parallel velocity is v k , EP [ c s ] = .
0, while typical EP temperature is T EP [ T e ( s = )] = . T EP [ keV ] = .
70. Here, c s = p T e , max / m d is the sound speed with m d being the deuteriummass, T e , max = .
15 keV. Such kind of normalization for v k , EP and T EP is used further in this work.4 a) (b) (c) FIG. 2: Thermal species density (Fig. 2a) and temperature (Fig. 2b) profiles (blue and red lines),and corresponding typical EP profiles (green lines). Distribution functions of the thermal andenergetic ions in velocity space are shown in Fig. 2c.The EPs are modelled by means of the two shifted Maxwellians to qualitatively investigate theEGAM behaviour in the AUG plasma configuration. The presence of two bumps in the velocityspace is explained by the necessity to avoid the input of an additional momentum into the plasmasystem, which might have some effect on the EGAM dynamics.In the nonlinear (NL) electrostatic simulation, performed in Section III to study the EGAM5p-chirping, the EPs have a realistic density profile, shown in Fig. 2a. In other sections, the EPsare described by an axisymmetric Gaussian distribution function: F , EP ( s ) ∼ exp (cid:18) − ( s − s EP ) σ EP (cid:19) , (3)and are localised at a radial point s EP = .
50, with s = p ψ / ψ edge being a radial coordinate, where ψ is the poloidal flux coordinate. Such a simplification has been done to have more freedom inthe variation of the EP parameters. The radial width of the EP Gaussian σ EP , used in this work, isequal to 0 .
10 to have a localised EP beam. An EP concentration is defined as n EP / n e , where n EP , n e are the EP and electron densities, averaged in volume, respectively.In the simulations, presented in this work, three species are taken into account: gyro-kineticthermal deuterium, gyro-kinetic energetic deuterium (EP), and electrons, either adiabatic (AE)or drift-kinetic (KE). Simulations with AE have been performed electrostatically, while the oneswith KE have been made electromagnetically. The experimental pressure ratio at the referencemagnetic surface has been adopted: β e = n e T e ( s = ) B / ( µ ) = . · − . (4)Electromagnetic simulations have the advantage of the absence of the unstable ω H -mode , ob-served in ES simulations with KEs. The KEs have a realistic deuterium/electron mass ratio m D / m e = has been used for the mitigation of the cancellationproblem . Real space in ES simulations has been discretized using the following grid param-eters: n s = , n χ = , n φ =
32, - for the radial, poloidal and toroidal directions respectively.Only n = | m | = [ , .., ] poloidal mode numbers have been taken into account, sincegenerally the geodesic acoustic modes have mainly low m numbers. The ES cases have been sim-ulated in a radial domain s = [ . , . ] with a time step dt [ ω − ci ] =
20 (where ω ci = eB / m d and e being the absolute value of the electron charge). The number of numerical markers for the thermaldeuterium ( N d ) and for the EPs ( N EP ) in these cases are N d = N EP = · . The EM cases havebeen simulated with n s = , dt [ ω − ci ] = s = [ . , . ] , taking N e = . · with the same N d and N EP . The radial domain has been reduced to make the EM simulationsmore stable by avoiding the abrupt increase of the safety factor at the edge, reducing in such away the restriction on the numerical time step. Since EGAMs are localised mainly in the core, thisreduction of the radial domain has a negligible influence on the EGAM dynamics.6 . Mode-Particle-Resonance diagnostic To localise velocity domains of the most intensive EGAM-particle interaction, the power bal-ance diagnostic, implemented in ORB5 , has been recently extended by keeping informationabout the wave-plasma interaction in velocity space. Using this diagnostic, a damping/growth rate γ of an ES wave can be calculated as γ = −
12 1 n w T w Z n w T w PE f d t , (5) P = ∑ sp P sp = − ∑ sp Z sp e Z d V d W sp f sp ˙R , sp · ∇∇∇ ( J , sp Φ ) , (6) E f = ∑ sp Z sp e Z d V d W sp f sp Φ . (7)where P is the energy transfer signal (heating rate), E f is the field energy, f sp ( r , p z , µ , t ) = F , sp ( r , p z , µ ) + δ f sp ( r , p z , µ , t ) is the species distribution function, Φ , J , sp Φ are the ES poten-tial perturbation and the gyroaveraged one, ˙R , sp is the species equations of motion, unperturbedby the field perturbations, Z sp e is the species charge. The time integration in Eq. 5 is performed onseveral wave periods n w T w . Finally, the integration in Eq. 6 and Eq. 7 is performed over the real V and species velocity W sp domains. For the equilibrium distribution functions F , sp , consideredin this work, the part of P related to F , sp is negligible. The energy transfer is normalized in thefollowing way: P = P [ W ] T e ( s = )[ J ] ω ci [ s − ] . (8) III. GK SIMULATION OF EGAM CHIRPING IN THE AUG DISCHARGE
To model the nonlinear evolution of the EGAM frequency, the EP parameters are taken consis-tent with previous works (Refs. 14 and 15), where the order of the mode frequency has been re-produced in linear GK simulations. The following EP velocity v k , EP = .
0, temperature T EP = . n EP / n e = .
095 are applied. In Fig. 3a, one can see a zoom of the ’long’branch, described previously, of the experimental EGAM up-chirping. In Fig. 3b, experimentaland numerical spectrograms normalized to an initial EGAM frequency are indicated. The bluecurve shows the mode up-chirping, obtained from a nonlinear ES simulation with adiabatic elec-trons. The numerical spectrogram has been measured by nonlinear fitting of zonal electric field indifferent time windows at the radial point of the mode localisation s = .
70. In the experimental7 a) (b) FIG. 3: A zoom of a ’long’ branch of the EGAM up-chirping from the experimental AUGspectrogram is shown in Fig. 3a. Numerical and experimental EGAM spectrograms,normalized to corresponding initial EGAM frequencies ( ∼
39 kHz and ∼
47 kHz respectively),are indicated in Fig. 3b. The numerical chirping is calculated by nonlinear fitting of zonal electricfield E at a radial point of the EGAM localisation (Appendix A).discharge, the mode was found in a radial domain s ∼ [ . , . ] according to soft X-ray emissiondata . The red curve corresponds to the experimental EGAM spectrogram, but with a squeezingtime scale 0 .
15, since the characteristic time of the numerical frequency change is shorter thanthe experimental one. As one can see from Fig. 3b, despite the two shifted Maxwellians, whichare used as a distribution function of the EPs, the GK model is able to simulate the EGAM rela-tive up-chirping in this AUG discharge. In Appendix A, the chirping is calculated at other radialpositions, and some convergence tests of this computation are presented. On the other hand, toperform quantitative comparison with the experiment of absolute chirping and of the time scales,one should use an EP distribution function closer to the experimental one such as a slow-downdistribution function with a pitch-angle dependence. Since the EGAM is a strongly driven mode,it significantly depends on EP parameters. Moreover, the considered simulation is performed elec-8 a) (b) FIG. 4: Dependence of the EGAM frequency (Fig. 4a) and linear growth rate (Fig. 4b) on thetemperature and velocity shift of the EP beam.trostatically without taking into account the dynamics of the drift-kinetic electrons. As it is shownin Section VI, electrons can significantly modify the EGAM behaviour and, as a result, might havesome effect on the mode chirping as well.
IV. INFLUENCE OF EP PARAMETERS ON THE EGAM DYNAMICS
To investigate the plasma heating by the EGAMs, we start from the study of the mode behaviourdependence on EP parameters. First of all, the linear growth rate is strongly modified by EPtemperature variations, as one can see in Fig. 4b, while the EGAM frequency remains practicallythe same and does not depend on T EP . As it was already shown in different works (Refs. 34 and35), the EPs can excite an EGAM corresponding to a damped mode in the absence of the EPs. Byreducing the EP parallel velocity by a factor around of two, an EGAM with half the frequency ofthe original is excited (compare blue points and green crosses in Fig. 4a). It should be noted thatsince the mode excited by a smaller EP velocity is less unstable, it was necessary to increase theEP concentration from n EP / n e ≈ .
01 to n EP / n e ≈ . m can be analytically estimated fromthe equation v ( m ) k , res = qR ω GAM m . (9)By reducing the EP velocity v k , EP one can excite an EGAM which is driven by higher order GAM-EP resonances. In Fig. 5a, one can see that by injecting an EP beam with v k , EP = . a) (b) FIG. 5: Localisation of the EGAM interaction with the EPs (Fig. 5a) and with thermal deuterium(Fig. 5b). Here, two simulations are considered: with v k , EP = . , T EP = . v k , EP = . , T EP = .
25 (red ones). The solid lines indicate the mode-species energy transfer P sp , summed on the perpendicular velocity and averaged on several EGAM periods in time. Thedashed lines depict the localisation of the species initial distribution functions.line), one obtains an EGAM driven by the low order resonance m = v k , EP = . m = m = m = v m = k , res , num ≈ . , v m = k , res , num ≈ . . (10)10 a) (b) FIG. 6: Time evolution of zonal electric field E in two cases with low (Fig. 6a) and high (Fig. 6b)EP velocities.A GAM frequency can be computed from a simulation without EPs: ω gam ( s = . ) = kHz , (11)to compare with corresponding EGAM frequencies ω EP . egam ( s = . ) = kHz , ω EP . egam ( s = . ) = kHz (12)Using these frequencies and a safety factor value q ( s = . ) = .
3, (E)GAM-plasma resonancepositions can be analytically estimated from Eq. 9: v EP . k , res , theor = qR ω EP . egam ≈ . , v EP . k , res , theor = qR ω EP . egam ≈ . , (13) v gam , m = k , res , theor ≈ . , v gam , m = k , res , theor ≈ . . (14)The difference between v EP . k , res , theor and v m = k , res , num might be explained by the fact that v m = k , res , num isestimated from a signal averaged in a whole space domain while v EP . k , res , theor is calculated by usinglocal plasma parameters. Moreover, the first order resonance v m = k , res , num has a significant widthalong parallel velocity. Since m = m = v m = k , res , num and v m = k , res , num respectively, it is reasonable to consider an EGAM as a mode that is driven by thefirst order resonance in case with v k , EP = .
0, and by the second order resonance in case with v k , EP = . m =
2) resonance in both cases,11ndependently of the EP parallel velocity (Fig. 5b). In other words, even when the EGAM is drivenby an EP beam with a high parallel speed, the EGAM is still damped by the thermal deuteriumthrough resonances at higher order. This can be simply explained by the fact that the position ofthese resonances are closer to the bulk of the thermal ion distribution function. Since the thermalion energy transfer signal (Fig. 5b) is positive, and the EP energy transfer signal (Fig. 5a) is mostlynegative, there is an energy flow from the EPs to the thermal ions establishing bulk plasma heating.
V. PLASMA HEATING BY EGAMS
As it has been already discussed previously in Refs. 3–6, the EGAMs might provide an addi-tional route for plasma heating by EPs. The thermal species can obtain energy from the EGAMdue to the Landau damping of the mode, which in turn is driven by the energetic particles throughthe inverse Landau damping. Note that this is also true for Alfvén modes, with the substantial dif-ference that the EP radial redistribution due to EGAMs is negligible with respect to Alfvén modes.Therefore, EGAMs represent a privileged mode for this plasma heating mechanism. This processis investigated here by varying the EP parameters in the AUG plasma configuration by performingnonlinear GK simulations firstly in the ES limit with adiabatic electrons.To calculate the amount of energy transferred from the mode to the thermal ions, the EGAM-species energy exchange signal P sp is integrated in time. As it is shown in Fig. 7, at the beginning,a rise of the heating rate due to the growth of the mode is observed. After a while, the modesaturates, the heating rate’s level decreases due to the relaxation of the EP distribution function, andin the deep saturated domain the energy transfer from the mode to the thermal ions is practicallysuppressed. It means that the total energy transferred to the bulk plasma (red points in Fig. 7)remains practically unchanged after the mode saturation. The amount of this energy dependson the EP parameters (Fig. 8a), which is directly reflected in the bulk ion temperature evolution(Fig. 8b).By varying the EP parameters in NL electrostatic simulations, one can observe a general corre-lation between the EGAM saturation levels and the amount of energy transferred from the modeto the bulk deuterium plasma (Fig. 9). However, for example, the case with v k , EP = . , T EP = . , n EP / n e = .
01 (the most right red star) has practically the same saturation level as the casewith v k , EP = . , T EP = . , n EP / n e = .
09 (the most right blue point). Having said that, weshould notice that the corresponding EGAM-bulk ion energy exchange for these cases is signif-12IG. 7: Time evolution of the heating rate (solid blue line), and the amount of energy transferredfrom the EGAM to the bulk deuterium plasma (red points). A case with n EP / n e = . , v k , EP = . , T EP = .
25 is considered here. The heating rate P D is normalizedaccording to Eq. 8. (a) (b) FIG. 8: Energy transferred from the EGAM to the bulk ions for different EP temperatures andparallel velocities (Fig. 8a). Evolution of the thermal ion temperature in the same cases (Fig. 8b).icantly different. This means that having the same mode level, one can achieve higher plasma13IG. 9: Dependence of the EGAM-thermal deuterium energy exchange on the EGAM saturationlevel in nonlinear ES simulations with adiabatic electrons. The saturation levels are calculated asa maximum value in time of r.m.s in space of the zonal electric field E . The electric field isnormalized to T e ( s = . )[ eV ] ω ci [ s − ] / c s [ m / s ] . Here, m spres indicates an order of a resonance,where the mode-species energy exchange occurs. The blue points correspond to the cases with v k , EP = . n EP / n e = . T EP = [ . , . , . , . ] , while the red stars correspond to thecases with v k , EP = . n EP / n e = . T EP = [ . , . , . , . ] .heating by varying EP parameters. In Appendix B, these two cases are analysed in detail. Here,in particular, increase of the EP density n EP / n e and lowering of the EP temperature T EP has led toa high energy transfer from the EPs to the mode. On the other hand, decrease of the EP velocityhas made the mode transfer higher part of its energy to the thermal plasma enhancing in such away the energy exchange between energetic and thermal ions, and keeping the EGAM amplitudeon the same level. 14 a) (b) FIG. 10: Influence of electron dynamics on the EGAM linear growth rate for different EPconcentrations (Fig. 10a). "AE" indicates ES simulations with adiabatic electrons, "KE" indicatesEM simulations with drift-kinetic electrons. Localisation of the EGAM interaction with theelectrons (case with n EP / n e = .
01, Fig. 10b). The black cone indicates the approximate positionof the boundary between the passing and trapped electrons of Eq. 15. The horizontal lines are theanalytical estimation of the EGAM-plasma resonances, calculated in Eq. 9 (for m = v k , EP = . , T EP = . VI. INFLUENCE OF ELECTRON DYNAMICS
In this section, the effect of the drift-kinetic electrons on the EGAM dynamics is examined,by considering electrons with a realistic mass m d / m e = v p − tr k = p εµ , (15)where ε is an inverse aspect ratio, the parallel velocity v k is normalized to the sound speed c s , andthe magnetic moment µ is normalised to m d c s / ( B ) . The energy exchange between the EGAMand trapped electrons occurs inside of this cone. The fact that the energy transfer resonancesare close to the passing-trapped boundary, means that the mode interacts mainly with the barelytrapped electrons with a bounce frequency close to that of the mode.15 a) (b) FIG. 11: Comparison between an electrostatic NL simulation assuming adiabatic electrons andan electromagnetic NL simulation with the drift-kinetic electrons. In both cases n EP / n e = . , v k , EP = . , T EP = . (a) (b) FIG. 12: Nonlinear ES and EM cases with n EP / n e = . , v k , EP = . , T EP = . VII. CONCLUSIONS
Energetic-particle-driven geodesic acoustic modes (EGAMs) have a significant influence onthe EP dynamics and on the plasma confinement. These modes provide an additional mechanismof the energy exchange between the energetic particles and the thermal plasma, enhancing directheating of the bulk ions. The significant progress done in the last decade in studying the EGAMnonlinear physics has been expanded here, by including kinetic electrons and the geometricaleffects of a realistic magnetic shape in the GK simulations.In this work, linear and nonlinear global GK simulations of the EGAMs in an ASDEX Up-grade discharge have been presented. While a corresponding tokamak configuration and thermalspecies profiles have been used here, the EPs have been represented by two shifted Maxwelliansin velocity space. Using such a model, an EGAM relative up-chirping close to the experimentalone has been obtained. It has demonstrated that the GK code is able to handle zonal structuresand anisotropic distribution functions properly in experimentally relevant cases which is a crucialingredient before going onto more complex plasma systems. On the other hand, to compare abso-lute mode chirping and characteristic time scales of the mode frequency evolution, it is necessaryto use an EP distribution function as close as possible to the experimental one.By varying the EP parameters, general correlation between the mode level and the EP-thermal17lasma energy flow has been revealed. Apart from that, it has been emphasized that the plasmaheating by EGAMs benefits from the presence of high order mode-particle resonances which areoften responsible for the EGAM-thermal species interaction. More precisely, it has been shownthat having the same mode levels, one can achieve a higher energy exchange between energeticand thermal species through an EGAM by adjusting EP parameters. Finally, it was indicated thatalthough the EGAM transfers most of its energy to the thermal ions, and not to the electrons, theelectron dynamics might significantly reduce the plasma heating by EGAMs by lowering the modeamplitude.Since the nonlinear GK simulations with drift-kinetic electrons have shown a significant influ-ence of the electron dynamics on the EGAM behaviour, it would be reasonable to study effect ofthe kinetic electrons on the mode chirping. The results obtained here will be useful in a futurework to perform investigation of the EGAM influence on the turbulent dynamics in the ASDEXUpgrade, that will be a next step in the code validation and an interesting example of a complexplasma system.
ACKNOWLEDGMENTS
The authors would like to acknowledge stimulating discussions with the Multi-scale Energeticparticle Transport in fusion devices (MET) team, and especially with Prof. Fulvio Zonca. Thiswork has been carried out within the framework of the EUROfusion Consortium and has receivedfunding from the Euratom research and training program 2014-2018 and 2019-2020 under grantagreement N o Appendix A CONVERGENCE TESTS
Convergence tests of some simulations discussed above are presented here. The main attentionis dedicated to the mode saturation levels and the EGAM chirping. First of all, the mode frequencyevolution considered in Sec. III is analysed. That computation has been done with dt = , n s = , N d = N EP = . · , which has a finer space grid (and, respectively, a higher number ofmarkers) than a typical electrostatic simulation described in Sec. II. The radial structure of the18ode (of the zonal electric field) is shown in Fig. 13a, where one can see that the mode is localisedaround a radial point s = .
70. In Fig. 13b, the EGAM chirping is investigated at different radialpositions near the mode localisation. There is some dependence of the frequency evolution on theradial position. The point s = .
70 is chosen since the chirping there is closer to the frequencyevolution of a signal averaged around the mode localisation (black markers in Fig. 13b). The modefrequency has been calculated by using a nonlinear fitting of the zonal electric field (Fig. 13c) to atest function ∼ cos ( ω t ) exp ( γ t ) (A.1)Time evolution of the signal is split on a set of time intervals of a length ∆ t with several EGAMoscillations. By performing the nonlinear fitting in every time domain, one can find a mode fre- (a) (b) (c) (d) FIG. 13: EGAM radial structure and its time evolution at s = .
70 (Figs. 13a and 13crespectively). Fig. 13b: EGAM chirping at different radial points, where every curve isnormalized to its own initial frequency. The black markers indicate frequency evolution of anaveraged zonal electric field in s = [ . , . ] interval. Here, ∆ t = .
11. Fig. 13d: frequencyevolution at s = .
70 with different time interval lengths ∆ t for the nonlinear fitting.19 a) (b) FIG. 14: Variation of the EP density profiles (Fig. 14a). Fig. 14b: corresponding mode chirpingcalculated at s = . ∆ t . The case with ∆ t = .
11 is chosen since the corresponding chirpinghas one of the smallest errorbars.By slightly varying the localisation of the energetic particles (Fig. 14a), we analyse how themode chirping changes due to the uncertainty in the definition of the EP density profile. InFig. 14b, one can see that the chirping does not disappear and remains more less stable withrespect to the indicated density variation.Convergence tests of some nonlinear ES simulations are shown in Fig. 15. Computation of theEGAM chirping with different numerical parameters are presented in Fig. 15a, where one can seethat the simulation is converged. As an additional test, we present here a convergence test (withrespect to a number of EP markers N EP ) of ES nonlinear simulations with different EP velocities(Fig. 15b and Fig. 15c). In Sec. V, all simulations have N EP = e
7. Increase of N EP by a factor of10 does not change T D growth rate and T D saturation level.Convergence tests of the nonlinear EM case with drift-kinetic electrons, described in Sec. VI,are presented in Fig. 16, where one can see that the EM simulation is converged.20 a) (b) (c) FIG. 15: EGAM chirping in simulations with different numerical parameters is shown inFig. 15a. Thermal plasma heating with different number of EP markers for two cases: case v k , EP = . , T EP = .
15 (Fig. 15b), case v k , EP = . , T EP = .
40 (Fig. 15c).21 a) (b) FIG. 16: Convergence tests for the nonlinear EM case with drift-kinetic electrons, described inSec. VI.22 a) (b) (c) (d) FIG. 17: Time evolution of energy transfer signals between an EGAM and thermal plasma(Fig. 17a), and between the mode and energetic particles (Fig. 17c). These signals are integratedin a real and velocity spaces. Localisations of the EGAM-thermal deuterium resonances in casewith v k , EP = . v k , EP = .
5) shown in Fig. 17b are obtained from the corresponding energytransfer signals integrated along a perpendicular velocity and time averaged in a blue (red) timeintervals. The same procedure is applied to find the mode - EP resonances in Fig. 17d.
Appendix B DETAILED ENERGY TRANSFER ANALYSIS OF EGAM DYNAMICS
It has been mentioned in Sec. V that energetic particles of significantly different energies(parallel velocities) can drive EGAMs with similar saturation levels, but with a different effi-ciency in the energy transfer from the EPs to thermal plasma. Here, we are going to considerthis process in detail on an example of two simulations. In the first one, an energetic beam has v k , EP = . , T EP = . , n EP / n e = .
09, while in the second simulation, the EPs are representedby a v k , EP = . , T EP = . , n EP / n e = .
01 beam. The time evolution of energy transfer signalsbetween EGAM and thermal and energetic ions are shown in Fig. 17a and Fig. 17c, respectively.23 a) (b) (c) (d) FIG. 18: EGAM radial structures are shown in Fig. 18a (case with v k , EP = .
0) and Fig. 18c (casewith v k , EP = . s = . v k , EP = .
5, is less effective. To compensate that, the EP concentration has been increased forthe v k , EP = . v k , EP = . v k , EP = . v k , EP = .
5, the mode grows slower, and theflattening of the EP distribution function develops slower as well (Fig. 19) than with v k , EP = . a) (b) FIG. 19: Evolution in time of EP distribution functions for the case v k , EP = . , T EP = . v k , EP = . , T EP = .
40 (Fig. 19b). (a) (b)
FIG. 20: EGAM-electron resonance in a nonlinear EM case (Fig. 20a). The horizontal linesindicate an estimated position of the mode-particle m = m = m = s = .
47 having a linear frequency ω egam [ ω ci ] ≈ · − . v EP = .
0, described in Sec. VI. In Fig. 20a, the EGAM-electron resonances are localised near anestimated electron passing-trapped boundary (Eq. 15), which is consistent to the linear computa-tion shown in Fig. 10b. Apart from that, as it has been observed in the ES cases discussed above,the mode here is driven through the m = m = EFERENCES C. Boswell, H. Berk, D. Borba, T. Johnson, S. Pinches, S. Sharapov,Observation and explanation of the jet n=0 chirping mode, Physics Letters A 358 (2) (2006)154 – 158 (2006). doi:https://doi.org/10.1016/j.physleta.2006.05.030 .URL G. Y. Fu, Energetic-particle-induced geodesic acoustic mode, Phys. Rev. Lett. 101 (2008)185002 (Oct 2008). doi:10.1103/PhysRevLett.101.185002 .URL https://link.aps.org/doi/10.1103/PhysRevLett.101.185002 M. Sasaki, K. Itoh, S.-I. Itoh, Energy channeling from energetic particles to bulk ions via beam-driven geodesic acoustic modes—GAM channeling,Plasma Physics and Controlled Fusion 53 (8) (2011) 085017 (jun 2011). doi:10.1088/0741-3335/53/8/085017 .URL https://doi.org/10.1088%2F0741-3335%2F53%2F8%2F085017 M. Osakabe, T. Ido, K. Ogawa, A. Shimizu, M. Yokoyama, R. Seki, C. Suzuki, M. Isobe, K. Toi,D. A. Spong, K. Nagaoka, Y. Takeiri, H. Igami, T. Seki, K. Nagasaki, LHD experiment group,Indication of bulk-ion heating by energetic particle driven geodesic acoustic modes on lhd, in:25th IAEA fusion energy conference, 2014 (2014).URL https://inis.iaea.org/search/search.aspx?orig_q=RN:47070985 K. Toi, K. Ogawa, T. Ido, S. Morito, S. Ohdachi, N.A. Pablant, K. Tanaka, H. Funaba,T. Oishi, M. Osakabe, A. Shmizu, H. Tsuchiya, K.Y. Watanabe, LHD Experiment Group,Observation of non-collisional bulk ion heating by energetic ion driven geodesic acoustic modes in lhd,in: Proc. 16th IAEA Technical Meeting on Energetic Particles in Magnetic Confinement Sys-tems — Theory of Plasma Instabilities, Shizuoka, Japan, 2019 (2019).URL https://nucleus.iaea.org/sites/fusionportal/Shared%20Documents/EPPI%202019%20Material.pdf H. Wang, Y. Todo, M. Oasakabe, T. Ido, Y. Suzuki,Simulation of energetic particle driven geodesic acoustic modes and the energy channeling in the large helical device plasmas,Nuclear Fusion 59 (9) (2019) 096041 (aug 2019). doi:10.1088/1741-4326/ab26e5 .URL https://doi.org/10.1088%2F1741-4326%2Fab26e5 D. Zarzoso, P. Migliano, V. Grandgirard, G. Latu, C. Passeron,Nonlinear interaction between energetic particles and turbulence in gyro-kinetic simulations and impact on turbulence properties,Nuclear Fusion 57 (7) (2017) 072011 (jun 2017). doi:10.1088/1741-4326/aa7351 .URL https://doi.org/10.1088%2F1741-4326%2Faa7351 M. Sasaki, K. Itoh, K. Hallatschek, N. Kasuya, M. Lesur, Y. Kosuga, S.-I. Itoh,Enhancement and suppression of turbulence by energetic-particle-driven geodesic acoustic modes,Scientific Reports 7 (2017) 16767 (2017). doi:10.1038/s41598-017-17011-y .URL https://doi.org/10.1038/s41598-017-17011-y Z. Qiu, I. Chavdarovski, A. Biancalani, J. Cao, On zero frequency zonal flow and second harmonic generation by finite amplitude energetic particle induced geodesic acoustic mode,Physics of Plasmas 24 (7) (2017) 072509 (2017). arXiv:https://doi.org/10.1063/1.4993053 , doi:10.1063/1.4993053 .URL https://doi.org/10.1063/1.4993053 Z. Qiu, L. Chen, F. Zonca, Kinetic theory of geodesic acoustic modes in toroidal plasmas: a brief review,Plasma Science and Technology 20 (9) (2018) 094004 (jul 2018). doi:10.1088/2058-6272/aab4f0 .URL https://doi.org/10.1088%2F2058-6272%2Faab4f0 L. Chen, Z. Qiu, F. Zonca, Short wavelength geodesic acoustic mode excitation by energetic particles,Physics of Plasmas 25 (1) (2018) 014505 (2018). arXiv:https://doi.org/10.1063/1.5003142 , doi:10.1063/1.5003142 .URL https://doi.org/10.1063/1.5003142 H. S. Zhang, Z. Lin, Trapped electron damping of geodesic acoustic mode, Physics of Plas-mas 17 (7) (2010) 072502 (2010). arXiv:https://doi.org/10.1063/1.3447879 , doi:10.1063/1.3447879 .URL https://doi.org/10.1063/1.3447879 G. Merlo, S. Brunner, Z. Huang, S. Coda, T. Görler, L. Villard,A. B. Navarro, J. Dominski, M. Fontana, F. Jenko, L. Porte, D. Told,Investigating the radial structure of axisymmetric fluctuations in the TCV tokamak with local and global gyrokinetic GENE simulations,Plasma Physics and Controlled Fusion 60 (3) (2018) 034003 (jan 2018). doi:10.1088/1361-6587/aaa2dc .URL https://doi.org/10.1088%2F1361-6587%2Faaa2dc I. Novikau, A. Biancalani, A. Bottino, A. Di Siena,P. Lauber, E. Poli, L. Villard, N. Ohana, S. Briguglio,Implementation of energy transfer technique in orb5 to study collisionless wave-particle interactions in phase-space,Comp. Phys. Comm. (2019). doi:10.1016/j.cpc.2019.107032 .URL A. Di Siena, A. Biancalani, T. Görler, H. Doerk, I. Novikau, P. Lauber, A. Bottino, E. Poli, The27SDEX Upgrade Team, Effect of elongation on energetic particle-induced geodesic acousticmode, Nucl. Fusion 58 (2018) 106014 (2018). doi:10.1088/1741-4326/aad51d . S. Jolliet, A. Bottino, P. Angelino, R. Hatzky, T. Tran, B. Mcmillan, O. Sauter, K. Ap-pert, Y. Idomura, L. Villard, A global collisionless pic code in magnetic coordinates,Computer Physics Communications 177 (5) (2007) 409 – 425 (2007). doi:https://doi.org/10.1016/j.cpc.2007.04.006 .URL A. Bottino, E. Sonnendruecker, Monte carlo particle-in-cell methods for the simulation of thevlasov–maxwell gyrokinetic equations, Journal of Plasma Physics 81 (5) (2015) 435810501(2015). doi:10.1017/S0022377815000574 . E. Lanti, N. Ohana, N. Tronko, T. Hayward-Schneider, A. Bottino, B. McMil-lan, A. Mishchenko, A. Scheinberg, A. Biancalani, P. Angelino, S. Brun-ner, J. Dominski, P. Donnel, C. Gheller, R. Hatzky, A. Jocksch, S. Jolliet,Z. Lu, J. M. Collar, I. Novikau, E. Sonnendrücker, T. Vernay, L. Villard,Orb5: A global electromagnetic gyrokinetic code using the pic approach in toroidal geometry,Computer Physics Communications (2019) 107072 (2019). doi:https://doi.org/10.1016/j.cpc.2019.107072 .URL A. Biancalani, I. Chavdarovski, Z. Qiu, A. Bottino, D. Del Sarto, A. Ghizzo, Ö. D. Gür-can, P. Morel, I. Novikau, Saturation of energetic-particle-driven geodesic acoustic modes dueto wave–particle nonlinearity, Journal of Plasma Physics 83 (6) (2017) 725830602 (2017). doi:10.1017/S0022377817000976 . A. Biancalani, N. Carlevaro, A. Bottino, G. Montani, Z. Qiu, Nonlinear velocity re-distribution caused by energetic-particle-driven geodesic acoustic modes, mapped withthe beam-plasma system, Journal of Plasma Physics 84 (6) (2018) 725840602 (2018). doi:10.1017/S002237781800123X . H. Berk, C. Boswell, D. Borba, A. Figueiredo, T. Johnson, M. Nave, S. Pinches, S. Sharapov,JET EFDA contributors, Explanation of the JETn= 0 chirping mode, Nuclear Fusion 46 (10)(2006) S888–S897 (sep 2006). doi:10.1088/0029-5515/46/10/s04 .URL https://doi.org/10.1088%2F0029-5515%2F46%2F10%2Fs04 H. Berk, C. Boswell, D. Borba, A. Figueiredo, T. Johnson, M. Nave, S. Pinches, S. Shara-pov, JET EFDA contributors, Explanation of the JETn= 0 chirping mode, erratum, Nuclear Fu-28ion 50 (4) (2010) 049802 (apr 2010). doi:10.1088/0029-5515/50/4/049802 .URL https://doi.org/10.1088%2F0029-5515%2F50%2F4%2F049802 D. Zarzoso, D. del Castillo-Negrete, D. Escande, Y. Sarazin,X. Garbet, V. Grandgirard, C. Passeron, G. Latu, S. Benkadda,Particle transport due to energetic-particle-driven geodesic acoustic modes, Nuclear Fusion58 (10) (2018) 106030 (aug 2018). doi:10.1088/1741-4326/aad785 .URL https://doi.org/10.1088%2F1741-4326%2Faad785 P. Lauber, AUG test case description (2015).URL P. Lauber, B. Geiger, G. Papp, G. Por, L. Guimarais, P.ZS. Poloskei, V. Igochine,M. Maraschek, G. Pokol, T. Hayward-SCHNEIDER, Z. Lu, X. Wang, A. Bot-tino, F. Palermo, I. Novikau, A. Biancalani, G. Conway, THE ASDEX UPGRADETEAM, THE EUROFUSION ENABLING RESEARCH ’NAT’ and ’NLED’ TEAMS.,Strongly non-linear energetic particle dynamics in asdex upgrade scenarios with core impurity accumulation,in: Proc. 27th IAEA Fusion Energy Conference - IAEA, Gandhinagar, India, 2018 (2018).URL https://conferences.iaea.org/indico/event/151/papers/6094/files/5080-IAEA_proceedings_Lauber_v3.pdf L. Horváth, G. Papp, P. Lauber, G. Por, A. Gude, V. Igochine, B. Geiger,M. Maraschek, L. Guimarais, V. Nikolaeva, G. Pokol, The ASDEX Upgrade Team,Experimental investigation of the radial structure of energetic particle driven modes, Nu-clear Fusion 56 (11) (2016) 112003 (jul 2016). doi:10.1088/0029-5515/56/11/112003 .URL https://doi.org/10.1088%2F0029-5515%2F56%2F11%2F112003 F. Vannini, A. Biancalani, T. Hayward-Schneider, P. Lauber, A. Mishchenko, I. Novikau, E. Poli,the ASDEX Upgrade team., submitted to Physics of Plasmas (2019). [link].URL https://arxiv.org/abs/1910.14489 F. Jenko, W. Dorland, M. Kotschenreuther, B. N. Rogers, Electron temperature gradient driventurbulence, Phys. Plasmas 7 (2000) 1904–1910 (May 2000). doi:10.1063/1.874014 . H. Luetjens, A. Bondeson, O. Sauter, The chease code for toroidal mhd equilibria,Computer Physics Communications 97 (3) (1996) 219 – 260 (1996). doi:https://doi.org/10.1016/0010-4655(96)00046-X .URL W. Lee, Gyrokinetic particle simulation model, Journal of Computational Physics 72 (1) (1987)243 – 269 (1987). doi:https://doi.org/10.1016/0021-9991(87)90080-5 .29RL A. Mishchenko, A. Bottino, R. Hatzky, E. Sonnendruecker, R. Kleiber, A. Koenies,Mitigation of the cancellation problem in the gyrokinetic particle-in-cell simulations of global electromagnetic modes,Physics of Plasmas 24 (8) (2017) 081206 (2017). arXiv:https://doi.org/10.1063/1.4997540 , doi:10.1063/1.4997540 .URL https://doi.org/10.1063/1.4997540 R. Hatzky, A. Koenies, A. Mishchenko, Electromagnetic gyrokinetic pic simulation with an adjustable control variates method,Journal of Computational Physics 225 (1) (2007) 568 – 590 (2007). doi:https://doi.org/10.1016/j.jcp.2006.12.019 .URL N. Tronko, A. Bottino, E. Sonnendruecker, Second order gyrokinetic theory for particle-in-cell codes,Physics of Plasmas 23 (8) (2016) 082505 (2016). arXiv:https://doi.org/10.1063/1.4960039 , doi:10.1063/1.4960039 .URL https://doi.org/10.1063/1.4960039 J.-B. Girardo, D. Zarzoso, R. Dumont, X. Garbet, Y. Sarazin, S. Sharapov,Relation between energetic and standard geodesic acoustic modes, Physics of Plas-mas 21 (9) (2014) 092507 (2014). arXiv:https://doi.org/10.1063/1.4895479 , doi:10.1063/1.4895479 .URL https://doi.org/10.1063/1.4895479 D. Zarzoso, A. Biancalani, A. Bottino, P. Lauber, E. Poli, J.-B. Girardo, X. Garbet, R. Dumont,Analytic dispersion relation of energetic particle driven geodesic acoustic modes and simulations with NEMORB,Nuclear Fusion 54 (10) (2014) 103006 (sep 2014). doi:10.1088/0029-5515/54/10/103006 .URL https://doi.org/10.1088%2F0029-5515%2F54%2F10%2F103006https://doi.org/10.1088%2F0029-5515%2F54%2F10%2F103006