Impact of supra-thermal particles on plasma performance at ASDEX Upgrade with GENE-Tango simulations

This paper presents global gyrokinetic simulations on the transport time scale of an ASDEX Upgrade H-mode discharge showing a pronounced peaking of the on-axis ion temperature profiles. Leveraging the newly developed GENE-Tango tool, which combines the global gyrokinetic code GENE with the transport solver Tango, we investigate the impact of energetic particles and electromagnetic effects on the improved plasma performance observed in the experimental discharge. Our results reveal that a striking agreement between the GENE-Tango simulations and the experimental measurements can be achieved only when energetic particles and electromagnetic effects are simultaneously retained in the modeling. In contrast, when these are neglected we observed a significant underestimation of the on-axis ion temperature, aligning with profiles computed using TGLF-ASTRA. The peaking in the ion temperature profile observed in the simulations can be attributed to the effective suppression of turbulence by high-frequency electromagnetic modes, likely Kinetic Ballooning Modes/Alfvén eigenmodes. These modes play a critical role in enhancing zonal flow activity and shearing rate levels which thus lead to a localized increase in the temperature gradient. However, it is crucial to maintain these modes at a state of marginal stability or weak instability to prevent energetic particle turbulence destabilization. Otherwise, the result would be a flattening of all the thermal profiles. Interestingly, we found that global GENE-Tango simulations are required to model correctly the linear dynamics of these high-frequency modes. Additionally, global simulations demonstrate greater tolerance than flux-tube simulations for marginal instability of these high frequency modes while maintaining power balance agreement.


Introduction
To support the quest for sustainable and efficient fusion energy production, it is crucial to develop sophisticated tools that can accurately capture the main physical effects occurring in plasma experiments.These tools will enable reliable predictions of plasma profiles, facilitating the investigation of various plasma regimes and assessing optimal conditions for maximizing fusion performance while avoiding unexpected confinement degradation.
In recent years, there has been significant progress in developing reduced turbulence models that efficiently capture the essential features of plasma micro-instabilities and turbulence.Notable examples of such models include TGLF (trappedgyro-Landau-fluid) [1,2] and QuaLiKiz [3][4][5][6].These models offer the advantage of fast computation times, enabling turbulent flux calculations-based on saturation rules-within seconds on a single CPU.They have been successfully integrated with transport codes like TGLF-ASTRA [7,8], TGLF-TGYRO [9], and QuaLiKiz-JINTRAC [10].However, despite the ability of these reduced turbulence models to accurately reproduce experimental profiles across a wide range of plasma conditions, they still lack the capability to capture nonlinear energetic particle effects on core plasma turbulence [11][12][13][14].They should be able to capture some linear fast ion effects [15].As a result, they tend to underestimate the on-axis ion temperature in regimes where energetic particles play a significant role.This limitation might be particularly severe in future burning plasma experiments, impacting both the initial ramp-up phase, which relies on auxiliary systems to reach the optimal conditions for fusion reactions, as well as the subsequent operational phase where heating is solely reliant on MeV alpha particles.
Nonlinear gyrokinetic codes have proven effective in capturing energetic particle effects on turbulence in various plasma regimes [16][17][18][19][20], offering valuable insights into the underlying physical mechanisms that lead to observed ionscale turbulence suppression in experiments.However, these codes are computationally significantly more demanding than quasi-linear models, despite extensive work to improve code performances, e.g.via GPU porting [21].As a result, gyrokinetic simulations are often limited to single flux-tubes, and even in the few global simulations conducted so far, plasma profiles are not allowed to freely evolve due to the combined effect of turbulence and external sources, but are kept fixed by using artificial Krook-like operators.Recent advancements have successfully addressed these limitations [22][23][24][25].A notable example is the coupling of the global gyrokinetic code GENE with the transport solver Tango [23,25].This coupling takes advantage of the distinct time-scale separation between microscopic and macroscopic physics, resulting in a substantial reduction in computational costs allowing, -for the first time-high-fidelity simulations that extend up to the transport time scale.In this coupling, the turbulence code determines the turbulence fluxes, while the transport code extracts transport coefficients based on the plasma profiles and evolves the plasma profiles themselves considering sources and turbulent fluxes.This approach ensures that the heavily computational gyrokinetic code is no longer constrained by the micro-turbulence time scale, which is several orders of magnitude shorter than the macroscopic scale at which background profiles evolve.Currently, there is an ongoing effort [26] to investigate potential differences between stand-alone gyrokinetic simulations up to the transport time-scales and gyrokinetic codes coupled to transport codes to speed-up the calculations, such as GENE-Tango.However, due to the immense computational cost of stand-alone gyrokinetic codes up to the transport time scale, this work has been confined to scenarios involving adiabatic electrons and simplified geometries, neglecting any energetic particle species.Future endeavors aim to broaden this investigation to include kinetic electrons.
In our study, we present the first GENE-Tango simulations of an ASDEX Upgrade H-mode discharge with large fast ion content, demonstrating a substantial improvement in plasma performance attributed to energetic particles.Remarkably, we find that only when the energetic particle species is consistently included in the radially global nonlinear electromagnetic GENE-Tango simulations, we observe an excellent agreement with the experimental plasma profiles.Removing the fast ions from the model leads to a significant flattening of the ion temperature profile, aligning with profiles computed using TGLF-ASTRA, which significantly under-predicts the on-axis ion temperature.Detailed analysis reveals that the thermal ion temperature gradient (ITG) increases markedly in correspondence with Kinetic Ballooning Modes (KBMs)/Alfvén eigenmodes (AEs).These modes enhance zonal flow activity and significantly enhance shearing rate levels compared to GENE-Tango simulations without energetic particles, resulting in turbulence suppression consistent with previous flux-tube simulations.When these modes become more unstable during the GENE-Tango simulation, turbulent fluxes increase for all species and Tango reacts relaxing thermal profiles.
Additionally, we observe that while the global GENE simulations can tolerate marginal instability of these modes and still achieve power balance agreement, flux-tube simulations impose stricter constraints on the linear stability of these modes.As a consequence, flux-tube simulations are likely to undergo a more pronounced flattening of plasma pressure profiles compared to global simulations to reduce the drive of these modes until they reach a state of marginal stability.These findings suggest that a global gyrokinetic treatment might be necessary to accurately model the effects of energetic particles on turbulence in regimes where energetic particles and electromagnetic effects are expected to play a non-negligible role on plasma turbulence.
This paper is organized as follows.The ASDEX Upgrade discharge is described in section 2. While the GENE-Tango coupling is detailed in section 3, we outline the numerical setup employed for the GENE-Tango simulations in section 4. In section 5, we present the steady-state profiles obtained from TGLF-ASTRA, highlighting a significant under-prediction, particularly of the on-axis ion temperature.Building upon this, section 6 shows the plasma profiles resulting from various GENE-Tango simulations that retain or neglect energetic particles and electromagnetic effects in GENE.Notably, we demonstrate that the only case that successfully matches the experimental plasma profiles is when both fast particles and electromagnetic effects are retained in the GENE-Tango calculations.To shed light on the physical mechanisms driving the observed on-axis peaking of the temperature profile in the GENE-Tango simulations with fast ions and electromagnetic effects, we provide a detailed analysis in section 7. Additionally, stability analyses on the final steady-state GENE-Tango profiles that successfully reproduce the experimental results are conducted in section 8.The impact of strongly unstable high frequency electromagnetic modes on the thermal profiles is discussed in section 9. Section 10 presents a comparison between flux-tube and global GENE simulations for the steady-state GENE-Tango profiles.Finally, conclusions are drawn in section 11.

