Energetic particle effects as an explanation for the low frequencies of Alfvén modes in the DIII-D tokamak

During beam injection in the DIII-D tokamak, modes with lower frequencies than expected for toroidicity induced Alfvén eigenmodes (TAEs) are often observed. The experimental `TAE' frequency is often ≈0.8 of the nominal theoretical frequency of the TAE, fTAE, while the typical frequency of beta induced Alfvén eigenmodes (BAEs) is (0.2-0.4)fTAE. An analysis is presented of an unstable discharge with a high n stability code, HINST, that includes the effect of energetic ions on mode frequency. The analysis shows that the experimental `TAE' and `BAE' could be resonant branches of the TAE and the kinetic ballooning mode, respectively.


Introduction
The toroidicity induced Alfvén eigenmode (TAE) [1] has been studied extensively in tokamaks [2]. Theoretically, the nominal frequency of the ideal MHD TAE branch in a low beta plasma is f TAE = v A /4πqR, where v A is the Alfvén speed, q the safety factor and R the major radius. During intense neutral beam injection into DIII-D plasmas, an instability with a frequency close to f TAE was observed [3], as well as another mode with a frequency about half as large that was dubbed the beta induced Alfvén eigenmode (BAE) [4]. In this article, we call these two experimentally observed instabilities the 'TAE' and the 'BAE' (irrespective of their correct theoretical identification). In their most virulent form, both instabilities cause large losses of fast ions [5].
Subsequent analysis of the 'TAE' found some points of agreement with TAE theory and some points of disagreement. For controlled parameter scans, the mode frequency scales linearly with toroidal field [6] and inversely with the square root of density √ n e [3], as expected for Alfvén waves. The observed stability threshold agrees to within a factor of 2 with calculations that include 'radiative damping' to the kinetic Alfvén wave [7] and other calculations that include mode coupling to drift and kinetic Alfvén waves [8]. On the other hand, comparisons of the mode structure with theory show significant discrepancies [9,10]. Calculations that include kinetic effects agree best with the measurements but other effects, such as modifications of the eigenfunction caused by the energetic ion population, may be required to explain the data.
To date, theory has been less successful in explaining the 'BAE'. Ideal MHD analysis predicted an eigenmode with predominately Alfvénic polarization in a low frequency gap in the MHD continuum caused by plasma pressure [11], which is why the instability was named the BAE in Ref. [4]. Later, three alternative theoretical identifications were proposed: a kinetic ballooning mode (KBM) [12], a mode that propagates at the ion thermal speed [13] and an Alfvén mode whose structure and frequency are modified by the fast ion population. The fast ion modified instability was called an energetic particle mode (EPM) in Refs [14][15][16] and a resonant TAE (RTAE) in Ref. [17]. In fact, two branches of EPMs were discussed in Refs [14,16]: a higher frequency branch related to the TAE and a lower frequency branch related to the KBM. (Ref. [17] discussed only the RTAE branch.) The frequencies of these EPMs lie in the continuum of ideal MHD, where strong damping is anticipated. The continuum damping is minimized near the accumulation points of the shear Alfvén continuous spectrum [18]. For the resonant TAE, the competition between fast ion drive and damping often results in a mode with a frequency near the bottom of the toroidicity induced gap in the Alfvén continuum. For the resonant KBM, the damping is minimized for ω/ω A → 0 (ω A is the Alfvén frequency). Experimental evidence for the existence of both of these  [22] and interferometer [23] measurements, and carbon density (multiplied by 6 (dashed curve)) measured with charge exchange recombination (CER) spectroscopy [24]. (b) Deuterium ion temperature (solid curve) inferred by TRANSP [25]  branches during high energy neutral beam heating was found in JT-60U [19].
In Ref. [20], an extensive database of frequency measurements of DIII-D 'BAEs' was compared with four simple analytical frequency scalings: an Alfvén eigenmode (f ∝ v A ), a KBM (f ∝ ω * i ), a mode that progagates at the ion thermal speed (f ∝ v i ) and an EPM (f proportional to the beam ion circulation frequency). None of these simple scalings fit all of the data. Reference [20] suggests that two factors contribute to this failure. First, the assumption that all the observations in the database are of the same instability is probably erroneous. Second, the simple analytical frequency scalings neglect many of the parametric dependences found in numerical studies, in particular for the EPM. Evidently, numerical calculations of the expected frequency are required on a case by case basis to establish definitively the correct identification of the experimental 'BAE'.
This article reports a comparison of EPM theory with DIII-D data for a representative case. A non-perturbative fully kinetic code HINST [21], which stands for high n stability code, is used to compute the expected frequencies of the RTAE and resonant KBM in a typical DIII-D plasma with slightly negative shear and with unstable 'TAE' and 'BAE' activity. The experimental 'TAE' actually has a frequency just below the toroidicity induced gap in the MHD continuum. As we will show, the code successfully finds unstable modes with frequencies that are close to the measured frequencies.
The calculated radial position and the toroidal mode numbers n of the unstable modes are also close to the experimental observations. The code predicts that the higher frequency RTAE is more unstable than the lower frequency resonant KBM, but, in the experiment, the 'BAE' seems to be slightly more unstable than the 'TAE.'  Frequency spectrum for an O mode reflectometer channel with a cut-off frequency corresponding to n e 3.1 × 10 13 cm −3 , which occurs inside the density pedestal. Contours of the logarithm of the amplitude are shown. The experimental 'BAE' has frequencies in the laboratory frame of ∼100 kHz, while the experimental 'TAE' has frequencies ≈200 kHz. The dashed line approximates the expected laboratory frequency of the TAE, v A /4πqR + nf rot , evaluated for the n = 5 mode near r/a 0.6.

