Non-linear dynamics of toroidicity-induced Alfvén eigenmodes on the National Spherical Torus Experiment

The National Spherical Torus Experiment (NSTX, (Ono et al 2000 Nucl. Fusion 40 557)) routinely operates with neutral beam injection as the primary system for heating and current drive. The resulting fast ion population is super-Alfvénic, with velocities 1 < vfast/vAlfven < 5. This provides a strong drive for toroidicity-induced Alfvén eigenmodes (TAEs). As the discharge evolves, the fast ion population builds up and TAEs exhibit increasing bursts in amplitude and down-chirps in frequency, which eventually lead to a so-called TAE avalanche. Avalanches cause large (≲30%) fast ion losses over ∼1 ms, as inferred from the neutron rate. The increased fast ion losses correlate with a stronger activity in the TAE band. In addition, it is shown that a n = 1 mode with frequency well below the TAE gap appears in the Fourier spectrum of magnetic fluctuations as a result of non-linear mode coupling between TAEs during avalanche events. The non-linear coupling between modes, which leads to enhanced fast ion transport during avalanches, is investigated.


Introduction
The interaction between fast ions and multiple toroidicityinduced Alfvén eigenmodes (TAEs [1]) is believed to represent one of the main loss mechanisms for fast ions in ITER [2]. Major consequences of enhanced losses are a decreased fusion efficiency and possible damage to in-vessel structures. It is therefore important to understand this phenomenon in present devices in order to limit, or possibly avoid, its deleterious effects in future fusion reactors [3]. A scenario of particular concern for ITER contemplates the coupling of multiple TAE resonances resulting from their overlap in phase space [4]. If the coupling is sufficiently strong, broad regions of the fast ion distribution function provide enough free energy to sustain an explosive growth of the modes [5]. In general, the modification of the fast ion phase space, which can be further enhanced by the resonance overlap, has been invoked to explain a variety of phenomena observed in TAE experiments. These include mode splitting into multiple frequencies [6], or periodic variations of the mode frequency (chirps) [7,8]. As a result of such non-linear TAE dynamics, an enhancement of fast ion loss or redistribution (in both velocity and real space) is predicted. Whether fast ions are actually expelled or not from the core plasma depends in part on the radial extension of the mode structure. Significant losses are expected for global, radially extended TAEs with low toroidal mode number, whereas more localized TAEs with a larger toroidal mode number are more likely to induce redistribution. The latter case is contemplated for the high-field, large major radius plasmas of ITER [9]. It should be noted that, although it is not directly affecting the overall confinement of fast ions, redistribution itself can be deleterious in a future reactor. First, modifications of the fast ion distribution may result in unwanted variations of the non-inductive current profile [10]. Second, a relaxed fast ion pressure can impact the stability of other modes, such as kink modes [11], sawteeth [12] and resistive wall modes [13], thus affecting the overall plasma stability.
The most challenging issue related to fast ion transport for future devices is arguably represented by events that happen on very short time scales, 1 ms, for which no possibility of external control (e.g. through actuators such as (NBs) or rf waves) currently exists. Although the typical scenarios on the National Spherical Torus Experiment (NSTX) [14] are considerably different from those expected in ITER and future fusion devices (cf sections 2 and 3), they provide a good test case to improve the present understanding of the basic physics of bursting TAEs and induced fast ion loss or redistribution. Bursts of TAE activity with duration ∼1 ms are commonly observed in NB-heated plasmas of NSTX. A prompt depletion of up to 30% of the confined fast ion population correlates with single events [15][16][17]. The TAE bursts are thought to be associated with an enhanced coupling between unstable modes. The experimental results discussed in this paper are consistent with the occurrence of a quadratic coupling between TAEs with consecutive toroidal numbers during their bursty, explosive growth.
The paper is organized as follows. The experimental scenario used for TAE studies on NSTX is summarized in section 2. The general observations on TAE dynamics are presented in section 3. The evidence of coupling between TAEs, as well as the conditions leading to coupling, is discussed in section 4. Section 5 summarizes the main finding of this work and concludes the paper.