ASDEX Upgrade discharge #39230
The plasma discharge under investigation in this study corresponds to the ASDEX Upgrade pulse #39230.This is a deuterium plasma with an on-axis magnetic field strength of approximately B 0 = 2.5 T. The plasma current is I p ≈ 0.8 MA, and the safety factor at the flux surface containing 95% of the poloidal flux is q 95 ≈ 5.The heating power is comprised of P NBI ≈ 7.5 MW of Neutral Beam Injection (NBI) with perpendicular beams and P ECRH ≈ 1.6 MW of Electron Cyclotron Resonance Heating (ECRH).Figure 1 illustrates the time evolution of various parameters, including the heating powers, onaxis plasma temperatures and density.Notably, the application of NBI is closely associated with a significant rise in ion temperature.This effect is most prominent in the time range of t = [2.5-3.2]s, i.e. at the maximum NBI power in the discharge, where the ion temperature reaches its peak values at the center, with relatively smaller impacts on electron temperature and density.Consequently, our analysis focuses on a specific time slice at t = 2.7 s, during which the maximum external NBI power is employed, resulting in a pronounced effect on the on-axis ion temperature.The electron temperature and density measurements are obtained using IDA [27],  while the ion temperature and plasma toroidal rotation are determined through IDI [28].
Considering the substantial NBI external heating power applied in this plasma discharge, there exists a non-negligible fraction of supra-thermal ions.In figure 2, we show the temperature and density profiles of these energetic particle species, computed via RABBIT [29] using the experimental profiles.

Description of the GENE-Tango loop
The numerical simulations presented in this paper were performed using the newly developed tool GENE-Tango, which combines the gyrokinetic code GENE [30,31] with the transport solver Tango [32].This coupling enables reliable profile predictions with significant speed-ups compared to standalone flux-driven gyrokinetic simulations on the transport timescale [25].The underlying scheme follows the approach employed in reduced turbulence models and transport codes, where turbulence and transport phenomena are simulated only on their natural timescales [22].
The coupling starts with a global GENE simulation using an initial guess for the temperature and density profiles, running until saturation.Afterwards, the turbulent fluxes (averaged over the saturated phase) are extracted and transferred to Tango.Tango, in turn, evaluates new plasma (temperature and density for each species) profiles consistent with the given turbulence levels and experimental sources, evolving temperatures and densities for each species.These updated profiles are then fed back into GENE, restarting from the previous turbulent phase and running until saturation.Given the restart from the previous iteration, the system quickly converges to a new saturated phase within a short time domain (approximately t = 150c s /a for the simulations performed in this paper, where c s represents the sound speed defined in section 4 and a the minor radius).This iterative process continues until the turbulent fluxes match the volume-averaged energy and particle injected by external sources.Achieving a steady-state solution is followed by further verification by performing a stand-alone GENE simulation with the final plasma profiles computed by Tango.The coupled GENE-Tango simulation is considered fully converged only when this standalone GENE simulation yields the same fluxes as the last GENE-Tango iteration.The coupling scheme is illustrated in figure 3.
The use of GENE-Tango significantly reduces the computational expense by several orders of magnitude, making first-principle core plasma simulations over the confinement time feasible, even for devices such as ITER, using currently available computing resources [25].The validity of GENE-Tango has been demonstrated through validation studies on multiple ASDEX Upgrade discharges with moderate external heating powers, where the effects of energetic particles on plasma profiles were negligible [25].Recently, the application of GENE-Tango has been extended to optimized stellarator devices [33].However, the results presented in this paper are based only on radially global GENE-Tango analyses.A radially global treatment is indeed required to describe energetic particle dynamics and energetic-particle-driven modes correctly [34].
In the numerical simulations conducted for this study, the NBI heating power deposited on electrons and ions, as well as the fast deuterium pressure, were calculated using RABBIT [29].The ECRH power was computed using TORBEAM [35,36], while the ohmic power density was determined from ASTRA [7] and remained fixed throughout the GENE-Tango loop.The power density of the collisional energy exchange between ions and electrons was self-consistently computed by Tango and varies at each iteration.Additionally, the radiated power density, which takes into account bremsstrahlung radiation, was calculated using the experimental effective charge profile and the tungsten density.While the effective charge profile and the tungsten density were kept constant during the GENE-Tango loop, the radiated power density is consistently adjusted at each iterations accounting for variations in the electron temperature profile.