Experiment: Description of the instability
The case selected for comparison is that of a double null divertor deuterium plasma that is heated by 9.5 MW of 76 keV deuterium neutral beams. At the time of interest, the safety factor profile is weakly reversed ( Fig. 1(d)) and the plasma has entered an ELM-free H mode, so the stored energy, neutron rate and density are rapidly increasing. As shown in Fig. 1, the central ion temperature is approximately 10 keV and there is an edge pedestal in the density and temperature profiles. Because the density is still relatively low and the confinement is high, the classically expected beam pressure ( Fig. 1(c)) is a significant fraction of the total plasma pressure. Phenomenologically, 'BAEs' are often observed in plasmas with large beam ion pressures [20].
Instabilities are observed by reflectometer channels and by Mirnov coils.
The spectrograph of the centremost reflectometer channel with a cut-off frequency corresponding to n e ≈ 3.1×10 13 cm −3 is shown in Fig. 2  Mirnov coil that is situated about 45 • below the outer midplane is shown in Fig. 2(a). Three types of coherent magnetic activity are observed. Low frequency (<20 kHz) low n (n = 1-2) modes appear intermittently. Regular bursts of an n = 5 'BAE' with a frequency of 90-100 kHz are also seen ( Fig. 2(b)). In some cases, coherent n = 3, 4 and 6 modes also appear in this frequency band. Bursts of 'TAE' activity occur between 150 and 250 kHz ( Fig. 2(c)). Virtually all these bursts have a dominant n = 4 or n = 5 mode with a frequency between 170 and 200 kHz, but additional peaks above 200 kHz appear occasionally. The 'BAE' at about 95 kHz and some of the 'TAEs' at about 180 kHz are observed on several reflectometer channels but, because of the weak density gradient in the plasma interior, it is not possible to reconstruct the spatial eigenfunction from the available data.
In previous work, the appearance of Alfvén modes generally correlates with reductions in the volume average neutron rate [20]. In this discharge, the measured neutron rate is 90-100% of the classically expected rate predicted by the code TRANSP [25] so global losses of beam ions are small. Any correlation between MHD bursts and reductions in neutron rate is undetectable ( 1%), so it is reasonable to assume that the beam pressure profile is similar to the one computed by TRANSP.
Comparison of the experimental mode frequency f lab (which is measured in the laboratory frame) with the predicted frequency (which is computed in the plasma frame) is complicated by a Doppler shift, which causes the inferred experimental frequency in the plasma frame f pl to be a function of the assumed radial position of the mode [29]. For strongly rotating plasmas, the measured frequency in the plasma frame is approximately f pl f lab − nf rot , where f rot is the toroidal rotation frequency of the plasma. If a cluster of toroidal modes appears and if one assumes that the different toroidal modes  are excited near the same minor radius, it is possible to infer f pl uniquely from the measured values of f lab and n [29]. In this discharge, multiple 'BAE' and 'TAE' peaks occur on some bursts (Fig. 2). In Fig. 3, the inferred frequencies f pl are plotted versus position for these 'BAE' and 'TAE' peaks. Since the frequency of the TAE in the plasma frame is expected to be independent of mode number, we assume that the modes are localized where the Doppler corrected frequencies coincide. The measurements imply that both modes are excited at a normalized radius near 0.55, that the 'BAE' frequency is f pl = 14 ± 9 kHz and that the 'TAE' frequency is f pl = 108 ± 9 kHz.
The dominant contribution to the uncertainty in f pl is associated with uncertainty in the radial location of the mode δr. The uncertainty associated with other quantities is smaller. For the ensemble of bursts near 1190 ms, the variation in f pl due to uncertainties in f lab and f rot at fixed radius is only 2.9 kHz for the 'BAE' and 5.5 kHz for the 'TAE'. (The uncertainty in n is negligible.) However, the uncertainty in radius also effects f pl through the Doppler shift correction because of shear in the toroidal rotation profile f rot . If we assume δr = 0.04, then δf pl 7 kHz from this effect alone for an n = 5 mode; adding this uncertainty in quadrature with the others gives an estimated uncertainty of δf pl 9 kHz. There is also a systematic uncertainty of ∼10 kHz associated with the assumption that different toroidal mode numbers have the same frequency in the plasma frame.