Experimental scenario for TAE studies on NSTX
NSTX operates at toroidal field 3.5-5.5 kG, with typical density (3-10)×10 19 m −3 and temperature T e ≈ T i 1.5 keV. NB injection is the primary system for heating and current drive, with a maximum available power P NB = 7 MW. The range of injection energies is E inj = 60-90 keV. The resulting fast ion population is super-Alfvénic with velocities 1 < v fast /v Alfven < 5. Fast ions provide the drive for a variety of Alfvénic instabilities, including TAEs [15]. The latter have a typical toroidal mode number of up to n = 8 and frequency 60 < f < 250 kHz.
An example of the plasma scenario used for TAE studies on NSTX is shown in figure 1 for a L-mode deuterium plasma with 1.5 MW of injected NB power. The reversed-shear safety factor profile, q(R), evolves in time. Its minimum, q min , decreases from 4 to ≈2 during the time of interest ( figure 1(b)). Since no direct measurement of q(R) through the motional Stark effect (MSE) is available for this discharge, the safety factor is reconstructed a posteriori from similar discharges where MSE was available through the Grad-Shafranov equilibrium code LRDFIT [18].
The spectrum of magnetic fluctuations, measured at the low-field side of the vacuum vessel by Mirnov coils, is shown in figure 2(a). Quasi-periodic variations of the mode amplitude correlate with frequency variations <10 kHz for t < 310 ms. Later in time, the modes exhibit larger bursts in amplitude and frequency down-chirps >10 kHz, which eventually lead to a so-called TAE avalanche [17]. Avalanches cause fast ion losses of up to ∼30% over ∼1 ms [16], as inferred from the neutron rate, cf figure 2(b). (The neutron rate being mostly determined by beam-plasma reactions, the relative drop in neutron rate is approximately equal to the fractional reduction in the fast ion population). The drop in the neutron rate as a function of the mode amplitude, integrated over 15 < f < 250 kHz, is illustrated in figure 2(c) for a set of similar discharges. As expected, the losses increase with the mode amplitude [17,19]. In addition, activity with n = 1 and f ∼ 25 kHz (i.e. well below the TAE gap), along with a weaker n = 2 harmonic at twice that frequency, is detected during strong TAE bursts. This fluctuation has a duration comparable to that of the avalanches. The temporal correlation with TAE bursts suggests that the n = 1 fluctuation may have a role in the fast ion loss process, as discussed in the next sections.

General characteristics of TAE dynamics
For the type of discharges investigated in this work, TAEs exhibit common, general features that do not vary substantially from shot to shot. A first operational distinction can be made between the behaviour of each single mode and that of all the unstable TAEs observed at a given time as an ensemble. A typical result for the quasi-periodic frequency and amplitude variations of single TAEs in the bursting/chirping phase is shown in figure 3(a). During each burst the mode amplitude grows by an order of magnitude, with a relative change in frequency 10%, then decays to the noise level. A new burst follows shortly after, with the initial frequency sweeping down as the amplitude reaches its maximum and then decays. The cycle is repeated on time scales 1 ms. In fact, the relationship between amplitude and frequency during the bursts (figure 3(b)) is reminiscent of a limit cycle for nonlinear harmonic oscillators, such as the forced Van der Pol oscillator [20]. This observation may provide the basis for the interpretation of the bursting/chirping single-mode TAE regime observed experimentally, in parallel with other models developed in recent years (cf [21] and references therein).
Other salient features involve the entire set of TAEs, for instance the almost constant frequency separation between peaks with consecutive ns. This suggests that the modes share a common frequency in the plasma frame, f TAE 0 , such that [22,23] f TAE lab,n = f TAE 0 + n f TAE