Numerical setup and grid resolution
The GENE simulations are performed employing the realistic ion-to-electron mass ratio.External flow shear and collisions are retained and the latter modeled using a linearized Landau-Boltzmann collision operator with energy and momentum conserving terms [37].Neoclassical fluxes are expected to be negligible for the analyzed discharge and are therefore not included in the calculations [38,39].For simplicity, each species has been modeled with a Maxwellian distribution function.Specifically, for the energetic particle species, this involves a Maxwellian distribution with the zeroth and second-order moments equal to the fast ion density and energy as computed by RABBIT.In the near future, the simplification of an equivalent Maxwellian distribution will be relaxed.While the GENE-Tango loop allows us to evolve the energetic particle profiles through Rabbit, we have kept the energetic particle profiles fixed at the ones computed for the experimental plasma profiles for each GENE-Tango simulations.This approach facilitates the comparison between the various simulations performed that cover electrostatic and electromagnetic regimes considering energetic particles.
The use of a more realistic distribution functions for modeling the energetic particle species is not expected to yield significant differences in the numerical results [40].While GENE-Tango has recently been coupled with CHEASE [41] to enable the self-consistent evolution of magnetic geometry as the plasma pressure profiles evolve within the Tango iterations, for simplicity, we have kept the magnetic equilibrium fixed to the one reconstructed via IDE [42] for the ASDEX Upgrade discharge #39230 at the time t = 2.7 s.
In this paper, we perform GENE-Tango simulations using only the radially global version of GENE.These simulations employ a grid resolution of (n x × n ky × n z ) = (225 × 36 × 32) in the radial (x), bi-normal (y), and field-aligned (z) directions.The velocity space grids consist of (n v × n µ ) = (48 × 32) points, with v /v th,s = [−3.5,3.5] and µB 0 /T s = [0, 12] for each species.Here, v denotes the velocity component parallel to the background magnetic field, and µ is the magnetic moment.These extended velocity ranges can capture for changes in pressure profiles within the radial box.The radial domain covered in the radially global GENE-Tango simulations is ρ tor = [0.025,0.7].Here, ρ tor represents the radial coordinate based on the toroidal flux Φ, calculated as ρ tor = √ Φ /Φ LCFS , where Φ LCFS represents the value of the toroidal flux at the last closed flux surface.
In GENE, v is normalized to the thermal velocity of each species s, v th,s = (2T s /m s ) 1/2 , where T s is the temperature and m s the mass; and the magnetic moment is normalized to B 0 , the on-axis magnetic field and T s .In the electrostatic simulations, we set the minimum toroidal mode number n 0,min = 2, while in the electromagnetic simulations, we used n 0,min = 1.The toroidal mode number n is discretized as n = n 0,min • j, where j ranges from 0 to n ky − 1.In the electromagnetic global simulations B fluctuations are neglected.
The global GENE simulations are performed in gradientdriven mode, using Krook-type particle and heat operators to maintain the plasma profiles close to the reference profiles provided by Tango during each iteration [31].To achieve this, the Krook heat coefficient (γ k ) and particle coefficient (γ p ) are set to γ k = 0.01c s /a and γ p = 0.01c s /a, respectively.Here, c s represents the sound speed, calculated as (T e /m i ) 1/2 , where T e is the electron temperature at the reference radial position, and m i is the bulk ion mass in proton units.In order to mitigate the possible impact of unresolved scales and instabilities such as electron temperature gradient driven modes, a numerical fourth-order hyperdiffusion method is employed to dampen fluctuations [43].
Additionally, buffer regions covering 10% of the GENE simulation's radial domain are implemented.These buffers serve the purpose of suppressing fluctuations to zero near the boundaries of the domain and ensuring compatibility with the Dirichlet boundary conditions.Within these regions, a Krook operator with an amplitude of γ b = 1.0c s /a is applied.
Tango employs different strategies to extend turbulent fluxes and sources in the GENE buffers, accounting for the undesired damping effects on plasma fluctuations.Within the inner buffer region, the physical sources are intentionally set to zero, with a minimal impact on the volume integral of the injected sources due to the reduced plasma volume near the magnetic axis.At the outer buffer, we employ extrapolation regions in Tango, with distinct radial domains chosen for temperatures and density.These regions are carefully selected in areas that remain unaffected by the Krook operator.Within these specified domains, Tango performs interpolation of the GENE turbulent fluxes.As we move beyond this region towards the outer boundary, the turbulent fluxes are replaced by a linear extrapolation until reaching the end of the GENE-Tango grid.We refer the reader to references [23,25,32] for further details.
To enhance numerical stability, Tango incorporates relaxation coefficients for plasma pressure (α p ) and turbulent fluxes (α q ), facilitating effective averaging over previous iterations and mitigating statistical variance arising from turbulence simulation [23,25,32].In this study we used α p = 0.1 and α q = 0.3.Furthermore, Neumann boundary conditions are applied for the inner boundary, while Dirichlet boundary conditions are used for the outer boundary in Tango.

Underestimation of temperature profiles with TGLF-ASTRA
We begin our analysis of the ASDEX Upgrade discharge #39230 at t = 2.7 s using TGLF-ASTRA covering the same GENE-Tango radial domain (ρ tor = [0.025− 0.7]).In this analysis, we employed the SAT2 saturation rule for TGLF [44].The steady-state profiles obtained from this analysis are shown in figure 4 and compared with the experimental profiles.While we observe a good agreement for the electron density, we note a significant underestimation of the temperature profiles, especially concerning the ion temperature predicted by TGLF-ASTRA, in comparison to the experimental measurements.This finding is consistent with previous observations indicating that TGLF-ASTRA generally fails to capture the positive effects of supra-thermal particles on turbulent transport often linked to nonlinear electromagnetic effects, resulting in larger turbulent fluxes and, thus flatter profiles [11][12][13][14].It is worth mentioning that recently, wave-particle resonant effects between energetic particles and ion-scale micro-instabilities have been successfully reproduced with TGLF simulations [15].As fast ions primarily impact the thermal ion turbulent fluxes, this discrepancy leads to a noticeable decrease in the on-axis T i , with a minor effect on the electron temperature while still maintaining qualitative agreement with the density profile.