Theory: Alfvén spectrum analysis
The Alfvén stability of this discharge at 1190 ms (Fig. 2) is analysed with the HINST code [21], which is a nonperturbative fully kinetic code. It is able to reproduce both the RTAE and resonant KBM branches with a drive from  Fig. 2; the solid circles at the intersection of these curves represent the approximate frequencies and radial locations of the modes in the plasma frame. The central Alfvén frequency is ω A0 = 1.23 × 10 6 rad/s. fast particles. Calculations of mode drive and damping include bulk plasma and fast particle finite Larmor radius (FLR) effects. Radiative damping supported by trapped electron collisional effects and ion Landau damping are also incorporated. Even though HINST robustly finds solutions with high toroidal n numbers that have radially localized mode structures, it can be used for medium n to low n modes in the local version, i.e. without resolving the two dimensional (2-D) mode structure. HINST employs a shooting technique to find the mode frequency, growth rate and one dimensional (1-D) mode structure in ballooning co-ordinates. Note that the global HINST 2-D solution requires radial localization of the mode and high toroidal mode numbers n.
HINST uses the s-α model for plasma equilibrium [30], which assumes an isotropic plasma pressure. Since RTAEs have a ballooning structure similar to that of TAEs, i.e. with the maximum of the mode amplitude at the low magnetic field side of the plasma, the local equilibrium can be approximated as isotropic. Empirically, the 'BAE' is destabilized by circulating beam ions [20], so the beam ion distribution function is approximated as a slowing down distribution consisting solely of passing ions.
The results of the analysis are shown in Fig. 3. The code finds two unstable branches. The higher frequency RTAE corresponds to the experimental 'TAE', while the lower frequency resonant KBM corresponds to the experimental 'BAE'. For both branches, the maximum theoretical growth rate occurs for r/a 0.6, which is close to the mode location inferred from the experimental Doppler shift analysis.
Consider first the results of the local (1-D) calculations for the RTAE, which are represented in Fig. 3(a) by crosses. In the region of low fast particle pressure gradient, the frequency of the local solution is inside the toroidicity induced gap of the Alfvén continuum, as expected for an ordinary TAE [1]; however, the mode is weakly damped (γ < 0) in these regions. On the other hand, where the gradients of the fast particle and the background pressure are both strong near r/a = 0.6, the mode is unstable (γ > 0) and the local frequency ω drops into the lower continuum of an ideal MHD. Both the maximum of γ and a local minimum of ω occur in this strong gradient region, which is also in the region of weak magnetic shear. The observed mode is apparently destabilized in this region of maximum drive near r/a = 0.6. Since the frequency of the local solution changes significantly on a short radial scale in this region, we must consider the 2-D structure of the eigenfunction and the associated modifications to the mode frequency. An estimate of the global (2-D) frequency can be obtained from the Wentzel-Kramers-Brillouin (WKB) formalism [31,32]. We need to apply the following quantization condition: where θ k is the ballooning variable used in the ballooning mode theory to describe the mode radial envelope and is a parameter in the local eigenmode equation, k is the radial mode number and the integral is taken along the constant frequency of the local solution ω = ω(r, θ k ). The path of integration in Eq. (1) depends on the function ω(r, θ k ), and can have socalled open or closed trajectories in (r, θ k ) space. For the open trajectory case, the function ω(r, θ k ) has local minima; the 1-D result shown in Fig. 3 is of this type. Further calculations can be simplified by expanding the numerically obtained local frequency ω near its minima, which is also near the most unstable region, i.e. r 0 /a = 0.6: where the frequencies are normalized to the central Alfvén frequency ω A0 and 0 , X and Y are constants. The discharge under consideration has a q profile that, in the vicinity of its minima, can be approximated with good accuracy by the simple formula q = q min +8.1(r−r qmin ) 2 , where q min = 1.7 and r qmin = 0.53. For the experimentally observed n = 5, the local frequency is fitted to the parabolic dependence in minor radius and in θ k using the least squares method, which results in 0 = 0.412, X = 27.4 and Y = 0.0188. With these parameters, the global frequency can easily be derived from the quantization condition given by Eq. (1), which gives for the first radial mode (i.e. k = 1) 1 = 0.50 and for the second radial mode (k = 2) 2 = 0.58. Note that there is no third mode since the 'potential well' described by the dependence given by Eq. (2) is shallow in the direction of θ k , so that the second mode already has a frequency close to the edge of the well. The other result from the application of Eq. (1) is the mode location, which is limited by two turning points r 1 = 0.543 and r 2 = 0.655 for the mode with 1 = 0.5.
The predictions for other toroidal mode numbers for the RTAE (Fig. 3) are qualitatively similar to the ones shown in Fig. 3 for n = 5. The code predicts that the RTAE for the lowest toroidal mode number should be more unstable than those for higher mode numbers; the predicted frequency ω is also lower. (Note that global calculations may alter the stability prediction.) Inclusion of 2-D effects raises the predicted frequency, as it did for the n = 5 mode. For example, at n = 7, we obtain 0 = 0.47, 1 = 0.55, with localization between r 1 = 0.548 and r 2 = 0.65.
Next, consider the local frequency of the resonant KBM branch (the chain curve in Fig. 3(a)). For this branch, the frequency changes gradually with position. This means that the global mode does not form a closed trajectory in (r, θ k ) space. Since n 1, this means (Eq. (1)) that the KBM is more localized and that the local frequency is a good approximation to the true frequency. As shown in Fig. 3, both the frequency and the growth rate of the resonant KBM depend sensitively on beam pressure. In JT-60U, the frequency of this mode gradually rose into the TAE gap as the pressure of the plasma and fast particles increased [19]. This stronger variation of the frequency is due to the higher characteristic transit frequency of the negative neutral beam ions in JT-60U, which had an injection energy of ≈360 keV (≈80 keV is used in DIII-D). As with the RTAE, low toroidal mode numbers are more unstable than higher mode numbers (Fig. 3).