Doppler
(1) where f TAE lab,n is the mode frequency for the toroidal mode number n as measured in the laboratory frame and f TAE Doppler the Doppler shift caused by plasma rotation. In practice, a Gaussian fit of the peaks in the fast Fourier transform (FFT) spectra from 1.25 ms time windows is first used to determine the centre frequency and frequency spread of the different modes. Then, a linear fit of f TAE lab,n as a function of n, based on equation (1)   is also observed in f TAE 0 . The jump is followed by a qualitative transition in the mode dynamics from quasi-stationary to a bursting/chirping regime. As a speculation, we observe that the mode structure and its (poloidal) harmonic composition may be changing, for example because of the evolution of the q-profile, and the different mixture of the poloidal harmonics of each mode favours its bursting character. In addition, the fast ion population also increases, as deduced by the rise in the neutron rate (figure 2(b)), indicating that the mode's drive may also be different for the two phases. Although the change in regime is observed in all discharges, not all of them show this sudden frequency variation.
As a further step, the radius R TAE at which the plasma rotation f rot (from charge-exchange recombination spectroscopy measurements of carbon rotation [24]) matches f TAE Doppler is calculated and compared with other quantities. The latter include the location of the steepest fast ion density gradient, obtained through fast ion D-alpha measurements [25,26] of the fast ion profile, and of the maximum rotation shear (figure 5). By using the condition f rot (R TAE ) = f TAE Doppler , the diamagnetic contribution to f TAE Doppler [22] is 1 kHz and is neglected. Also, the differential rotation between impurity and main ions [27] is neglected and the measured f rot is used as rotation of the main ion species. The radii R TAE and of the steepest fast ion gradient are approximately equal and show a similar temporal evolution [28]. The same position also corresponds to the location of largest rotation shear. Although the high shear does not affect the macroscopic mode dynamics (see [28]), it implies that small changes in R TAE may lead to relatively large variations of f TAE lab,n due to different Doppler shifts experienced by the modes, see below. The coincidence of the spatial locations of the modes' Doppler shift and drive and of the maximum rotation shear is presumably caused by the dependence of the fast ion and rotation profiles on the same source, namely the NB deposition profile. This suggests that NB injection, as opposed to other heating schemes, is intrinsically favouring the destabilization of multiple modes at similar radii by sustaining a steep fast ion gradient in a limited spatial region, thus paving the way for mutual interactions between modes. The implication for the projection of these results to future reactors, mostly envisaging rf waves as heating and current-drive tools (e.g. ITER), is twofold. First, the occurrence (and nature) of non-linear coupling between TAEs destabilized by alpha particles may differ considerably from what is observed in those present devices in which fast ions originate from NB injection. Second, reliable predictions must be based upon a robust knowledge of the evolution of the fast ion population, including its source terms. The small error bars in the fit of f TAE 0 and f TAE Doppler indicate that equation (1) describes well the frequency evolution over time scales 1 ms. The mechanism through which TAEs with different ns lock on a similar frequency is not yet clear. At present, it can be excluded that the modes are harmonics of a fundamental TAE mode or that they couple through a persistent (non-bursting) low-frequency mode, as was observed in [22]. The sub-millisecond dynamic is more complex. Each mode usually chirps independent of the others. The largest bursts represent an exception, in that most of the modes seem to lock onto a similar dynamic. For short time windows, the procedure illustrated above to find f TAE 0 and f TAE Doppler gives more scattered results ( figure 4). Nevertheless, two features emerge for the largest bursts: (i) f TAE Doppler increases by a few kilohertzs, implying that R TAE moves inwards (up-hill along the toroidal rotation profile); (ii) the frequency down-chirp in the laboratory frame is almost entirely ascribed to a downchirp of f TAE 0 , see figure 6. A more complete description of the TAE dynamic during a burst requires a different analysis, as explained in the following section.