Global GENE-Tango simulations
We apply the global version of GENE interfaced with Tango to simulate this ASDEX Upgrade discharge, employing the numerical setup described in section 4. The plasma profiles were initialized to the ones reconstructed via IDA for electron temperature and density, and IDI for toroidal rotation and ion temperature.While the toroidal rotation remained fixed during the GENE-Tango iterations, the other profiles were selfconsistently adjusted due to the combined impact of heating   and particle sources.When included in the modeling, the temperature and density profiles of energetic particles were computed using RABBIT.Additionally, the thermal ion density is kept equal to the electron density profile to ensure quasineutrality in the simulations without energetic particles.
Four simulations were performed: (i) electrostatic without fast particles, (ii) electromagnetic without fast particles, (iii) electrostatic with fast particles, and (iv) electromagnetic with fast particles.In each case, GENE-Tango was run until the global turbulent fluxes computed by GENE matched the volume integral of the injected sources at each location, e.g. as shown in figure 5 for the electromagnetic simulation with fast particles.The final steady-state profiles are illustrated in figure 6 for each of the different cases and compared with the experimental profiles.
Figure 6 shows a particularly good agreement between GENE-Tango and experimental profiles for the electron temperature and plasma density.However, when focusing on ion temperature, it is found that achieving an agreement with experimental data requires the inclusion of both fast particles and electromagnetic effects in the modeling.This suggests that the joint influence of these factors is crucial in accurately reproducing the observed T i profiles.It further highlights the importance of considering the complex interplay between fast particles and electromagnetic effects for a realistic description of plasma behavior in plasma regimes with high external heating power and energetic particle concentration.
The analysis also reveals that the magnetic equilibrium alone, characterized by features such as reversed shear and rational surfaces, is insufficient in explaining the T i peaking, since it was kept fixed in all of the different simulations (figure 7).
Electrostatic simulations, regardless of the inclusion of fast ions, predict excessive turbulent transport.As a consequence, the central T i remains at a relatively constant value of T i,0 = 4 keV, with minimal effects on T e and plasma density.On the other hand, the inclusion of electromagnetic effects without fast ions leads to a mild peaking of T i .Although there is some improvement in reproducing the experimental profiles, the influence on the on-axis T i value remains limited.

Physics behind the T i peaking
The central role played by electromagnetic effects and fast particles on the ion temperature profile peaking becomes evident when looking at the radial profile of the logarithmic gradient of the thermal ion temperature and the energetic particle heat flux.This is shown in figure 7, together with the results obtained from the electromagnetic GENE-Tango simulations without fast particles (as a reference).
Figure 7 illustrates a clear correlation between the regions experiencing a significant rise in the logarithmic gradient of the thermal ion temperature and the locations where the energetic particle heat flux increases.This occurs at the rational surface q = 1 and at s = 0.These locations correspond to the excitation of high frequency KBMs/AEs (as discussed in section 8) [45].This indicates that the presence of these modes is likely to be linked to the turbulence suppression of the thermal ion turbulent fluxes, resulting in the peaking of the thermal ion temperature profiles.
The presence of these modes is particularly important in the nonlinear GENE simulations to strongly suppress ionscale turbulent transport and increase the logarithmic gradient of the thermal ion temperature profiles.This effect is illustrated in figure 8 where the time evolution of the radial profile of the flux-surface averaged radial electric field is shown for the electromagnetic simulation with supra-thermal particles.Interestingly, we observe the generation of zonal layers in the proximity of the radial location where the energetic particle heat fluxes are larger in the nonlinear GENE simulations.However, no such behavior is observed in the zonal currents (not shown here), as they do not exhibit a noticeable increase at the rational surface q = 1 or at the location s = 0. Additionally, figure 8 compares the time-averaged (over the saturated nonlinear phase) radial electric field computed in electromagnetic simulations with and without energetic particles during the saturated nonlinear phase.It is evident that the simulation without fast particles lacks the layer at s = 0 of the radial electric field, which is driven by the KBM/AE that is absent in the simulation without fast particles.
Additionally, it is worth noting that while the electromagnetic simulation neglecting fast particles does indicate an increase in the radial electric field around q = 1, this effect is confined to a localized region.In contrast, in the simulation incorporating fast particles, the radial electric field exhibits a similar magnitude at q = 1, but with the notable difference of a more extensive zonal structure that spans a broader radial range from ρ tor = [0.2-0.3].This extended zonal structure may account for the larger amplification of the logarithmic temperature gradient in the thermal ions near the rational surface q = 1 when fast ions are present, as compared to the case without fast ions.
These results are consistent with previous findings observed in flux-tube simulations, showing that fast particles play a significant role in destabilizing marginally stable modes through nonlinear interactions [18,46].In particular, as proposed by [47], the mechanism leading to the nonlinear excitation is likely the inverse Compton scattering by energetic ions.These modes are shown to lead to an energy redistribution from the thermal modes to the fast ion-driven modes, initially depleting the energy content of ion-scale turbulence.When the drive from fast ions is sufficiently large, it was observed that these fast particle modes interact with zonal flows, leading to a substantial increase in shearing rate levels, ultimately resulting in a strong suppression of heat and particle fluxes.The impact of energetic particle modes in enhancing zonal flow and zonal currents has also been observed in global gyrokinetic [48,49] simulations and predicted by analytic theory [47,50,51].
Finally, figure 7 also shows that the other rational surfaces do not significantly impact the fast ion heat flux.This can be attributed to the relatively low fast ion density and temperature at ρ tor > 0.3, which strongly reduce the drive of these modes.We observe that the simulations without fast ions-despite having same geometry-do not show clear improvements of the logarithmic gradient of the thermal ion temperature profile at s = 0 and q = 1.

