First principles of modelling the stabilization of microturbulence by fast ions

PAPER • OPEN ACCESS First principles of modelling the stabilization of microturbulence by fast ions To cite this article: G.J. Wilkie et al 2018 Nucl. Fusion 58 082024 View the article online for updates and enhancements. Related content Turbulence in high-beta ASDEX upgrade advanced scenarios H. Doerk, A. Bock, A. Di Siena et al. - Global anomalous transport of ICRH- and NBI-heated fast ions G J Wilkie, I Pusztai, I Abel et al. - Ion-scale turbulence in MAST: anomalous transport, subcritical transitions, and comparison to BES measurements F van Wyk, E G Highcock, A R Field et al. - Recent citations Predictive multi-channel flux-driven modelling to optimise ICRH tungsten control and fusion performance in JET F.J. Casson et al - Impact of fast ions on a trapped-electron- mode dominated plasma in a JT-60U hybrid scenario S. Mazzi et al - Impact of fast ions on density peaking in JET: fluid and gyrokinetic modeling F Eriksson et al - This content was downloaded from IP address 87.218.223.151 on 15/07/2020 at 16:30 1 © EURATOM 2018 Printed in the UK 1. Fast ion stabilization of ITG microturbulence Microturbulence fundamentally limits the confinement time in current and future tokamak experiments [2]. The very gradients that are required to achieve high central densities and temperatures also provide a source of free energy. This free energy drives transport of particles, momentum, and energy, usually far in excess of collisional transport. The ion temper ature gradient (ITG) mode, for instance, limits the core temper ature of tokamaks [3–5]. This turbulence occurs on the scale of the thermal ion Larmor radius ρi . Gyrokinetics is the reduction of the Fokker–Planck kinetic equation that rigorously handles electromagnetic fields whose fluctuations vary on spatial scales similar to ρi , but on timescales much slower than the gyro-frequency Ωi [6–8]. A multitude of computational tools have been developed to solve the nonlinear gyrokinetic equation[9–12]. These tools take as inputs the equilibrium magnetic field, density, and temperature profiles, and predict the ensuing turbulent fluctuations and the associated transport fluxes. By using experimentally-determined profiles and comparing the computed fluxes to experimentally-inferred fluxes (usually from powerbalance calculations), and through more detailed comparisons of turbulence characteristics, the validity of the gyrokinetic approach can be verified [13]. This has been widely demonstrated [14–17]. However, due to the scarcity of computational resources, these matching exercises of necessity entail the use of simulations that neglect various parts of the complex physics of the experimental setup. It was discovered during one of these Nuclear Fusion First principles of modelling the stabilization of microturbulence by fast ions G.J.Wilkie1, A.Iantchenko1,2, I.G.Abel1, E.Highcock1, I.Pusztai1 and JET Contributors1,a 1 Department of Physics, Chalmers University of Technology, Gothenburg, SE-41296, Sweden 2 Swiss Plasma Center, Ecole Polytechnique Federale de Lausanne, Lausanne Switzerland E-mail: [email protected] Received 21 December 2017, revised 22 February 2018 Accepted for publication 16 March 2018 Published 29 June 2018 Abstract The observation that fast ions stabilize ion-temperature-gradient-driven microturbulence has profound implications for future fusion reactors. It is also important in optimizing the performance of present-day devices. In this work, we examine in detail the phenomenology of fast ion stabilization and present a reduced model which describes this effect. This model is derived from the high-energy limit of the gyrokinetic equationand extends the existing ‘dilution’ model to account for nontrivial fast ion kinetics. Our model provides a physicallytransparent explanation for the observed stabilization and makes several key qualitative predictions. Firstly, that different classes of fast ions, depending on their radial density or temperature variation, have different stabilizing properties. Secondly, that zonal flows are an important ingredient in this effect precisely because the fast ion zonal response is negligible. Finally, that in the limit of highly-energetic fast ions, their response approaches that of the ‘dilution’ model; in particular, alpha particles are expected to have little, if any, stabilizing effect on plasma turbulence. We support these conclusions through detailed linear and nonlinear gyrokinetic simulations. Keywords: turbulence, fast ions, stabilization, microturbulence, gyrokinetics, simulation (Some figuresmay appear in colour only in the online journal) G.J. Wilkie etal First principles of modelling the stabilization of microturbulence by fast ions Printed in the UK 082024 NUFUAU © EURATOM 2018 58 Nucl. Fusion NF 10.1088/1741-4326/aab727 Paper 8 Nuclear Fusion IOP Original content from this work may be used under the terms of the Creative Commons Attribution 3.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. International Atomic Energy Agency a Litaudon etal [1] 2018 1741-4326 1741-4326/18/082024+15$33.00 https://doi.org/10.1088/1741-4326/aab727 Nucl. Fusion 58 (2018) 082024 (15pp) G.J. Wilkie etal 2 exercises, in which experiments on the Joint European Torus [18] were analysed, that effects from non-thermal minority ions had to be accounted for. Indeed, the predictions of the gyrokinetic codes in the absence of such ions suggested far greater transport than was observed [19]. This effect—the fast ion stabilization of ITG turbulence—is the subject of this paper. This effect is important because the presence of energetic ions are essential to sustain fusion-relevant bulk temperatures. Two external heating methods are used to supply the plasma with high-energy ions: neutral beam injection (NBI) [20] and ion cyclotron resonance heating (ICRH) [21]. As well as the fast ions resulting from external heating, another class of fast ions, energetic alpha particles, are generated from the fusion reaction itself. In the future, it is anticipated that alpha particles will provide the majority of the heating, and thus, that the plasma will be energetically self-sustaining (i.e. a burning plasma). In the local gyrokinetic simulations used herein, the practical differences between these several classes of fast ions (as well as, of course, their mass and charge) are in their distribution in phase space. As the main medium of heat injection, fast ions are a critical component of the fusion plasma and thus worthy of study in their own right. However, above and beyond this, it has been observed that a large fast ion population can improve the plasma confinement [22], leading to elevated core temperatures and/or densities. It has been shown in gyrokinetic simulations [19] that this effect is non-trivial. It is important to understand how and why this stabilization occurs, par ticularly when extrapolating to future net-power-producing fusion reactors, where alpha particles provide most of the heating. Although, like externally driven fast ions, alpha particles will contribute significantly to the plasma pressure in reactors, their density is very small ( nα0.01ne ) compared to NBI and ICRH ions. It is therefore not obvious that they will affect turbulence in a similar way. Owing to the potentially profound implications of this and other effects of fast ions on microturbulence, there has been widespread study of the topic. Follow-up works [23, 24] further characterize the stabilization and examine cases both where it is weak and where it is strong. Other research has sought to determine which properties of the fast ions contribute most strongly to their effects; these properties are typically categorized as dilution (of the main ions), non-Maxwellian velocity space distributions, and electromagnetic effects. Tardini etal and Holland etal [25, 26] studied the effect that dilution has on suppressing ITG turbulence, with [27] presenting a reduced model for this effect (expanded further in this work). The active response of a hot Maxwellian impurity was investigated analytically [28] and in simulations [29]. The effects of non-Maxwellian fast ions on turbulence were initially studied using isotropic fast ion velocity-space distributions in [30] and was further generalized for anisotropic fast ions in [31]. The electromagnetic effects of fast ions consist a rich topic, primarily focused on the destabilization of Alfvén eigenmodes (AEs) and energetic particle modes (EPMs) [32, 33]. The theory of how electromagnetic fluctuations interact with ITG turbulence is continuing to be developed [34–39]. In the present work, an analytic model for the effect of fast ions on ITG turbulence (in contrast to modes which are driven by fast ions themselves) is presented. We seek to explain the dominant mechanism for fast ioninduced stabilization of ITG turbulence from first principles, focusing on a case where the effect was very strong: discharge 73224 of the Joint European Torus (JET) [40]. The results of a comprehensive linear study are presented in section2, in which it is examined to what extent the stabilization can be characterized as a change in the linear growth rate. Much about the phenomenon will be learned by doing this because it allows a wider variety of high-resolution simulations to be performed than would be unfeasible with nonlinear turbulence simulations. Then, in section 3, the basic problem is simplified further by approximating the fast ion distribution function in the energetic limit, which yields an approximate analytic solution to the gyrokinetic equation. With the most important elements distilled, this simplified model is inserted into Maxwell’s equations, in which an effective parameter model, applicable to linear and nonlinear simulations alike, makes itself evident. The electrostatic and electromagnetic effects are modelled respectively as effective modifications to the temper ature ratio τ=Ti/Te and βe = 8 π neTe / B2 , parameters to which microturbulence is known to be sensitive. This model is then benchmarked against nonlinear gyrokinetic simulations and elaborated upon further. 2. Characterizing the effect of fast ions on the linear ITG mode We begin by considering the simplified linear gyrokinetic system and ask whether, and how, the presence of fast ions affects the growth of unstable ITG modes that give rise to microturbulence. We accomplish this by analyzing a case in which their effect is particularly strong: JET discharge 73224 [19]. We will find that the presence of fast ions has a nontrivial effect on the ITG mode beyond their global effects on the plasma and their dilution of turbulence driven by the thermal ions3. In this sectionmany of the results and analyses presented are a summary of the more extended linear analysis expounded in [41]. When the fluctuations are small enough to treat linearly, ITG may cause the electromagnetic fields to grow exponentially. Thus we find, for example, that the fluctuating electric potential φ∝e−i˜ωt , where ˜ω=ω+iγ , ω is the frequency (by convention, positive for waves propagating in the ion diamagnetic direction), and γ (if positive) is the growth rate. In this section, we repeatedly calculate this growth rate because, although eventually the fields grow until the nonlinear interaction between modes competes with the linear physics [42] and the system reaches a saturated turbulent state, the linear growth rate has a direct impact on the strength of the saturated turbulence [43]. Indeed, it is often used for the estimation 3 The ions which are approximately Maxwellian close to the electron temperature which make up most of the positive charge and drive the ITG mode will alternatively be referred to as ‘bulk’, ‘thermal’, or ‘main’ ions throughout this work. Nucl. Fusion 58 (2018) 082024 G.J. Wilkie etal 3 of saturated field amplitudes in reduced transport models [44, 45]. Now in order to study the ‘effect of fast ions on the linear growth rate’, it is necessary to somewhat artificially isolate particular effects because, without generating fast ions, the discharge would be fundamentally different. In other words, since fast ions are used to heat the plasma, drive the current, and act as a particle source, it is not immediately clear how to compare cases ‘with’ and ‘without’ fast ions. In this section, we define this comparison as being between a baseline case that includes fast ions and a hypothetical equivalent case in which fast ions do not participate in the turbulence, but in which global properties of the main ions and electrons are nevertheless the same, except for constraints imposed locally by the absence of fast ion charge. To be effective at heating, the fast ion pressure ought to be comparable to the thermal pressure. Therefore an additional effect on the magnetic geometry is expected. Our observations indicate that fast ions have little effect on the flux surface shape, but have a non-trivial impact on the safety factor, magnetic shear, and Shafranov shift, which are all known to play a significant role in microturbulence. For more details on the stabilizing role of the fast ion pressure gradient, the reader is directed to [22, 24, 41]. Henceforth, in order to to isolate the effect of fast ions on the local ITG mode, the equilibrium plasma parameters, including the safety factor, magnetic shear, and flux surface shape, will remain fixed. 2.1. Gyrokinetic framework and baseline case Apart from these global effects, there remains a nontrivial effect of fast ions on the ITG mode as manifest in local gyrokinetics. The simulations presented in this work were performed with the gs2 code [9, 46], which solves the gyrokinetic equation[6]: ∂h s ∂t+vb·∇hs+vD·∇hs+vχ·∇hs−C[hs ] = − Zse∂χR ∂t ∂F0s ∂E− vχ ·∇ F0s (1) for the perturbed distribution function for several isotropic species s: δ fs=Zseφ∂ F 0s ∂E +hs. While we will focus on simulations of isotropic fast ions in this work, we will occasionally discuss the implications of anisotropy. In this case, the gyrokinetic equation is the same, but there is an additional contribution to δfs from the possible μ-dependence of F0s. The mass and charge of the species are ms and Zse respectively. The fluctuating fields φ and A (the vector magnetic potential) are represented by a scalar electromagnetic potential χ≡φ−v·A/c . The equilibrium distribution of species s is F 0s E, µ , σ  , and hs is the non-adiabatic part of the fluctuating distribution function which does not depend on gyrophase ϑ. The equilibrium magnetic field has magnitude B and points in the direction of the unit vector b≡B/B . The parallel velocity is defined as v =σ  2[ E− µB(θ)] /ms, with energy E , the exactly-conserved magnetic moment μ, and sign σ=±1 . For a Maxwellian species with a temperature Ts, the thermal speed vts ≡2Ts / ms , and for a non-Maxwellian species, this represents a characteristic speed defined using an effective temperature T ∗ f≡− nf  (∂F0f/∂ E )d3v −1 . The characteristic Larmor radius is given by ρs≡vts/Ωs , where Ωs≡ZseB/msc . The tokamak minor radius a provides an approximate length scale on which the equilibrium and fluctuations parallel to the magnetic field vary, and ρ∗≡ρi/a1 . The quantity vD is the magnetic drift velocity of the guiding center and vχ≡(c/B)b×∇χR includes the E×B drift along with the drifts and streaming associated with the fluctuating magnetic field. The gyro-average at fixed guiding center R is denoted by  φ R = 2 π 0 φ( r )dϑ/2π. For simplicity, gradients of equilibrium plasma flows are ignored and the equationis solved in the frame rotating with the plasma. A conservative linearized Fokker–Planck collision operator C [47, 48] is employed to model collisions of ions and electrons, and of each species with themselves (collisions between ions of different species are omitted). Throughout this work, we will focus on JET discharge 73224 around the flux surface with a half-width of r = 0.375a (where a is the half-width of the last closed flux surface). The nominal parameters for this baseline case are based on those from [49] and are listed in tablesA1 and A2. The resolution for the linear simulations of this sectionare as follows: there are 58 grid points along the field line per poloidal turn (parametrized by poloidal angle θ), and the parallel domain extends to θ=±22π . The velocity space resolution is 36 grid points in energy and 46 in pitch angle. The time step is 0.03 a/vti, run until the complex frequency is conv erged to within a factor of 2×10−4 . There are two species of fast ions present in these simulations, each represented as a hightemperature Maxwellian. The local pressure gradient, used to rescale the local geometrical parameters according to the Miller prescription, is calculated consistently depending on the local fast ion pressure gradient. Unless otherwise stated, all simulations include perpendicular magnetic fluctuations ( A ), with compressive fluctuations ( δB ) artificially disabled. This is justified by the relatively low βe , and by the growth rates shown in figure1. Figure 1. Comparing ITG/KBM growth rates with fast ions with (green pentagons) and without (orange squares) compressive magnetic fluctuations (case B0 of tableA2) at kyρi=0.57 . Nucl. Fusion 58 (2018) 082024 G.J. Wilkie etal 4 At a critical βe , the kinetic ballooning mode (KBM) becomes unstable and the growth rate significantly increases. Only then is there a discernible difference when including δB . As long as the mode is ITG-like, increasing βe is stabilizing and only the fluctuations of A need to be considered. For a more detailed treatment on the effect of compressive fluctuations on ITG modes, see [38]. 2.2. Fast ion-induced stabilization of the ITG mode: general observations The effect of fast ions on the ITG mode growth rate and frequency in this case is shown in figure2. The baseline case has both NBI fast deuterium and ICRH fast helium-3, with both types together consisting 20% of the positive charge. Then, the growth rates were recalculated for a case with the different types of fast ions removed. A significant increase in the ITG mode growth rate is observed, with a modest change in the frequency spectrum. Finally, both fast ions were removed to generate a case ‘without fast ions’. This shows that the effect of NBI is relatively small compared to ICRH. While the growth rate is sensitive to the presence of ICRH fast ions, the ITG mode frequency and eigenfunction are not significantly altered by the presence of fast ions, as was also observed in [49]. A ‘mixing-length’ estimate for the bulk ion heat flux based on the linear physics is qi∼ γ/ k2 ⊥ . Using this estimate, our results suggest about a factor of 2 increase in the turbulence amplitude when fast ions are removed. This alone does not account for the strong nonlinear effect in [19] and later in figure8. These latter results showed about a factor of 10 increase in the thermal ion heat flux. It is in this sense that we say the effect of fast ions on turbulence cannot be explained in the context of linear gyrokinetics. The decrease in the ITG mode growth rate when fast ions are included as a kinetic species, and the subsequent (disproportionate) decrease in the turbulence amplitude, is the subject of the remainder of this paper. We begin with the simplest explanation for this phenomenon: that only the indirect effects of fast ions are responsible for the stabilization. Namely, we examine cases where the only effects of fast ions considered are their effects on the bulk ion and electron densities. 2.3. Dilution On the temporal and spatial scales of interest, the plasma remains quasineutral:  s Z s en s= 0. (2) When positively-charged impurities are included, ni=ne , and this means that the electromagnetic fields have proportionally less response to the bulk ions. Since the bulk thermal ions are responsible for the instability of the ITG mode, this can lead to a reduction in the vigor of the turbulence. This effect is known as dilution. For the same reason, dilution of the ions implies an enhancement of electron-driven microinstabilities, such as the trapped-electron mode [50] and the electrontemperature-gradient mode (ETG) [46]. In contrast to thermal high-Z impurities, which are known to stabilize ETG [51], the adiabatic response of singly-charged fast ions actually reduces Z( ETG ) eff ≡i Z2 i (n i /n e )(T e /T i) due to their high temperature. Taking the radial derivative of equation(2) gives:  s Zse∂ n s ∂r≡−  s Zse n s Lns = 0. (3) In a local simulation, this provides a second independent constraint that must be satisfied. This represents dilution on neighboring flux surfaces, which influences the relationship between the electron and thermal ion density gradients. This also has a nontrivial effect on ITG turbulence and is included when we speak of ‘dilution’. When results are thereby labelled, it means that fast ion fluctuations are not included in the gyrokinetic simulations, but their equilibrium effect on the density of the bulk species and on the magnetic drift (through their contribution to ∂β/∂r ) are included. In other words, we keep the effects of F0f and ∇F0f , but not δff . Note that, in the high-energy limit (explored in section3), this is equivalent to the adiabatic approximation, where δff≈hf≈0 . In constructing artificial cases lacking fast ions (in order to isolate and study their effect on the ITG mode and resulting turbulence), there is a choice to be made regarding whether the electron or thermal ion density (and density gradients) are changed to maintain quasineutrality. Fusion products and energetic tails heated from the bulk ion population deplete the thermal ion population. On the other hand, the physical origin of injected and minority-heated fast ions are such that they are accompanied by excess electrons. The difference between these choices for the kyρi=0.57 ITG mode is shown in figure3. As can be seen, this choice happens to have a relatively small effect on our results, although this may not be true in general. In this section, we choose the convention that it is the electron properties that are changed. Figure 2. ITG growth rate and frequency spectra with different concentrations and types of fast ions. Orange circles show the case without any fast ions (case BXe of tableA2), black/green pentagons show a case with only the NBI/ICRH fast ions with Zfnf=0.06/0.14 (cases BN/BI, respectively), and violet triangles show the baseline case (B0) with both NBI ( Zfnf=0.06 ) and ICRH ( Zfnf=0.14 ). Nucl. Fusion 58 (2018) 082024 G.J. Wilkie etal 5 2.4. Differentiating classes of fast ions The different types of fast ions are characterized by whether their strong radial variation is one of particle density (‘NBI- like’) or energy density (‘ICRH-like’). It has previously been shown that these different classes of fast ions respond differently to turbulence when passive [27, 52]. It is worth examining to what extent their effect on turbulence differs compared to their respective dilution effects. This contrast is shown in figures4(a) and (b). The violet spectra each include one of the respective fast ion species, and the dashed black spectra are the cases without the fast ions, for which the electron density and density gradients are changed to maintain quasineutrality. Then, cases identical to the violet (‘with fast ions’) were run, except that fast ions were not included as a kinetic species; only their effect on the equilibrium ne and Lne. These ‘dilution’ cases are shown in green. We see that, relative to dilution, NBI-type fast ions are actually destabilizing, whereas ICRH-type fast ions are more strongly stabilizing than dilution only. Therefore, the classification of fast ions, whether LTf Lnf (NBI) or Lnf LTf (ICRH), is critical to predicting and understanding their stabilizing effect. Note that alpha particles are considered NBI-like, but instead of being artificially injected, they are produced by the fusion reaction. Because the fusion source is a strong function of radius, the alpha particle density gradient is sharp. Now the fact that the two cases with fast ions differ from their respective ‘dilution only’ cases demonstrates that fast ions play a non-trivial role besides mere dilution. However, when both types are present, the (significant) kinetic effects of fast ions may cancel out and one can be left with the illusion that dilution is the dominant effect. We saw from figure4 that, relative to dilution, the NBI- like fast ions (strong density gradient) are destabilizing, while the ICRH-like fast ions (strong temperature gradient) are stabilizing. The parameter η f ≡T f(r)/n f(r)=L nf /L Tf characterizes the relative strengths of the gradients. Since the stabilization is clearly a function of this parameter, its effect is examined in figure5. Here, the temperature gradient of fast deuterium is modified from the baseline case and the growth rate is compared to the case of pure dilution (black lines). Note that the dilutive effect on the bulk plasma density is not affected when changing the temperature gradient. For both electrostatic and electromagnetic cases, the threshold in the value of ηf for when the fast ions become more stabilizing than simple dilution is around 0.7–1.0. This threshold will be explained with an analytic model in section3.4 after a model perturbed distribution function is obtained. The effect of fast ions in local gyrokinetic turbulence simulations is often ascribed to electromagnetic effects. Here, some qualitative differences between different kinds of fast ions at different values of βe are catalogued. This is complicated by the fact that βe itself has its own electromagnetic stabilizing effect, and separating which changes to the growth rate are due to the electromagnetic stabilization and which are due to fast ions is not trivial. Therefore, scans in βe for the growth rate of the kyρi=0.57 mode, using several different fast ion parameters, are shown in figure6. In all these cases, increasing βe is stabilizing. However, increasing the Figure 3. ITG mode growth rates at kyρi=0.57 , as a function of the driving temperature gradient. In one case, we assume that the thermal ions replace the removed fast ions (orange circles—case BXi of tableA2). In the other, electrons are removed along with fast ions (green pentagons—case BXe). Figure 4. ITG mode growth rate spectra when NBI and ICRH fast are included (violet circles), compared to to dilution (green pentagons) and no fast ions (black dashed). The black dashed line in both cases is case BXe of tableA2. Part (a) shows cases BN and BND, whereas part (b) shows cases BI and BID. Figure 5. Scan in fast ion η=T(r)/n(r) for electrostatic (green circles) and electromagnetic (violet triangles) ITG simulations relative to electrostatic/electromagnetic dilution (dashed/solid black line, respectively) for case B0 of tableA2 at kyρi=0.57 . Nucl. Fusion 58 (2018) 082024 G.J. Wilkie etal 6 fast ion pressure gradient (either from the density gradient or the temperature gradient) does not universally decrease the ITG mode growth rate. In particular, note the case of NBI- like fast ions in figure6(b) (a/LT,f = 0) at low βe : increasing a/Ln,f is actually net-destabilizing electrostatically. Dilution cases are not shown in figure6, but these still have a lower growth rate than their corresponding high-a/Lnf cases, and is approximately equivalent in the a/Lnf = 0, a/LTf = 0 case, as expected. Therefore, the stabilization effect of fast ions and dilution itself is somewhat sensitive to βe [53]. The perturbed fast ion parallel current is small, but the direct effect of fast ions on the electrostatic potential could couple to A , causing the behavior shown in figure6. We have explored several potential explanations for the fast ion stabilization including: a modified linear growth rate, dilution of bulk ions, changes to magnetic geometry, and magnetic fluctuations. Although these explanations are physically motivated and in aggregate have a non-trivial effect, they are found lacking in describing the full order-of-magnitude stabilization observed in JET discharge 73224. A more fundamental treatment, valid in fully nonlinear turbulence, is thereby motivated. In the next section, a simplified model for the gyrokinetic equation, valid in the high-energy limit, will be derived and this will later be used to develop a reduced model to explain the qualitative effects presented in this section, in addition to the effect of fast ions in nonlinear simulations. 3. Reduced model for fast ion effect on microturbulence In the previous sectionit was demonstrated that the effect of fast ions in stabilizing the ITG mode can go beyond dilution and depends on the details of the equilibrium fast ion phase space distribution. In this section, we examine the leadingorder behavior of the gyrokinetic equationin the high-energy limit. This subsidiary expansion yields an analytic solution that applies rigorously to fast ions with a strong radial dependence. Then, this model for hf is reduced further with additional assumptions based on the ITG mode structure. When this approximate distribution function is inserted into Maxwell’s equations, this leads to a physically-transparent effective parameter model. We then discuss this model and benchmark it against linear and nonlinear gyrokinetic simulations. 3.1. Energetic expansion of the gyrokinetic equation In order to obtain a nonlinear model for the effects of fast ions, we will directly expand the gyrokinetic equation, making use of the high-energy nature of the fast ions. Indeed, it is common to perform such expansions for electrons, taking advantage of their small mass. Commonly-used models include: driftkinetic [54–56], fluid [57, 58], bounce-averaged [59], or adiabatic electron [60] models. The adiabatic electron model, for example, approximates their contribution to the field equationsas proportional to the electrostatic potential. In this section, an analogous model, applicable to energetic ions, will be developed. In the energetic limit of ≡vti/vtf 1 , the gyrokinetic equation(1) reads: vb·∇h f +v D ·∇h f =−v χ ·∇F 0f . (4) In arriving at equation (4), ion-scale microturbulence is in mind so that the bulk ions are what set the temporal and spatial scales of the turbulent fluctuations so that ∂hf/∂t∼ωhf and |∇hf|∼hf/ρi . The ∂hf/∂t term, the nonlinear term, and the first term on the right hand side of equation(1) are smaller than the magnetic drift term by factors of −2 , −3 / 2 , and −3 , respectively. The radial gradients of the fast ion equilibrium are ordered to be strong such that |∇ F 0f|∼O − 3 / 2 F 0f /a  . To leading order, only the magnetic drift term survives and we recover the adiabatic approximation: hf≈0 , which is equivalent to ‘dilution’ in this limit. To find nontrivial effects, we therefore wish to find a solution for hf correct up to O() , which is why the other terms are retained in equation(4) even though they are formally smaller than the magnetic drift term by one power of ε. As is customary, an eikonal representation is chosen for the fluctuations [61] so that χ=ˆχ(θ)eiS , where b·∇S=0 and the ballooning angle θ∈(−∞,∞) has been chosen as the coordinate along the magnetic field line. Hats will be dropped henceforth on h and χ where there is no ambiguity. The perpend icular spatial dependence of the fluctuations is embedded in S, which depends on the poloidal magnetic flux ψ (which labels the flux surface and is a useful ‘radial’ coordinate) and a field line label α such that B=∇ψ×∇α . We perform a Fourier transform in the plane spanned by ∇ψ and ∇α so that the gyroaverage is represented by the Bessel function J0s≡J0(k⊥v⊥/Ωs) . Now equation(4) becomes: v (b ·∇ θ) ∂h f ∂θ +i(vD ·∇ S)hf=ic ∂S ∂α ∂F 0f ∂ψ J0fχ . (5) This can be solved directly for hf. For the integrating factor, it will be convenient to define: z(θ)≡ θ θ0 ω Ddθ v  (b ·∇ θ) , (6) Figure 6. Scan of growth rate versus βe for different fast ion species: (a) ICRH-like, and (b) NBI-like, each with several different gradients of temperature and density, respectively. These are iterations of cases BI0 and BN0 of tableA2 with the fast ion gradients adjusted as indicated. For Case BN0, electron density and density gradient are likewise adjusted to maintain quasineutrality. Nucl. Fusion 58 (2018) 082024 G.J. Wilkie etal 7 where ωD(θ)≡vD·∇S . The primed variables denote which θ is the independent variable such that ω D≡ωD(θ) (similarly for v,b·∇θ , χ, J0, and z). In equation (6), the lower limit of integration, θ0 , is defined to be a point where ωD vanishes such that θ0<θ if v>0 , and θ0>θ if v<0 . Equation(5) can then be rewritten as ∂ ∂θ  eizhf  =eizic ∂S ∂α ∂F 0f ∂ψ J 0f χ v b·∇ θ , (7) whose solution is: h f=ic∂ F 0f ∂ψ ∂S ∂α θ − σ ∞ J 0fχei(z−z)1 (b ·∇ θ) dθ v  . (8) This solution is largely inspired by that of [34] and can alternatively be derived directly therefrom with the appropriate approximations. The equilibrium distribution can be factored out because F 0f= F 0f E, µ , σ  is not a function of θ in these coordinates. Equation(8) is the full solution of the gyrokinetic equationfor fast ions correct to O(vti/vtf ) . It is a complete linear model in the sense that it is the furthest one can take the (vti/vtf ) expansion linearly; the formally next-largest term that would appear in equation(4) is the electromagnetic nonlinear term −v / Bb×∇AR . Without further analysis, we can read off one important implication of this solution for hf: there is no zonal reponse of the fast ions. This follows from the fact that the zonal modes are defined as those without any α variation (in the context of gyrokinetic simulations, this is often written as ky = 0). Hence, for such modes ∂S/∂α =0 , and so hf(ky=0)=0 . Zonal flows (which arise from zonal fluctuations of φ) are well known to be important in the nonlinear saturation of ITG turbulence. Our approximate solution shows that fast ions have no direct impact on the zonal fields. To demonstrate the robustness of this result, simulations have been performed in which hf(ky=0)=0 is artificially enforced. There, no discernible difference in the nonlinear fluxes was found; see figure8. This result is not without consequence. The saturated level of turbulence is determined by a balance between ‘drive’ and zonal flows, which interact nonlinearly. Fast ions only directly affect the former and not the latter. Decreasing the drive will allow stronger zonal flows [53] and an overall damping of the turbulence: moreso than if all modes (including the zonal modes) were directly damped. A similar mechanism applies to an adiabatic electron response, which also vanishes for the zonal mode [60], and this makes turbulence remarkably sensitive to the ion-electron temperature ratio [62]. This is a possible explanation for the nonlinear enhancement of the fast ion stabilization. Zonal flows play an important role in the story, but it is precisely because fast ions do not have a zonal response. 3.2. Simplifying the model The rigorous solution, equation (8), still retains too much physics to be a useful reduced model. We can make further approximations which allow us to write hf∝χ . When moments of hf are thereby taken, the proportionality factors become response functions, which we will find very useful in interpreting the contribution of fast ions to the fluctuating electromagnetic fields. Integrate equation (8) directly by parts to obtain the approximate solution: h f=c∂S ∂α ∂F0f ∂ψ  J0fχ ωD −  θ −σ∞ ei(z−z)∂ ∂θ  J 0χ ω D  dθ  ≈ c∂S ∂α ∂F0f ∂ψ 1 ωD J0fχ, (9) where the approximation is made because the second term is smaller by one power of ε. A difficulty with this approximation is that ωD vanishes at specific θ points. We thus need to assume that χ vanishes sufficiently rapidly away from θ=0 that the contributions from such resonances are negligible due to the smallness of χ. This is referred to as the ‘strongly ballooning’ or ‘outboard mid-plane localization’ approximation. It is valid to the extent that the ITG mode structure peaks at the outboard midplane, which is typically the case. In practice, this is equivalent to approximating the integral  θ −∞ ≈ θ θ0+∆θ , stopping just short of the nearest resonance point θ0 . Alternatively, the solution (9) could also be found by ignoring the (formally small) parallel streaming term in equation(4), promoting ∇F0f by an additional order in ε, and solving for hf algebraically. The approximation in equation(9) is supported by the fact that it exhibits the correct energy dependence that we expect from the usual scalings of energetic particle transport in electrostatic and electromagnetic turbulence [30, 63–65]. The leading term in equation(9), which we will use for the fast ion response to fluctuating fields, does not contribute to the turbulent transport because it is in-phase with χ. The correction from the second term results in a nontrivial phase factor, which does contribute to turbulent transport with the correct energy dependence. 3.3. Model fast ion response function In this section, the model distribution function of the previous sectionis reduced further to remove the θ dependence, taking θ=0 as the only relevant location for ωD and J0f. In this case, the vj  moments of hf become related to numerical parameters Rjf, which can be transparently interpreted in the gyrokinetic field equations. First, consider the relevant Maxwell’s equations in the gyrokinetic limit: sZsδns=0 (Poisson’s equation/quasineutrality) and ∇2 ⊥A = − ( 4 π/ c ) s δ j s (the parallel component of Ampere’s law with δj s ≡Z s e δ f s rvd3v being the contribution of each species to the perturbed parallel cur rent). Assume the plasma consists of a thermal ion species, electrons, and fast ions. Considering the low mass and high thermal speed of electrons, their contribution to the perturbed current dominates over the ions. At low βe , the electron contrib ution to the perturbed charge density is approximately their adiabatic response to the electrostatic potential: Nucl. Fusion 58 (2018) 082024 G.J. Wilkie etal 8 δne=(nee/Te)  φ − φ  ψ  , where ψ denotes the flux surface average. The field equationsbecome: φ Z2 ie2ni Ti +e2ne Te +Z 2 fe 2 nf T∗ f  −nee2 Te φψ =Zie  J0ihid3v+Zfe  J0fhfd3v , (10) ck2 ⊥ 4π A= − e  J0ehevd3v+Zfe  J0fhfvd3v . (11) In equation(11), it was assumed that the thermal ion contribution to the parallel current is small compared to that of the electrons or the fast ions. Equation (11) is not inconsistent with employing the adiabatic electron approximation in equation(10). This is because the leading adiabatic behavior of electrons, even if dominant at low βe , does not contribute to the parallel current, while non-adiabatic corrections do contribute. With singly-charged bulk thermal ions and a temperature ratio τ≡Ti/Te , the non-zonal components of quasineutrality (we have already established that the zonal fast ion contribution is negligible) can be written: φ ni ne +τ+Z2 f nf ne Ti T∗ f  =Ti ene  J0ihid3v+Zf  J0fhfd3v . (12) To capture the gyrokinetic effect of fast ions in response to the fluctuating fields, let us define the response functions: Rjf =Zf cTi en e ∂S ∂α  v v ti j 1 ω D0 J2 0  k⊥0v⊥ Ω f0 ∂F0f ∂ψ d3v, (13) where ωD0≡ωD(θ=0) (similarly, k⊥0 and Ωf0 take their values at θ=0 ). Here, the mid-plane localization assumption of the previous section was taken further to assume that vD takes its values at the outboard midplane. This assumption is also applied to v , b·∇θ , and Ωf , which is equivalent to taking the large aspect ratio approximation and ignoring trapped particle effects. Up until now, the magnetic geometry and equilibrium fast ion distribution are general. However, for the sake of straightforward parametrization, let us make the further assumptions of circular geometry and Maxwellian fast ions: Rjf ≈Z2 f Ti Tf nf ne R 2Lnf v vti j1+ηf  (v/vtf )2−(3/2)  v2 +v2 ⊥/2/v2 tf ×J2 0  k⊥0v⊥ Ωf0  e−v2/v2 tf π3/2v3 tf d3v. (14) This is the form of the response function that will be used henceforth. Note that R1f vanishes by odd symmetry in v when F0f is isotropic. This symmetry is important for the model that follows because the latter depends on the fast ions not coupling the electrostatic and electromagnetic field equations(this coupling is precisely R1f). Even if this symmetry is violated, one could still argue that the fast ion current is small (which removes their relevance in Ampere’s law) and, when βe is small, vA/cφ , which means that the prefactor on R1f is small in quasineutrality, equation(10). The response functions R0f and R2f are shown in figure7. Note that k⊥ρf and ηf consist the only nontrivial parameter dependency of the response functions; all other parameters appear as prefactors in equation(14). When |R0f|1 , this means that dilution is the dominant electrostatic fast ion effect by definition. This allows us to make an important conclusion for alpha particles: unless accompanied by an unphysically strong radial gradient, the prefactor of nα/Tα in equation(14) makes the electrostatic kinetic response of alpha particles very weak, in agreement with the results of [30]. We will find that the electromagnetic response R2f has a relatively weak effect at thermal ion scales, even though it comes with an additional prefactor of Tf/Ti . 3.4. Effective parameter model Having derived a simple fast ion response function in equation(14), we proceed to interpret it in the context of the gyrokinetic field equations. This is done by generalizing the model presented in [27] in a way that goes beyond dilution and electrostatic turbulence. Consider that hi is proportional to ni, but otherwise does not depend on the equilibrium ion density, except through the calculation of φ in equation (12). After multiplying equation(12) by ne/ni , it is seen that the same φ is ensured to be calculated from hi for all cases where the quantity ( ne/ni)  (ni/ne)+τ+Z2 f(nf/ne)  Ti/T∗ f − R0f  is constant. Equate this quantity to a case without fast ions where ni=ne , but with an artificially defined τeff . This defines an effective temperature ratio which mimics effect of fast ions on microturbulence: τ eff =ne ni  Ti Te +Z2 f nf ne Ti T∗ f −R0f . (15) This is a generalization of the dilution model presented in [27], which did not include the kinetic effect of fast ions Figure 7. (a) Electrostatic and (b) electromagnetic response functions for fast ions as functions of k⊥ρf . Select values of ηf=0, 1, 5 are shown in cyan, blue, and black respectively. The parameters that appear as prefactors in equation(14) have been normalized out. Note that R2f is multiplied by an additional factor of Tf/Ti . Nucl. Fusion 58 (2018) 082024 G.J. Wilkie etal 15 [57] LinZ. and ChenL. 2001 A fluid-kinetic hybrid electron model for electromagnetic simulations Phys. Plasmas 81447 [58] ChenY. and ParkerS. 2001 A gyrokinetic ion zero electron inertia fluid electron model for turbulence simulations Phys.Plasmas 8441 [59] GangF.Y. and DiamondP.H. 1990 A nonlinear bounce-kinetic equationfor trapped electrons Phys. Fluids B 22976 [60] DorlandW. and HammettG.W. 1993 Gyrofluid turbulence models with kinetic effects Phys. Fluids B 5812 [61] AntonsenT.M. and LaneB. 1980 Kinetic equationsfor low frequency instabilities in inhomogeneous plasmas Phys. Fluids 231205 [62] PettyC.C., WadeM.R., KinseyJ.E., GroebnerR.J., LuceT.C. and StaeblerG.M. 1999 Dependence of heat and particle transport on the ratio of the ion and electron temperatures Phys. Rev. Lett. 833661 [63] HauffT., PueschelM.J., DannertT. and JenkoF. 2009 Electrostatic and magnetic transport of energetic ions in turbulent plasmas Phys. Rev. Lett. 102 [64] ZhangW., DecykV., HolodI., XiaoY., LinZ. and ChenL. 2010 Scalings of energetic particle transport by ion temperature gradient microturbulence Phys. Plasmas 17055902 [65] PueschelM., JenkoF., SchnellerM., HauffT., GünterS. and TardiniG. 2012 Anomalous diffusion of energetic particles: connecting experiment and simulations Nucl. Fusion 52103018 [66] CourantR., FriedrichsK. and LewyH. 1967 On the partial difference equationsof mathematical physics IBM J. Res. Dev. 11 215 [67] TegneredD., OberparleiterM., StrandP. and NordmanH. 2017 Impact of a hollow density profile on turbulent particle fluxes: gyrokinetic and fluid simulations Phys. Plasmas 24072303 [68] PorkolabM. etal 2012 Transport and turbulence studies in the linear ohmic confinement regime in Alcator C-Mod Plasma Phys. Control. Fusion 54124029 [69] PusztaiI., MollénA., FülöpT. and CandyJ. 2013 Turbulent transport of impurities and their effect on energy confinement Plasma Phys. Control. Fusion 55074012 [70] EnneverP. etal (Alcator C-Mod Team) 2015 The effects of dilution on turbulence and transport in C-Mod ohmic plasmas and comparisons with gyrokinetic simulations Phys. Plasmas 22072507 [71] FuG.Y. 2008 Energetic-particle-induced geodesic acoustic mode Phys. Rev. Lett. 101 Nucl. Fusion 58 (2018) 082024