Verification of the conditions for mode coupling
The presence of mode-mode interactions is first investigated through bicoherence [29] which quantifies the degree of correlation between triplets of modes via a phase-weighted average of the complex Fourier bispectrum, normalized to the spectral amplitude. The bicoherence calculated before a TAE burst (figure 7(a)) is vanishing for small amplitude oscillations, indicating a weak or null coupling between modes which, nevertheless, still exhibit a bursting/chirping character. Conversely, values of up to 75% are obtained during the burst for frequencies corresponding to the TAEs ( figure 7(b)) and to the modes at f < 60 kHz and f > 160 kHz. Because of its limitations in temporal resolution and the lack of statistics for events of short duration, Fourier analysis has limited use to study the evolution of the mode coupling process during single bursts (hundreds of microseconds). An analysis in the time domain is instead used here. The analysis is performed as follows. Signals are band-pass filtered around the mode frequencies to obtain the magnetic fluctuations,ṡ n (t) for each mode number n = 1, 2, . . .. The sum and difference of signals from sensors toroidally displaced 180 • apart are used to further separate modes with even or odd mode number. The B-field fluctuation signals, s n (t), are calculated by integratingṡ n (t) via software. The evolution of frequency and amplitude of the modes, f n and A n , are then obtained from the peak-to-peak amplitude and maxima separation of the resulting (sinusoidal) signals s n , each representing a mode with a specific frequency and toroidal mode number. Figure 8 shows the results of such analysis for a TAE burst at t ∼ 328.2 ms from the same discharge shown in figure 2(a). Before the burst only the primary TAEs are detected, with dominant n = 2, 3, 4. The frequency of the modes varies in time, but with no obvious correlation up to ≈328 ms. Then, their amplitude increases ( figure 8(c)) and other modes are clearly emerging at higher and lower frequencies, filling up the frequency spectrum up to 250 kHz. A progression in the toroidal mode number with frequency is observed within each group of modes (primaries, low and high frequency), with a roughly constant separation between subsequent n's ( figure 8(b)). The latter is consistent with the frequency of the lowest frequency n = 1 mode. In general, the different temporal evolution of the mode amplitudes and the measured frequency/wavenumber lead to exclude that the high-frequency modes are simple harmonics of the primary TAEs. From a simple model based on bilinear interactions between pairs of modes, the resulting fluctuation can be reconstructed asṡ with the constraints between mode numbers and real frequencies given by The right-hand side in equation (2) is filtered around the frequency f n3 and the dot indicates a time derivative. Note that s n2 is the complex conjugate of the true signal for difference interactions. c(n 1 , n 2 ) represents a coupling coefficient. For finite damping or growth of the modes, c(n 1 , n 2 ) is a complex quantity [30] whose phase, along with the relative phase between the modes involved in the coupling, determines whether an explosive instability can occur [31]. For the measured signals, s n → {s n } ( {. . .} is the real part). In practice, all phase terms, including the relative phase between modes and the phase of the coupling coefficient, will result in a phase, φ, between left-and right-hand sides of equation (2). In general, it can be expected that random, rapid changes in the phase destroy the phase coherence between modes and inhibit an efficient coupling. According to a simplified model with negligible net damping or growth (i.e. a real coupling coefficient), unbounded solutions for the modes' amplitude are found for φ = π (φ = 0) [30], corresponding to ideal difference (sum) interactions between modes. In this case, values φ = 0, π will therefore indicate a mismatch between either frequencies or wavenumbers. Although some of the experimental results are consistent with this simplified model, see below, its applicability still needs to be verified in detail. In the following, only the more general condition of temporal phase coherence is retained as necessary condition for a significant mode coupling.
The measured and reconstructedṡ n=1 (t) signals, shown in figures 9(a) and (b), are consistent with this model for the pairs n = (4, 3) and n = (3, 2) and a phase shift of 180 • . The n = 1 oscillation amplitude grows at the same rate as the (4,3) and (3,2) bilinear terms, then decays when the amplitude of (at least one of) the pump TAEs vanishes, or the frequency match condition is no longer satisfied. A similar temporal evolution, but with a phase shift φ ∼ 90 • , is observed for fluctuations resulting from sum interactions between the dominant modes. The roughly exponential decay of the n = 1 fluctuation is characterized by a time constant ≈175 µs (see the inset in figure 9(b)), which corresponds to an effective damping rate of 4-5% in the absence of the pump primary modes. The frequency f 1 overlaps with the difference f 4 − f 3 and, later in time, with f 3 − f 2 (figure 9(c)). When that happens, the mode amplitude first increases rapidly, then decays again as the condition is no longer satisfied. The frequency match for mode-mode coupling is therefore transiently verified, as shown in figure 9(c) for the reconstructed n = 1 amplitude according to equation (2). Figure 9(d) shows the displacement measured through reflectometry [32] for the dominant modes at the time of maximum amplitude. Within the uncertainties due to the reduced spatial coverage, the modes peak at approximately the same location, as expected for efficient coupling. It should be emphasized that these results from equation (2), obtained using the frequencies measured in the laboratory frame, are valid only if equation (1) is verified. If each mode experiences a different Doppler shift, the mode frequencies in the plasma frame, inferred from the actual mode locations and individual Doppler shift frequencies, must be used instead.
Additional information on the the spatial relationship between mode amplitudes is available from the entire set of 11 Mirnov coils distributed toroidally over 360 • . By fitting the signal from the coils against the toroidal angle, the toroidal mode structures can be qualitatively visualized along with their relative (toroidal) phase at each time. The mode amplitude is represented in a polar plot as the deviation from the unitary radius circle (the dotted line in figure 10) of the fit as a function of the toroidal angle for each specific mode number. All mode amplitudes are rescaled by the same factor. The example given in figure 8(d) and detailed in the polar plots in figure 10 clearly shows that the n = 1 perturbation is generated by the overlap of positive and negative lobes of the n = 4, 3 modes, producing regions of constructive and destructive interference between mode pairs with a definite n = 1 structure. For this mechanism to be efficient, three-wave matching conditions must be maintained for a sufficiently long time, of the order of tens of wave cycles of the primary (pump) modes, as observed from the experiments.