Linear and nonlinear stability analyses
Given the high relevance of the electromagnetic modes destabilized by fast particles on the T i temperature peaking observed in GENE-Tango, it is crucial to characterize their linear and nonlinear stability.To this aim we analyze the electrostatic and electromagnetic (A ) potentials obtained from the nonlinear stand-alone GENE simulations retaining fast particles and electromagnetic fluctuations.The results, shown in figure 9, are obtained using the GENE-Tango steady-state plasma profiles.Figure 9 reveals that the mode with toroidal mode numbers n = 7 at ρ tor = 0.13-which corresponds to the location of s = 0-exhibits the largest amplitude.Additionally, significant contributions are also observed from modes with n = [1, 2] at ρ tor = 0.23, located at the rational surface q = 1.for the electrons at the rational surface q = 1, obtained in the electromagnetic GENE stand-alone simulation with energetic particles using the steady-state GENE-Tango profiles.
These modes have a considerable impact on the turbulent heat fluxes at these specific locations for each plasma species.This is shown in figure 10, which displays the ion heat flux spectra (electron heat flux is negligible at s = 0) at the location of s = 0 (left) and the electron spectra at q = 1 (right).In particular, the turbulent fluxes at s = 0 are dominated by the mode at n = 7 and its harmonics.Notably, figure 10 shows a particularly large electromagnetic contribution (flatter) to the turbulent fluxes at these specific locations of interest for the modes n = 1 and n = 7, thereby highlighting the fundamental electromagnetic nature of these modes.In order to provide a more comprehensive characterization of these modes, we present in figure 11 the frequency spectra of the electrostatic potential for different toroidal mode numbers obtained from the nonlinear electromagnetic radially global GENE simulation using the steadystate GENE-Tango profiles taken at the location of s = 0 (left) and q = 1 (right).Figure 11 reveals that at the location s = 0, the electrostatic potential spectrum is predominantly influenced by high-frequency modes at n = [7,14], with only a minor contribution from the low-frequency ITG branch.As we move to the rational surface at q = 1, the significance of the low-frequency branch increases.However, it is important to note that a non-negligible contribution from a high-frequency mode with ω ≈ 150 kHz, occurring within the range of toroidal mode numbers n = [1-7], is still observed.also showcases slices of the electrostatic potential spectra for the toroidal mode numbers n = 7 and n = 1, respectively, at the locations of s = 0 and q = 1, revealing more clearly that the electromagnetic mode observed is characterized by a highfrequency component of approximately ω ∼ [150-200] kHz.
To provide further insights into the nature of this highfrequency mode, we perform radially global linear GENE simulations using the final steady-state GENE-Tango profiles corresponding to the case retaining electromagnetic effects and fast particles.The linear results are shown in figure 12, respectively for the growth rate and frequency varying the plasma beta.This has been done by artificially re-scaling the plasma beta of each species by constant value, keeping the same shape as the nominal beta profile (that measures β e = 0.013 at mid radius for the steady-state profiles).While the mode at n = 7 is found to be linearly unstable for the steady-state GENE-Tango profiles and gets marginally stable at β e = 0.008, the modes at n = [1,2] observed in the nonlinear spectra of the electrostatic potential (see figure 9) are found to be linearly stable for the steady-state profile since at these toroidal mode numbers we only see an electrostatic mode at low frequency.This finding suggest that the electromagnetic high-frequency mode with n = 1 is linearly stable for the final GENE-Tango profile and possibly destabilized nonlinearly.
The poloidal cross section of the electrostatic potential for the toroidal mode number n = 7 at β e = 0.013 is shown in figure 13.We can observe that the two dominant poloidal mode numbers are m = 6 and m = 7.Additionally, it is worth mentioning that the mode is well localized at the location where s = 0 also in the linear simulations, which is consistent with the nonlinear results presented previously.Figure 13 displays the radial frequency spectrogram of the linear simulations performed for n = 7 at the nominal plasma beta β e = 0.013, accompanied by the Alfvén continuum computed using the  expression from [52].Notably, figure 13 highlights that the high-frequency mode observed in the GENE linear simulations is destabilized between two distinct Toroidal AE gaps, with a frequency exceeding that of a reversed shear AE [45,53,54].
Interestingly, this high frequency electromagnetic mode can be destabilized even in the absence of energetic particles.This is shown in figure 14, where linear growth rates and frequencies are compared for various toroidal mode numbers in simulations both with and without energetic particles.The thermal profiles were initialized to those obtained from the GENE-Tango simulation retaining simultaneously energetic particles and electromagnetic effects.Specifically, figure 14 indicates that the high-frequency mode at n = [6-7] persists even in the absence of fast ions, provided the plasma beta reaches sufficiently large values, akin to what was observed in the GENE-Tango simulations with fast ions and electromagnetic effects.
Nonetheless, there are notable distinctions in the high frequency mode behavior with and without the presence of energetic particles.This is shown in figure 15 when looking at growth rate and frequency dependence on the total plasma beta (obtained as the sum of thermal ion, electron and fast ion contributions).Specifically, we observe that the linear growth rates are reduced, and the mode frequency is higher in the absence of fast ions.Moreover, figure 15 shows that the threshold for destabilizing this high frequency mode shifts to higher beta values when fast ions are incorporated into the model.The substantial increase in growth rates for both the high frequency mode and ITG, observed at the same total plasma beta in the absence of energetic particles, bears significant consequences for turbulent fluxes.Specifically, when conducting two stand-alone GENE electromagnetic global simulations with and without energetic particles initializing the same profiles as the GENE-Tango electromagnetic simulation with energetic particles and electromagnetic effects, we notice large differences in the turbulence levels.The simulation without energetic particles displays notably enhanced turbulent fluxes in each channel, as shown in figure 16.These findings imply that in the absence of energetic particles, both ITGs and high frequency modes are more unstable compared to scenarios involving energetic particles, resulting in overall larger turbulent fluxes.This aligns with the analyses presented in the previous sections.
It is essential to note that the magnetic coil spectrogram analyses for this ASDEX Upgrade discharge do not show signatures of these high-frequency modes.This absence of the mode may be attributed to the weak linear drive observed in the GENE simulations for the steady-state plasma profiles.Specifically, the mode with n = 7 at the location s = 0 exhibits a very small linear growth rate, only slightly exceeding the ITG growth rate, as seen in figure 12.In addition, the other mode observed in figure 9 with n = 1 at the rational surface q = 1 remains linearly stable and is possibly destabilized nonlinearly.As a result, also this mode is characterized by a notably weak drive.Therefore, it is possible that the experimental measurements of these weak signals is particularly challenging.Future investigations will focus on designing similar plasma discharges at ASDEX Upgrade to evaluate the existence of these high-frequency modes through density fluctuation measurements.Unfortunately, such measurements were not available in the plasma core for the discharge analyzed in this paper.