Comparison between theory and experiment
The theoretical predictions of the previous section are in reasonable agreement with the experimental observations. Both branches are predicted to be unstable, and the maximum growth rate occurs near the experimentally inferred radial location of the mode. Qualitatively, the 'TAE' corresponds to the RTAE, while the 'BAE' corresponds to the resonant KBM. The quantitative comparisons below include the local rotation frequency and are made at r/a 0.57.
If the experimentally observed 'TAE' corresponds to the second radial mode, the calculated frequency 2 /2π is ≈4 kHz higher than the experimental value of f pl , which is well within the experimental uncertainty of ≈9 kHz. If the measured 'TAE' actually corresponds to the first radial mode, then the calculated frequency is ≈12 kHz lower than the measured value.
The calculated frequency of the resonant KBM is ≈14 kHz higher than the measured frequency of the 'BAE'.
It is likely that the instabilities cause some redistribution of the beam population (although any redistribution must be modest or it would cause detectable reductions in neutron rate). A lower beam pressure lowers the resonant KBM frequency and raises the RTAE frequency. If the measured 'TAE' corresponds to the first radial mode, a reduction in beam beta would bring both calculated frequencies into agreement with experiment. Another important parameter is the magnetic shear, which may alter the predictions when slightly changed. We conclude that the predicted frequency is consistent with the observed frequency within experimental uncertainties.
For the RTAE, the code predicts that n = 4 should be more unstable than higher values of n. Experimentally, the largest mode in the 150-250 kHz band has n = 3-5 at each burst, with n = 4 being most common. For the resonant KBM, the n = 4 is predicted to be most unstable. Experimentally, n = 5 is the most unstable mode. Note that, since the HINST code does not calculate the 2-D eigenfunction, the growth rate prediction becomes increasingly unreliable as n → 0.
The theory predicts that the RTAE is more unstable than the resonant KBM. Experimentally, both the 'TAE' and the 'BAE' are observed, a situation that only occurs when the instability thresholds are comparable. From the reflectometer and Mirnov coil measurements, the 'BAE' appears more unstable than the 'TAE', but the relative sensitivity of the measurements depends on the eigenfunctions, which differ theoretically [19], so it is possible that the 'TAE' is as strong as the 'BAE' at r/a 0.6.
In the light of the quantitative success of the EPM theory in explaining this discharge, a qualitative re-examination of the extended database of Ref. [20] is in order. One perplexing feature documented in Fig. 21 of Ref. [20] is the variability of the mode frequency in successive bursts. The mode frequency sometimes jumps between low frequencies (f pl < 0.3f TAE ) and high frequencies (f pl > 0.6f TAE ) on successive bursts, when only slight changes in plasma properties occur. Also, the frequency variability within frequency bands is larger than the estimated experimental uncertainty. Such rapid changes in frequency and stability properties are difficult to explain if the instabilities are normal modes of the background plasma. (For example, changes in rotation frequency occur on a transport timescale and can only account for changes 0.05f TAE .) On the other hand, fast ion transport can rapidly alter the fast ion distribution, which has a strong effect on both the growth rate and frequency of energetic particle modes. Thus, the jumps between bands may occur when the most unstable mode alternates between the resonant KBM and the RTAE, and the variability within a band may be associated with changes in the fast ion distribution function.
Reference [20] dismisses the KBM [12] as a plausible explanation for the 'BAE' because the data do not scale with the ion diamagnetic frequency ω * pi , because the temporal evolution of the mode frequency is often opposite to ω * pi and because instabilities occur in discharges that are far from the MHD ballooning mode stability boundary. None of these objections apply to the resonant KBM branch. Theoretically, the resonant KBM can occur even in plasmas with ω * pi = 0 [16]. Because the mode frequency and growth rate depend strongly on the fast ion population, the theoretical frequency differs significantly from ω * pi and instability can occur well below the ideal MHD threshold.
The original observation of the 'BAE' reported that f pl ∝ v A [4], which supported identification of the instability as an Alfvén eigenmode; however, the dependence on v A in the extended database is very weak [20]. Theoretical work by Zonca et al. [33] suggests the following relationship between the KBM and the BAE. The original theory of the KBM [12] assumed that the ion transit frequency ω ti = v th,i /qR is small compared with the ion diamagnetic frequency ω * pi . In this limit, the low frequency accumulation point of the Alfvén continuum occurs at ω * pi . The opposite limit, ω ti ω * pi , was considered in the original theory of the BAE [11,34], where the solution was shown to be localized inside a low frequency beta induced gap of the Alfvén continuum that exists due to core plasma compressibility. In this limit, the low frequency accumulation point of the Alfvén continuum occurs at a frequency of ∼ √ γ ω ti . In reality, the actual low frequency shear Alfvén wave is neither a pure KBM nor a pure BAE in plasmas with ω ti ∼ ω * pi . The original observation of a 'BAE' occurred in plasmas with ω ti 2ω * pi , so the pure BAE theory may be applicable. The hotter plasmas observed subsequently are closer to the conditions assumed in pure KBM theory. In HINST, ω ti is included in the nonadiabatic plasma response, which gives a Landau type damping but which cannot reproduce a separate BAE branch. For the plasma analysed here, ω ti 0.2ω * pi , thus in the HINST code the KBM approximation that the low frequency gap is determined by ω * pi alone is reasonable.
The parameters with the strongest effect on mode stability are the beam pressure and the angle of beam injection, with there being little dependence on other plasma parameters [20]. Theoretically, EPM stability is very sensitive to the details of the fast ion distribution and pressure.

Conclusion
Local analysis of a DIII-D discharge with 'BAE' and 'TAE' activity finds unstable EPMs with frequencies and mode numbers comparable to those found in experimental observations. This lends credence to the notions that the 'BAE' is a resonant KBM and that the 'TAE' is an RTAE.
In future work, more stringent tests of the theory will be possible. To verify the importance of the energetic particle population, the predictions of HINST can be compared with measured variations in frequency and growth rate associated with experimental modification of the beam ion distribution function. With a global non-perturbative code, more rigorous predictions of the most unstable n modes are possible. The predicted 2-D structure should also be compared with internal measurements of the eigenfunction.