Saturation of single toroidal number Alfvén modes

The results of numerical simulations are presented to illustrate the saturation mechanism of a single toroidal number Alfvén mode, driven unstable, in a tokamak plasma, by the resonant interaction with energetic ions. The effects of equilibrium geometry non-uniformities and finite mode radial width on the wave-particle nonlinear dynamics are discussed. Saturation occurs as the fast-ion density flattening produced by the radial flux associated to the resonant particles captured in the potential well of the Alfvén wave extends over the whole region where mode-particle power exchange can take place. The occurrence of two different saturation regimes is shown. In the first regime, dubbed resonance detuning, that region is limited by the resonance radial width (that is, the width of the region where the fast-ion resonance frequency matches the mode frequency). In the second regime, called radial decoupling, the power exchange region is limited by the mode radial width. In the former regime, the mode saturation amplitude scales quadratically with the growth rate; in the latter, it scales linearly. The occurrence of one or the other regime can be predicted on the basis of linear dynamics: in particular, the radial profile of the fast-ion resonance frequency and the mode structure. Here, we discuss how such properties can depend on the considered toroidal number and compare simulation results with the predictions obtained from a simplified nonlinear pendulum model.


Introduction
Plasma physics finds one of its most relevant applications in the field of controlled thermonuclear fusion research. A promising approach to the realisation of a fusion reactor is the so-called magnetic confinement, in which a hot plasma (temperatures of the order of those characterising stars' core) is contained by means of strong magnetic fields. Tokamak is the most advanced magnetic-confinement concept. In order to bring plasma to conditions in which fusion reactions can occur at a sufficiently large rate, several heating methods (based on radio-frequency waves and neutral beam injection) have to be added to the ohmic heating produced by the electric current flowing through the plasma itself. These methods, as well as fusion reactions, produce energetic particles (EPs) with velocities of the order of the Alfvén speed, the typical propagation velocity of shear Alfvén waves (SAWs). Such particles can then easily resonate with SAWs, with their pressure gradient acting as a freeenergy source for destabilising these waves. On their turn, Alfvénic fluctuations can cause an enhanced transport of EPs, affecting their confinement, preventing their thermalisation in the central region of the plasma and increasing thermal and particle loads on the tokamak wall. The assessment of EP confinement properties in next generation fusion experiments then strongly relies on understanding Alfvén mode dynamics, with regard both to linear stability properties (which modes are expected to be driven unstable) and the nonlinear saturation mechanisms (which saturation level is expected for the mode amplitude and which effects on the EP confinement).
In nonuniform tokamak plasmas, the SAW spectrum becomes a continuum. Poloidal asymmetries open 'gaps' in such continuum, allowing for weakly damped global modes, such as the toroidal Alfvén eigenmode [1] or the reversed shear Alfvén eigenmode [2], to exist. The resonant interaction with EPs can drive these modes 2. Overview of single mode problem Nonlinear wave-particle interactions for a single wave in one-dimensional (1D) systems was originally investigated by O'Neil and co-workers in 1960s [31][32][33]. In 1990s, the beam-plasma nonlinear problem was reconsidered by Berk and Breizman [23][24][25][37][38][39] and applied to the interpretation of the dynamics of the Alfvén eigenmodes driven by EPs near marginal stability. Such 'bump-on-tail' paradigm treated EP contribution as a perturbative one, i.e. the structure of Alfvénic fluctuations is determined by MHD equations (no kinetic contributions), and its width is much larger than the nonlinear radial excursion of resonant EP orbits. Recently, a 'fishbone' paradigm [22,27] has been introduced, which consider the response of EPs in the frame of a nonperturbative approach, in which it can modify the plasma dielectric response as well as mode structure and frequency. In this paradigm, equilibrium geometry and plasma non-uniformity play important roles and become more and more relevant with increasing drive. In the present work, we will show how these effects affect the nonlinear dynamics of the modes and their saturation amplitude.

1D collisionless beam-plasma system
The original works on nonlinear wave-particle dynamics in 1D systems [32] studied the behaviour of a suprathermal electron beam interacting with a thermal plasma in a strong axial magnetic field. Various processes such as Landau damping in a finite amplitude wave [31] and nonlinear wave-particle interactions [33] were understood for the first time within the 1D beam-plasma model.
In [32] a collisionless uniform 1D beam-plasma system is considered, with Maxwellian bulk electrons characterised by density n e and thermal velocity v th , much lower than the electron-beam drifting speed v d . The beam distribution function is assumed to be Lorentzian, with density n n b e  and velocity spread Dv b . The resulting total electron distribution is known as 'bump-on-tail' distribution. Thermal ions are treated as a fixed neutralising background. In such system, a plasma wave is driven unstable. The wave grows and saturates at a scalar-potential amplitude of order f~( ) n n mv e; b e 2 3 d 2 after that, it starts oscillating [33]. The reason for this is that, as the wave amplitude grows, more and more resonant electrons slosh back and forth in the potential well of the wave. The flux associated to such wave-trapped electron motion yields a temporary flattening of the velocity-space electron distribution function (with consequent mode saturation), followed by a reconstruction of the distribution-function gradient, a new flattening, and so on (with the corresponding mode amplitude oscillation). The oscillation frequency is related to a mean value of w B , where w B is the bounce frequency of resonant electrons in the potential well of the wave.
On longer time scales, the spread in w B values produces increasingly finer structures of the electron distribution function around the resonance. The net mode-particle energy exchange reduces, along with the mode amplitude oscillations. By this phenomenon of phase mixing, the system asymptotically reaches a steady state [31], characterised by a flattened coarse-grain particle distribution function [40] and a constant mode amplitude.

Berk-Breizman model: bump-on-tail paradigm
In a series of papers in 1990 [23][24][25], Berk and Breizman reconsidered the beam plasma bump-on-tail nonlinear problem and applied it to Alfvén Eigenmode dynamics. Their work adds, to the previous analysis by O'Neil [31,41], the treatment of particle interactions with a finite-amplitude wave in the presence of sources and collisions. Because of the associated dissipation, the system reaches saturation while the distribution function still presents a finite gradient, such that the residual drive balances the background dissipation. Berk and Breizman's analysis relies on several assumptions: (1) mode amplitude is low and the linear mode structure is fixed, but sampled by the perturbed particle orbits; (2) finite background dissipation is independent of the wave amplitude; (3) wave dispersiveness is set by the background plasma and is not affected by the beam (in bump-on-tail problem) or the EPs (in general tokamak plasma system). Under these assumptions, mode saturation level comes out to be reduced with respect to the dissipationless system [23,24].
Furthermore, the formation of propagating phase-space holes and clumps can occur, giving rise to adiabatic frequency chirping: the mode adjusts it in order to adapt to the resonance frequency characterising the instantaneous hole/clump phase space localisation [37][38][39]. The 'bump-on-tail' paradigm has been successful in explaining several phenomena observed in experiments [42][43][44][45][46].

Resonance detuning and radial decoupling
The importance of plasma non-uniformity has been pointed out by [22,27], in the frame of the so-called 'fishbone' paradigm. The basic point is that, for increasing drive and, hence, mode amplitude, particle orbits become able to explore the finite radial-width of the fields. At the same time, the resonance condition is satisfied, in this limit, in a wider radial region (whose width scales with the mode growth rate). The spatial region where the mode-particle power exchange can take place becomes than limited by the mode width rather than the resonance width. The consequences of this fact on the nonlinear dynamics and the saturation process of Alfvén modes in tokamaks have been investigated numerically [26] by evolving a suited set of test particles in the fields computed by self-consistent single toroidal mode number simulations. Provided that the mode frequency remains constant during the nonlinear evolution (condition well satisfied in the simulations considered in [26]), the system exhibits two invariants of the perturbed motion: besides the magnetic moment M, the quantity w = f C P nE, with ω being the mode frequency, f P the toroidal angular momentum, n the toroidal number and E the particle kinetic energy. The invariance of M and C implies that gradients of the distribution function along those directions do not play any role in the dynamics (they do not appear in the Vlasov equation). Then, cutting the phase space into slices orthogonal to the axes M and C, mode-particle power exchange can be described as the combination of different-slice contributions, each of them evolving in an independent way. In order to explore in great detail resonant-particle behaviour, test particles are selected by fixing the values of = M M 0 and = C C 0 in correspondence with a peak, in the phase space, of the linear-phase mode-particle power transfer. For each test particle, the information related to its evolution is collected at the equatorial-plane crossing times. Particles are referred to the coordinates Q f ( ) P , , with Θ being the wave-particle phase (e.g. see figure 1). During the whole linear phase, in which the field amplitude is so small that its effects on the particle orbits is negligible, f P is almost conserved. Then, during such phase, particle trajectories in the Q f ( ) P , plane essentially reduce to fixed points for = f f P P res , with f P res being the value at which the resonance condition for the considered set of particles, (with w res being the resonance frequency), is satisfied, while they correspond to drift along the Θ axis in the positive/negative direction, for f P greater/less than f P res . In the nonlinear phase, f P varies because of the mode-particle interaction (e.g., radialÉ B drift). Even particles that were initially resonant are brought out of resonance, getting non zero Q t d d and drifting in phase until the drift in f P is inverted. Particles that cross the = f f P P res line revert values of Q t d d as well. Thus, their orbits are bounded and they would properly close if the field amplitude were constant in time. This is true for particles born close to the resonance, while particles born with f P far from the resonance maintain drifting orbits, as they do not cross = f f P P res . Bounded orbits and drifting ones are divided by an instantaneous separatrix, whose f P width increases with increasing field amplitude. The formation of such structures in the Q f ( ) P , plane corresponds to a particle flux from the higher density side of the resonance to the lower density one, and a corresponding flattening of the resonant-particle density (as well as pressure) profile. The drive is consequently reduced. Saturation occurs when the flattening region extends by the whole region in which the mode-particle power exchange can occur. Such region will be limited both by the finite mode width and the finite resonance width. This is represented in figure 2, where the radial profiles have been considered, instead of f P profiles, taking into account that for given ( ) M C , 0 0 there is a one-to-one correspondence between f P and the (equatorial-plane) radial coordinate r: res L it then increases with the linear growth rate g L . Two opposite regimes are represented in figure 2 (left) and (right), respectively. In the former ('resonance detuning' regime), the most stringent constraint on power exchange is set by the resonance width; in the latter ('radial decoupling' regime), it is set by the mode width. The transition from the resonance detuning regime to the radial decoupling one occurs, ceteris paribus, as larger drive, narrower modes and/or flatter resonance frequency profiles are considered. In [26], examples of resonance detuning as well as radial decoupling have been shown, on the basis of numerical simulations performed by the HMGC code [35].
In the next section, we will show that each of the two regimes corresponds to a specific scaling of the saturation amplitude with the linear growth rate. This correspondence has been explained in [47] and [48] on the basis of a simple nonlinear pendulum model. An approximate analytical solution of that model connects the radial width of the resonant particle density flattening with the instantaneous mode amplitude and the linear growth rate. We will see how such solution allows us to predict the dependence of the saturation amplitude observed in the different regimes.
, at certain times in the linear (left) and nonlinear (right) phase. The position of each marker corresponds to the last equatorial-plane crossing of a test particle. Colours (blue or green) depend on the initial P f value of each test particle and do not change with time. The boundary between two colours is initially fixed at = f f P P res . The formation, in the nonlinear phase, of bounded orbits separated from drifting orbits by an instantaneous separatrix can be seen.

Simulations of bump-on-tail instability
We perform a numerical simulation of the bump-on-tail system by the df Vlasov-Poisson code PIC1D-PETSc [49]. The distribution function for a homogeneous electron beam-plasma system can be written as: with v th being the electron thermal velocity, v d the beam velocity and Dv b the beam thermal spread. The parameters and v d =5.0 are fixed. The length of the system is p = L k 2 with k=0.3. Periodic boundary conditions are used. We initially perturb the distribution function by a space dependent amount with  1  . In our simulations, time t is normalised to the inverse electron plasma frequency w - By increasing the beam density n b , the growth rate increases, as well as the saturation level, as shown in figure 3(a). The mode reaches a first saturation, and, after that, amplitude oscillations occur. Such oscillations are due to the particles trapped in the wave bouncing back and forth, and it has shorter period for larger saturation amplitude. Figure 3(b) shows that the saturation amplitude scales quadratically with the linear growth  rate g L . Such scaling is also important for our later discussion concerning the two different saturation mechanisms: namely, resonance detuning and radial decoupling.

XHMGC model
In the present paper, the extended version [36] of the nonlinear MHD-Gyrokinetic code HMGC [15,35] (XHMGC) has been used to simulate self-consistenly the BAE mode driven by anisotropic Maxwellian fast ions. In XHMGC, a shifted circular magnetic flux surfaces, low-β tokamak equilibrium is adopted. The fluid response of the thermal background plasma is described by a set of ( )  O 3 -reduced MHD equations [50] (with  being the inverse aspect ratio). Fast-ion and thermal-ion kinetic dynamics enter suche equations via the respective pressure tensors, which are computed by solving the Vlasov equation for each species in the drift-kinetic limit, by particle-in-cell techniques. Finite-Larmor-radius effects are ignored, while finite-orbit widths are taken into account in order to retain resonant wave-particle dynamics associated with guiding-center magnetic-curvature drift in toroidal geometry. Kinetic contributions are treated in a non perturbative way: pressure-tensor terms contribute to determine both structure and evolution of electromagnetic fields.

Test particle technique
In order to investigate in greater detail the resonant-particle behaviour, we can resort to test-particle diagnostics. Following [26], we select the set of test particles in the following way.
(1) We identify the coordinates, ( (2) We observe that the perturbed motion preserves a second invariant (besides the magnetic moment M), provided that the mode is characterised by a single mode number and a constant frequency; namely, the quantity C, defined, at the lowest order, as

Simulation parameters
The simulation parameters are the same as those used in [47]. A tokamak equilibrium is considered, characterised by aspect ratio = R a 10 0 and safety factor = + -( )( ) q q q q r a a 0 0 2 , with q 0 =1.9 and q a =2.3. Thermal ions are characterised by flat density and temperature radial profiles, such that their diamagnetic effect can then be ignored. Kinetic thermal-ion compressibility effects are instead retained, and they give raise to the formation of a kinetic thermal-ion gap in the shear-Alfvén continuum spectrum. As a consequence, a BAE mode exists and localised around » r a 0.5, i.e. the q=2 rational surface, in our current equilibrium. Kinetic thermal ions are characterised by an isotropic Maxwellian initial distribution function. The fast-ion initial distribution function is instead an anisotropic Maxwellian with a single pitch angle a = V E m cos 2 H . In the following, we will report the results of simulations related to modes characterised, respectively, by toroidal numbers n=2, 3 and 4. Different toroidal numbers correspond to different shapes of the shear-Alfvén continuum; then, both mode structure and mode frequency change. In principle, we could find that, for different n, different regions of the phase space yield the most relevant contribution to the destabilisation of the mode. Moreover, as the mode frequency does not change proportionally to the toroidal numbers, the relative weight of angular momentum and kinetic energy in the C expression changes with n. This means that, even the same couple ( ) M C , 0 0 would correspond to a different radial profile of the parallel velocity V and, then, of the resonance frequency. We will show that such differences affect the scaling of the mode saturation amplitude with the linear growth rate.

Different toroidal number n simulations
In this section we will analyze the effect of increasing the toroidal number on the scaling of the mode saturation amplitude. We will show that both the variation of the radial profile of the resonance frequency and the narrowing of the mode radial profile play a role.

Overview of mode evolution
In [47], the results of n=2 simulations were reported. It was shown that both co-passing and counter-passing fast ions drive the mode unstable via transit resonance, with very similar real frequencies. Growth rates increase with increasing fast-ion fraction. We start comparing the results obtained for a given n=2 simulation with those obtained for n=3 and n=4. In particular, we consider the case characterised by fast ion density = n n 0.0014 H i and pitch angle α such that a = cos 0.5 for co-passing fast ions, a =cos 0.5 for counterpassing fast ions. The other parameters are those given in section 3.4 and in [47].
The normalised mode structure of the scalar potential and the mode frequency spectrum are shown in figure 4 for each n. The unstable modes are all located near the shear-Alfvén continuum accumulation point for each n. Mode frequencies are very similar for co-passing fast ion driven modes and counter-passing ones for all n values. They increase as the toroidal mode number increases. The mode structure is narrower for larger toroidal mode number n. In figure 5, the power exchange between mode and particles is shown in the velocity space ( ) M V , . The power exchange is integrated over a toroidal shell centred around the mode (   r 0.2 0.8). Boundaries between the passing-particle region and the trapped-particle one at r=0.2 and r=0.8 are plotted. The power exchange is obviously maximum for passing particles along the pitch angle direction. No significant dependence on the toroidal number is observed for the value of magnetic moment where the power exchange is peaked. Figure 6 shows the nonlinear time evolution of the scalar potential amplitude (normalised to T e H H 0 , where T H0 being the on-axis fast ion temperature and e H the charge of fast ions) at the peak of the radial structure shown in figure 4 for each toroidal number n. Each case shows amplitude oscillations after the mode reaches the first saturation. For the considered case, all simulations show that the co-passing fast ion driven mode reaches higher saturation level than the one driven by counter-passing fast ions. For co-passing fast ion driven modes, the larger toroidal number yields larger growth rate. The largest saturation level occurs for n=2 mode, with the n=3 mode reaching close levels, though n=2 mode is characterised by much smaller growth rate. Counterpassing fast ion driven ones yield analogous results, apart from that related to the saturation level of the n=2 mode: in this case, the lower growth rate is not accompanied by a larger saturation level.

Scaling of mode saturation amplitude with the growth rate
Let us now examine how the mode saturation amplitude varies as the growth rate increases. The results obtained at different values of n n H i or different toroidal mode numbers are reported in figure 7, for modes driven by copassing fast ions (left) and counter-passing ones (right). A clear transition from quadratic to linear scaling is observed in the counter-passing fast ion case, for all the values of the toroidal mode number n. The growth rate value at which the transition occurs slightly increases with increasing n. At the same time, the corresponding mode saturation amplitude decreases. The situation is less defined in the co-passing fast ion case. In particular, plane, for the BAE driven by fast ions, as obtained from XHMGC simulations for different choices of the toroidal mode number n. The normalised radial structure of the dominant harmonic for the scalar potential is also reported for each case (dashed lines). Top frames refer to modes driven by co-passing fast ions; bottom frames, to modes driven by counter-passing fast ions. Here w A0 is the one-axis Alfvén speed.   no transition is observed for n=2, with the mode saturation amplitude scaling linearly over the whole considered range of n n H i . For n=3, although a quadratic scaling can not be revealed, a deviation from the linear scaling is observed when moving from the stronger to the weaker cases. It is still true that the transition growth rate value increases with n (for n=2 and n=3 we can only set upper limits to such value, with the limit relative to the n=2 case being surely lower than that pertaining to the n=3 case). Nothing can be said, instead, about the variation of the transition values of the mode saturation amplitude.
It has been demonstrated [26] that mode saturation occurs as the flattening of the resonant particle distribution function profile (which represent the free-energy source for mode instability) extends over the whole region where the mode-particle interaction can take place. This region is limited both by the finite radial structure of the mode (Dr mode ) and the finite radial extension where the resonance condition, w w g -| ( )|  r res L , is satisfied (Dr res ). Saturation will be reached as In [47,48], it has been argued, on the basis of a nonlinear pendulum model, that the radial width of the flattening region scales as flat L with A(t) being the instantaneous mode amplitude, at time t. Linearising the frequency mismatch as r res res , we get the following approximate expression for the radial width of the resonance: Taking into account that the mode width exhibits a negligible dependence on the growth rate, we see that the quadratic scaling for the mode saturation amplitude is obtained when the most stringent constraint is represented by the resonance width (resonance detuning regime); the linear scaling, when it is given by the finite mode width (radial decoupling regime). Transition from the former to the latter regime is expected for g~D | | ( ) S r , 1 0 In the present case, the mode structure has been inspected in figure 4. In section 3.5.3 we will analyse the radial structure of the resonance frequency, in order to explain the results shown in figure 7.

Radial resonance structure
As sketched in figure 2, we need to compare the resonance width with the radial mode structure given in  where w D is the precession frequency and w b is the transit frequency defined as In the above definition, σ is the sign of parallel velocity V, f D is the change in toroidal angle over the bounce time t b defined as t q q º ∮ḋ b andq is the safety factor integrated along the particle orbit s p q º ∮ ( ) q q 2 d . In the present cases, the dominant bounce harmonic k for both co-passing and counter-passing fast ions is k=1, as discussed in detail in [47]. Here, we show how the resonance radial profile changes for different toroidal numbers in the same equilibrium. By initialising test particles along the curves V(r) shown in figure 8, the corresponding resonance profiles are shown in figure 9. For passing particles, the precession frequency is so small that its contribution to the whole resonance frequency is in fact negligible. The red line represents the radial profile of the transit frequency w b , the green dashed line is the profile of the term proportional to  -(¯) nq m , which is the contribution of k V  . Adding these two values (k=1, as stated above) gives the resonance frequency radial profile w ( ) r res . The fast-ion transit frequency monotonically decreases along radius for both copassing and counter-passing fast ions. On the other hand, for the considered q profile, the term related to k V  monotonically increases from for co-passing fast ions; it decreases for counter-passing fast ions; and changes the sign when crossing = k 0  radial position which is corresponding to the q=2 rational surface in the current equilibrium. The final resonance profiles are given in figure 9 (bottom), showing that, co-passing fast ions have flatter resonance than counter-passing ones. For n=2 case, which is discussed in [47], co-passing fast ions have an extremely flat profile for inner values of the radius, so that the resonance width is much larger than the mode width for all the considered growth rate. For the counter-passing fast ions, the radial resonance profile is steeper, so that the resonance width for low growth rate is smaller than the mode width. As the growth rate increases, a transition from resonance detuning to radial decoupling is expected.
For the n=3 case, transit frequency radial profile is slightly flattened by the fact that the resonant particle velocity profile, V(r), is flatter. However, this is compensated by the steeper profile of the k V  term. The resonance frequency radial profile has then a larger gradient in the region around the mode localisation. Figure 10 shows the dependence of the resonance frequency steepness and the mode width on the toroidal number, for both co-passing and counter-passing ions. It allows us for interpreting the results shown in figure 7 and compare them with the predictions based on the nonlinear pendulum model. Inserting these results in equations (10) and (11), we find, for the transition values for growth rate and mode saturation amplitude the dependence on the toroidal number shown in figure 11 (note that the mode amplitude is determined apart from a proportionality constant). These predictions appear to be in qualitative agreement with most of the simulation results reported in figure 7. In particular, as far as the counter-passing fast ion case is concerned, the decrease of the transition mode saturation amplitude with n is predicted, as well as, the weak n dependence of the transition growth rate. The latter dependence on the positive n dependence of the steepness of the resonance frequency, which compensate the negative dependence of the mode width (see figure 10).
With regard to the co-passing fast ion case, figure 11 (left) predicts a positive n dependence for the transition growth rate. This is in agreement with the fact that our simulations show the transition from the resonance detuning regime to the radial decoupling one only for n=4, while for the lower n cases, only an upper limit to  the transition growth rate can be recovered, as a proper resonance detuning regime is not observed in the growth rate range considered (for n=3, however, a clear deviation from the linear scaling is observed for the weaker cases). This dependence is motivated by the strong positive increase of the resonance frequency steepness with n, which beats the opposite dependence of the mode width (see figure 10 and equation (10)). Figure 11 (right) predicts a positive dependence on the toroidal mode number for the transition saturation amplitude as well. The fact that no transition is observed in our simulations for the lower-n cases does not allow a comparison between model results and simulation ones.
Note that equation (8) is a quite rough approximation for the resonance width if the growth rate is not too small and the dependence of the resonance frequency on the radial distance from the resonant surface is not linear. A slightly better quantitative agreement can be obtained, for the transition values of the growth rate by computing Dr res directly from the full resonance frequency profile, as res res res L  Figure 12 compares the values of the resonance width obtained from equation (14) with the measured values of the mode width at different toroidal numbers. Figure 13 reports the transition growth rate values recovered from figure 12 and the corresponding saturation mode amplitudes, computed, consistently with equation (9), multiplying the transition growth rate by the mode width. We observe that, while no relevant difference is obtained for the transition saturation amplitude, the transition growth rate for the counter-passing fast ion case exhibits a weak positive n dependence. Moreover, for the co-passing fast ions, the transition growth rate values result smaller than those reported in figure 11. Both these results appear to be in better agreement with the results reported in figure 7.  (14), compared with the measured mode radial width, at different toroidal numbers, for co-passing (left) and counter-passing (right) fast ions. Figure 13. Values of the transition growth rate measured from figure 12 and the corresponding value of the saturation mode amplitude, computed multiplying the transition growth rate by the mode radial width.

Summay and conclusions
In this paper, the saturation mechanism is investigated for BAEs driven unstable by anisotropic fast-ion populations (co-passing or counter-passing ions) for different single toroidal numbers. Numerical simulations are performed by the hybrid MHD-particle code XHMGC. The fact that a single toroidal mode number is considered and that the mode frequency does not appreciably vary in time allows for the existence of a second constant of (the perturbed) motion, apart from the magnetic moment M: namely, the quantity w = f C P nE. Cutting the phase space into slices orthogonal to the axes M and C, the dynamics of each slice is independent on that of other slices. We can then look at the evolution of that slice (or those slices) where the mode-particle power transfer is maximum to explore the nonlinear evolution of the system. Once determined the relevant (M, C) slice, the resonance frequency depends only on the radial coordinate. In our case, the shape of the resonance frequency is mainly determined by the bounce (transit) frequency for passing particles. The resonance frequency profile is steeper for the counter-passing fast ions than for the co-passing ones due to the opposite sign of their k V  terms in the resonance condition. In both cases, the steepness increases with the toroidal mode number. Correspondingly, the width of the radial region where the resonance condition is satisfied decreases with n. The same happens to the radial mode width.
Taking into account that the mode-particle power exchange can take place only in the radial region where the mode amplitude is not negligible and the resonance condition is satisfied, two different saturation regimes can be distinguished, depending on the relative importance of these two characteristic widths. In the first regime, the most stringent constraint is given by the resonance condition (D < D r r res mode ). Saturation occurs when the fast-ion density flattening induced by the radial flux associated to resonant particles captured in the wave potential well extends over the whole region where the frequency matching is satisfied. This regime is called resonance detuning regime and is characterised by a quadratic scaling of the mode saturation amplitude with the linear growth rate. In the second regime, called radial decoupling regime, the relevant constraint comes from the finite mode width (D < D r r mode res ). Saturation is characterised by a density flattening extending over the whole region where the mode amplitude is not vanishing, and a linear scaling of such amplitude with the growth rate.
In the case considered in this paper, the radial profile of the resonance frequency for the co-passing fast ions is flatter than that for the counter-passing ones. This results in a much greater resonance width in the former case, such that, for low toroidal mode numbers, only the radial decoupling regime is observed in the explored growth-rate range. The transition from the resonance detuning regime to the radial decoupling one is appreciated, for co-passing ions, only for the larger toroidal number considered, n=4. Modes driven by counter-passing fast ions, instead, exhibit such transition even for the lower n values.
These results and the n-dependence of the growth-rate value at which the resonance-detuning/radialdecoupling transition occurs and the corresponding saturation mode amplitude, for both co-passing and counter-passing fast ions, compare fairly well with the prediction based on a simplified nonlinear pendulum model.
Note that the analysis presented in this paper does not encompass several aspects of Alfvén mode nonlinear dynamics, able to play an important role for different equilibria, like frequency chirping, large nonlinear modification of mode structure and synergic interaction of different toroidal harmonics, which will be discussed in later works.