Further destabilization of the KBM/AEs at q = 1
During the electromagnetic GENE-Tango convergence loop performed with energetic particles, in a few iterations Tango increased the plasma pressure at the point to destabilize modes at the rational surface q = 1 above their linear stability threshold.When this happened in the GENE-Tango loop (i.e. when the high frequency mode became strongly unstable), we observed a significant increase in the turbulent fluxes of the energetic particle species.As a result, there is an increase in turbulent fluxes for each thermal species.Figure 17 illustrates this, comparing the turbulent fluxes of energetic particles and thermal ions in two iterations: (i) when the mode at the rational surface q = 1 is marginally stable, and (ii) when Tango increases the pressure beyond the linear stability threshold,  leading to additional destabilization.These findings demonstrate the negative impact of these modes on plasma profiles and overall performance when left uncontrolled consistently with previous studies [46,48,49,55].
Figure 17 shows the impact of the increased fluxes on the thermal profiles, revealing how Tango responded by reducing the pressure gradients.This adjustment was intended to regulate and mitigate the destabilizing effects induced by energetic particles at q = 1 by reducing their drive.These observations indicate that the simulation aims to achieve a critical equilibrium where the modes at q = 1 are maintained marginally stable or weakly unstable, as these modes corresponds to an increase in the thermal turbulent fluxes.In other words, the profile stiffness in these regimes and regions is determined by the energetic particles if the total plasma beta and its radial gradient reach certain thresholds.

Comparison with flux-tube
In this section, we compare the GENE-Tango turbulent fluxes with the ones obtained running the flux-tube version of GENE at different locations using the steady-state profiles for the case retaining both fast ions and electromagnetic effects.Specifically, we seek to determine if the flux-tube simulations can appropriately represent the energetic particle stabilizing effects, leading to a matching of the volume integral of the injected sources at different locations, similar to what is observed in the global GENE-Tango simulations.This comparison is particularly relevant when considering computational cost of the simulations, especially when looking at future fusion reactors, where the value of 1/ρ * becomes very large.Under these conditions, local gyrokinetic simulations require significantly lower resolution to accurately represent the relevant bi-normal and radial scales together with the velocity space.However, it is important to note that fusion reactors will have fundamentally very different fast particles (alpha particles) compared to those with moderate energies generated through heating schemes.Therefore, the findings of this section should be applied specifically to fast particle effects of present devices.
The local results are obtained by running GENE fluxtube simulations at seven different radial locations, covering a range of ρ tor = [0.1-0.7].We used a grid resolution of (n kx × n ky × n z ) = (256 × 48 × 32) in the radial (x), bi-normal (y), and field-aligned (z) directions.The velocity space grids consist of (n v × n w ) = (32 × 16) points, covering the velocity range v /v th,s = [−3.0,3.0] and µB 0 /T s = [0, 9] for each species (s).Periodic boundary conditions are employed along the radial direction allowing us a Fourier decomposition along the radial direction as well.The plasma geometry and energetic particle profiles are kept the same as the ones used for the global GENE-Tango simulation with energetic particles and electromagnetic effects.
In figure 18, we present a comparison between the turbulent fluxes obtained from the flux-tube version of GENE and those from the global version (considering the last five GENE-Tango iterations), alongside the volume integral of the experimental sources.Strikingly, we find a significant discrepancy across all channels, with the flux-tube simulations greatly over-predicting turbulent fluxes for ρ tor > 0.2.Specifically, at ρ tor = 0.3, where the linearly unstable KBM/AE modes are present, there is a pronounced increase in the turbulent fluxes, particularly for the electrons (both heat and particles), and to a lesser extent, also for the thermal ions.
Interestingly, the mismatch between the global and fluxtube fluxes is not limited to the locations with high frequency modes but also extends to areas without these modes.This observation aligns with other analyses demonstrating that energetic particles can influence turbulent fluxes even in radial locations where locally fast ions play a less prominent role [34].Reference [34] shows differences on the turbulent fluxes qualitatively similar to the ones observed in this section in presence of unstable high frequency modes.In cases where a KBM/AE mode is destabilized in neighboring locations, global modifications of the underlying electric field can impact turbulent fluxes.Such global effects are not adequately captured in flux-tube models, where fluxes at each radial location are treated independently of neighboring positions.
It is also interesting to perform a linear stability analysis for each flux-tube selected for the nonlinear comparison.The resulting linear growth rates and frequencies are displayed in figure 19.As shown in this figure, we observe at ρ tor = 0.3 a linearly unstable high frequency with frequencies similar to those observed in the linear and nonlinear global GENE simulations.Surprisingly, another linearly unstable mode is found at the location ρ tor = 0.4 not present in the global simulations, but with a notably lower linear drive compared to the former case.This finding suggests that capturing the radial changes in plasma pressure and safety factor profiles is essential to reproduce the correct drive/damping mechanisms for these high frequency modes.Figure 19 reveals that the flux-tube approximation has a strong impact also on the linear properties of these high frequency modes.One notable limitation-among others-is the inability of local modes to correctly describe the interaction with the shear Alfvén wave continuum, which significantly impact the growth rate of the mode, thus strongly affecting in turn also the nonlinear dynamics that govern the interaction between high frequency modes, ITG turbulence, and zonal flows [45].
Notably, our results suggest that global GENE simulations exhibit greater tolerance for marginal instability of these modes, achieving power balance agreement even in the presence of such instabilities.On the other hand, flux-tube simulations impose stricter constraints on the linear stability of these modes.The flux-tube approach, in particular, shows a notable increase in overall injected heat and particle fluxes when a linearly unstable high frequency electromagnetic mode occurs, far exceeding the levels from external sources.These observations strongly advocate for the necessity of a global gyrokinetic treatment to accurately model the effects of energetic particles on turbulence.The limitations of flux-tube simulations hinder their ability to fully capture the intricate interplay between energetic particles and turbulence.Sometimes they yield qualitatively correct results [56].