Discussion
In this work it has been shown that pairs of TAEs can efficiently couple and generate lower and higher frequency perturbations. The matching conditions for mode-mode coupling are verified for frequency and wave-number of the primary and resulting modes. The coupling is inherently transient because of the different single-mode dynamics of the primary TAEs. As a consequence, random cross-phase variations between the primary modes on short time scales tend to invalidate the coupling conditions. However, when those conditions are met for hundreds of microseconds, i.e. for several tens of wave cycles of the primary modes, coupling can take place. When that happens the modes undergo a fast growth, which eventually lead to fast ion losses of O(10%) of the confined population. During this process the TAE amplitude grows by more than one order of magnitude. The effective rate of amplitude raise γ eff /2πf , estimated experimentally from the nearly exponential increase of the mode amplitude, is 10%, i.e. much larger than the typical linear growth rate of O(1%) estimated through the NOVA-K code [33] for TAEs on NSTX. The coupling between TAE modes is favoured by their accumulation around the same radius, which is plausibly caused by the steep fast ion density gradient measured at that location (cf figure 5). Similar conditions can be expected in fusion devices such as ITER. However, the higher n numbers predicted to be unstable in ITER [9] will have a narrower radial structure, which imposes more stringent conditions in terms of spatial overlap for an efficient coupling between modes. Nevertheless, even a more sporadic coupling between TAEs in ITER should be a concern for its uncontrollable effects on the fast ion population and associated quantities such as the NB current-drive profile.
The strongest mode resulting from bilinear interactions is a n = 1 oscillation, which mediates the non-linear coupling between the primary TAEs. This mechanism appears to be different from other mode coupling processes previously observed on NSTX [34], where a set of harmonics identified as fishbone modes concentrated the TAE modes into a slowly rotating wave packet. Unlike the case studied in this work, the low-frequency modes showed strong frequency chirps, characteristic of fishbone modes, which repeated themselves in pulses spaced by ≈5 ms. For the scenario discussed in this paper, the nature of the n = 1 (and the weaker n = 2) lowfrequency fluctuation remains an open question. Although the signal-to-noise of the reflectometer data is low (figure 9(d)), its structure peaks at the same location of that of the TAEs, around ≈132 cm. No data are available closer to the magnetic axis, R 120 cm. As a consequence, see equation (1), the mode frequency in the plasma frame is f n=1 Doppler . The more plausible interpretation is that f n=1 Doppler ≈ f TAE Doppler , i.e. the perturbation rotates rigidly with the TAE revolution frequency and is localized at similar radii as the primary TAEs. The spatial overlap between modes favours the coupling, and is in fact measured (cf figure 9(d)) through the reflectometer system.
To investigate the generality of these results for the n = 1 fluctuation, the analysis is extended to a broader data set consisting of data from five L-mode discharges with TAE activity. The presence of n = 1 activity during the larger TAE bursts is confirmed for these discharges. In contrast to the intuition, no clear dependence of the n = 1 amplitude upon the amplitude of the primary TAEs is observed. However, this can be understood by considering that multiple conditions, e.g. on the exact frequency and phase matching of the TAEs, are required for an efficient coupling. The relation between f TAE Doppler and the observed n = 1 frequency for the extended dataset is shown in figure 11. The results are based on the spatio-temporal analysis of the signals from the eleven Mirnov coils (see, for instance, figures 10(c) and (d). The signals from a moving time window of 1.2 ms are filtered in the frequency band f TAE Doppler ± 10 kHz, then a fit of the toroidal mode structure imposing n = 1 for all the experimental time points provides mode amplitude and frequency as a function of time. Three main conclusions can be drawn. First, a substantial correlation between the frequency of n = 1 activity and f TAE Doppler is observed when the n = 1 amplitude is above the noise level. This result is independent of the details of the band-pass filter used to isolate the n = 1 activity, as long as no other modes are present in the selected frequency range. Second, the spread in the data suggests that f n=1 lab in the rotating plasma frame is usually close, but not necessarily identical, to the TAE Doppler shift frequency. Statistically, it is found that f n=1 lab f TAE Doppler , i.e. f n=1 0 0. Third, the n = 1 amplitude is, at least qualitatively, larger when f n=1 lab ≈ f TAE Doppler . Based on this information, two possibilities are left open: (i) the n = 1 perturbation is a forced oscillation (or a quasi-mode) driven by the large amplitude primary TAEs, or (ii) it is an otherwise stable, kinklike plasma mode driven unstable by the coupling process and by the associated redistribution or loss of energetic ions. In this latter case, a significant coupling with the TAE modes can be expected only for those n = 1 eigenmodes which have a radial structure similar to that of the TAEs. Further experiments  with improved mode structure diagnostic capabilities (e.g. with enhanced spatial coverage of the reflectometer system and through beam emission spectroscopy [35]) and a direct comparison with self-consistent, non-linear codes such as M3D-K [36] are foreseen to elucidate this issue.