Conclusions
This paper presents the first gyrokinetic simulations of an ASDEX Upgrade discharge up to the transport time scale retaining energetic particles and electromagnetic effects.Specifically, the study focuses on the H-mode pulse #39230, which exhibits a pronounced peaking of the on-axis ion temperature.The numerical simulations are performed using the newly developed tool GENE-Tango, which integrates the radially global version of the gyrokinetic code GENE with the transport solver Tango.This coupling ensures reliable profile predictions-based on global gyrokinetic calculations-while significantly reducing computation time in comparison to standalone flux-driven simulations.
To investigate the influence of energetic particles on the observed improvement in plasma performance, we conducted four distinct GENE-Tango simulations: (i) electrostatic simulation without fast particles, (ii) electromagnetic simulation without fast particles, (iii) electrostatic simulation with fast particles, and (iv) electromagnetic simulation with fast particles.
We found that the inclusion of energetic particle species in the radially global nonlinear electromagnetic GENE-Tango simulations led to an excellent agreement with experimental plasma profiles.On the contrary, when fast ions and electromagnetic effects were simultaneously neglected in the GENE modeling we observed a significant flattening of the ion temperature profile, consistent with profiles computed using TGLF-ASTRA, which notably under-predicts the on-axis ion temperature.
Thorough analyses conducted in our study have revealed the critical role of high-frequency KBMs/AEs (ω ∼ [150-200] kHz) in suppressing ion-scale plasma turbulence.These modes play a critical role in the global GENE simulations enhancing zonal flow activity and increasing the shearing rate levels.These findings are consistent with recent experimental observations at DIII-D [57].
Interestingly, as the GENE-Tango simulation advanced, there were a few iterations when Tango raised the plasma pressure beyond a critical threshold leading to a destabilization of these high frequency modes.However, this led to a strong increase of the turbulent fluxes of each species within the GENE simulations.As a result, Tango responded by flattening the thermal profiles to achieve a match of the turbulence fluxes with the experimental power balance leadingpotentially-to a large confinement degradation.These findings suggest that the simulation seeks a delicate balance by keeping these modes marginally stable or only weakly unstable.
It is worth emphasizing that the magnetic coil spectrogram analyses of this ASDEX Upgrade discharge did not reveal any indications of the high-frequency modes observed in the GENE simulations.However, it is possible that the particularly weak linear drive of these quite central modes may impose a considerable challenge to detect their signatures through experimental measurements.Future investigations will focus on designing similar plasma discharges at ASDEX Upgrade to evaluate the existence of these high-frequency modes through density fluctuation measurements, not available in the plasma core for the discharge analyzed in this paper.
The mechanism responsible for the observed turbulence suppression bears qualitative similarity to that previously identified in flux-tube simulations where marginally stable energetic particle-driven modes were nonlinearly destabilized, leading to a depletion of the free energy of the ITG turbulence and a subsequent increase in the levels of zonal flows [18,46].However, we noted fundamental differences when comparing the global GENE-Tango simulations with results obtained from stand-alone (without Tango) simulations using the flux-tube code at various radial locations.Specifically, global GENE simulations exhibited a higher tolerance for marginal instability of these modes while maintaining power balance agreement.In contrast, flux-tube simulations imposed stricter constraints on the linear stability of these modes, limiting their ability to capture the full complexity of the interplay between energetic particles and turbulence.As a consequence, the plasma pressure profiles in flux-tube simulations would lead to a more pronounced flattening if coupled to a transport code compared to the one observed when running global gyrokinetic simulations to reduce the drive of these modes until they reached a state of marginal stability.These findings underscore the necessity of a global gyrokinetic treatment to accurately model the effects of energetic particles on turbulence, especially in scenarios where energetic particle modes are linearly unstable [56].

Figure 1 .
Figure 1.Time evolution of the (a) NBI and ECRH heating power; (b) experimental on-axis plasma density, electron and ion temperature; and (c) and ion temperature logarithmic gradients at two different radial locations for the ASDEX Upgrade discharge #39230.While the electron temperature and density are obtained using IDA, the ion temperature with IDI.

Figure 2 .
Figure 2. Radial profile of the energetic particle (a) temperature and (b) density profiles computed by RABBIT using the experimental plasma profiles for the thermal species for the ASDEX Upgrade discharge #39230 at t = 2.7 s.

Figure 3 .
Figure 3. Illustration of the GENE-Tango loop, encompassing CHEASE and RABBIT, for the self-consistent evolution of magnetic equilibrium, energetic particle temperature and density profiles, and heating profiles.Tango computes changes in the temperature and density profiles, while the gyrokinetic code GENE calculates the turbulent fluxes Qm corresponding to each pressure profile pm.The fluxes are then transferred back to Tango, which evolves the temperature and density profiles to the next macroscopic time step m + 1.This iterative process continues until the turbulent fluxes computed by GENE match the volume integral of the heat and particle sources at each radial location.

Figure 4 .
Figure 4. Comparison of the (a) ion, (b) electron temperatures and (c) density computed by TGLF-ASTRA (blue) and experimental measurements (black).The temperature and density profiles are obtained using Bayesian analyses.

Figure 5 .
Figure 5. Time-averaged radial profile of the (a) ion, (b) electron heat fluxes in MW and (c) particle flux in 1/s corresponding to the averaged last five GENE-Tango iterations (red) in the simulations considering energetic particles and electromgnetic effects.The shaded gray areas denote the buffer regions and the black circles the volume integral of the injected particle and heat sources.

Figure 6 .
Figure 6.Comparison of the (a) ion, (b) electron temperatures and (c) density computed by GENE-Tango for the four different cases considered: with both fast ions and electromagnetic effects (red), with fast ions but without electromagnetic effects (cyan), without fast ions but with electromagnetic effects (blue), without both fast ions and electromagnetic effects (magenta).Experimental measurements are represented by the black curve.

Figure 7 .
Figure 7. Radial profile of the time-averaged (a) energetic particle turbulent flux corresponding to the averaged last five GENE-Tango iterations (electromagnetic case), (b) safety factor profile reconstructed via IDE, (c) logarithmic thermal ion temperature gradient obtained in the electromagnetic GENE-Tango simulation with (red) and without (blue) energetic particles, and (d) ratio between the logarithmic thermal ion temperature gradient obtained in the electromagnetic simulations with and without energetic particles.The continuous vertical black lines denote the location of the rational surfaces q = [1, 3/2, 2], while the dotted black line represents the location of null shear s = 0.The logarithmic gradient is defined as ω Ti = −adlnT i /dρtor with a being the minor radius.

Figure 8 .
Figure 8.(a) Time evolution of the radial profile of the flux-surface averaged radial electric field obtained in the electromagnetic GENE stand-alone simulation with energetic particles using the steady-state GENE-Tango profiles.(b) Comparison of the time-averaged radial electric field obtained in the electromagnetic GENE-Tango simulation with (red) and without (blue) energetic particles.The continuous vertical black lines denote the location of the rational surfaces q = [1, 3/2, 2], while the dotted black line represents the location of zero shear s = 0.

Figure 9 .
Figure 9. Radial profile of the time-averaged (over the nonlinear saturated phase) of the (a) electrostatic and (b) electromagnetic (A ) potentials for different toroidal mode numbers obtained in the electromagnetic GENE stand-alone simulation with energetic particles using the steady-state GENE-Tango profiles.

Figure 10 .
Figure10.Nonlinear heat flux spectra (a) for the thermal ions at the zero-shear (s = 0) location and (b) for the electrons at the rational surface q = 1, obtained in the electromagnetic GENE stand-alone simulation with energetic particles using the steady-state GENE-Tango profiles.

Figure 11 .
Figure 11.Fourier spectra at different toroidal mode numbers of the electrostatic potential (a) at the zero-shear location and (b) at the rational surface q = 1, obtained in the electromagnetic GENE stand-alone simulation with energetic particles using the steady-state GENE-Tango profiles.Slices at the toroidal mode numbers n = 7 and n = 1 are shown in (c) and (d) respectively at the zero-shear location and at the rational surface q = 1.The electrostatic potential has been averaged over the field-aligned coordinate z over the nonlinear saturated phase.The amplitude of the signal is plotted on a logarithmic scale.

Figure 12 .
Figure 12.Contour plots of growth rates (a) and frequencies (b) obtained by radially global linear GENE electromagnetic simulations at different values of the toroidal mode number and βe (at the center of the GENE radial box).The nominal value of the plasma beta for the steady-state GENE-Tango profile is marked by the horizontal black line.Slices of growth rates (c) and frequencies (d) at different βe.The value of βe corresponding to the steady-state GENE-Tango profile is βe = 1.3% and depicted by the blue line in (c) and (d).

Figure 13 .
Figure 13.(a) Poloidal cross section of the electrostatic potential at n = 7 for the GENE global simulation using the steady-state GENE-Tango profile.The dashed black line delimits the location of s = 0; (b) radial frequency spectrogram for the electrostatic potential at n = 7 obtained from the GENE global simulation using the steady-state GENE-Tango profile.The solid black lines represent the Alfvén continuum for the poloidal mode numbers m = [6, 7].

Figure 14 .
Figure 14.Comparison of the linear growth rates (a) and frequencies (b) obtained by radially global linear GENE electromagnetic simulations at different values of the toroidal mode number for simulations retaining (blue) or neglecting (red) the energetic particle species.The thermal profiles were initialized to those obtained from the GENE-Tango simulation retaining simultaneously energetic particles and electromagnetic effects.

Figure 15 .
Figure 15.Comparison of the linear growth rates (a) and frequencies (b) obtained by radially global linear GENE electromagnetic simulations at different values of the total plasma beta (obtained as the sum of thermal ion, electron and fast ion contributions) at the toroidal mode number n = 7 for simulations retaining (blue) or neglecting (red) the energetic particle species.The thermal profiles were initialized to those obtained from the GENE-Tango simulation retaining simultaneously energetic particles and electromagnetic effects.

Figure 16 .
Figure 16.Comparison of the time-averaged radial profile of the (a) ion, (b) electron heat fluxes in MW and (c) particle flux in 1/s for simulations retaining (blue) or neglecting (red) the energetic particle species.The thermal profiles were initialized to those obtained from the GENE-Tango simulation retaining simultaneously energetic particles and electromagnetic effects.

Figure 17 .
Figure 17.Comparison of the (a) ion temperature profile, (b) time-averaged energetic particle and thermal ion (c) turbulent fluxes computed by GENE-Tango during an iteration with weakly (red) and strongly (blue) unstable high frequency modes.The continuous vertical black lines denote the location of the rational surfaces q = [1, 3/2, 2], while the dotted black line represents the location of null shear s = 0. Experimental measurements are depicted by the black curve for comparison with the simulated profiles.

Figure 18 .
Figure 18.Time-averaged radial profile of the (a) ion, (b) electron heat fluxes in MW, (c) particle flux in 1/s and (d) fast ion heat flux in MW corresponding to the averaged last five GENE-Tango iterations (red) and the flux-tube simulations at seven different radial locations (blue).The shaded gray areas denote the buffer regions and the black circles the volume integral of the injected particle and heat sources.

Figure 19 .
Figure 19.Contour plots of growth rates (a) and frequencies (b) obtained by flux-tube simulations at different radial locations for different values of the toroidal mode number using the steady-state GENE-Tango profiles.