# On the relationship between residual zonal flows and bump-on tail saturated instabilities

Xavier Garbet, P Donnel, C Ehrlacher, E Caschera, R Dumont, M Faganello, Virginie Grandgirard, Philippe Ghendrih, M Idouakass, Maxime Lesur, et al.

► To cite this version:

Xavier Garbet, P Donnel, C Ehrlacher, E Caschera, R Dumont, et al. On the relationship between residual zonal flows and bump-on tail saturated instabilities. 2016. hal-01356331

HAL Id: hal-01356331

https://hal.science/hal-01356331

Preprint submitted on 26 Aug 2016

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.
# On the relationship between residual zonal flows and bump-on tail saturated instabilities.

X. Garbet¹, P. Donnel¹, C. Ehrlacher¹, E. Caschera¹, R. Dumont¹, M. Faganello², V. Grandgirard¹, P. Ghendrih ¹, M. Idouakass², M. Lesur³, Y. Sarazin¹, D. Zarzoso²

¹CEA, IRFM, F-13108 St.Paul-lez-Durance cedex, France

² Aix-Marseille Université, CNRS, PIIM UMR 7345, Marseille, France

³ IJL, UMR 7198 CNRS, Université de Lorraine, 54500 Vandoeuvre-lès-Nancy, France

E-mail: xavier.garbet@cea.fr

**Abstract.** A connection is established between two classical problems: the non linear saturation of a bump-on tail instability in collisionless regime, and the decay of a zonal flow towards a finite amplitude residual. Reasons for this connection are given and commented.

## 1. Introduction

This paper establishes a connection between two well-known analytical results that are currently used for code verification. The first one is the Rosenbluth-Hinton theory of zonal flow decay [1], which predicts an evolution towards a finite asymptotic amplitude called residual flow. This prediction is routinely used to test and benchmark gyrokinetic codes. The second result addresses the saturated value of the bump-on tail instability in collisionless plasmas. It predicts a bounce frequency of deeply trapped particles that is proportional to the linear growth rate, or equivalently a saturated mode amplitude that goes like the square of the linear drive (see [2] for an overview). This result is an extension of a seminal paper by O'Neil that addresses Landau damping in a nonlinear regime [3] (see also the work by Mazitov [4]). It is used to test codes that solve 2D (x,v) Vlasov/Poisson problems. In fact MHD modes driven by energetic particles often exhibit a quadratic scaling of the mode amplitude with the linear growth rate when the drive is weak, in line with this prediction (cf. [5, 6] and references therein). It turns out that these two results, zonal flow residual and saturated level of bump on tail instability, are connected. The reason for this relationship, which is not entirely fortuitous, is detailed and commented in this paper.

The paper is organized as follows. Section 2 summarizes the O'Neil calculation of Landau damping in non linear regime and its extension to the bump on tail instability. The Rosenbluth-Hinton calculation of the time evolution of zonal flow is reminded in section 3. Finally the relationship between residual zonal flows and bump-on tail saturated instabilities is discussed in section 4. A conclusion follows.

## 2. Landau damping in non linear regime - extension to the bump on tail instability

The seminal O'Neil paper [3] aimed at studying Landau damping in collisionless regime in the non linear regime where particle trapping matters. The time evolution of the distribution
function and electric field were calculated. O'Neil realized that the asymptotic state can be captured by calculating the energy transfer from the wave to particles. This calculation was a source of inspiration for many subsequent works (e.g. [7, 8]). The O'Neil work was soon extended to the bump-on tail instability. The saturated state can be calculated by stating that the kinetic energy of resonant particles is fully converted to wave energy. The key point is then to calculate the distribution function. O'Neil prescription consists in using a phase mixing ansatz. The later states that the relaxed distribution function is the average of the unperturbed one along trajectories, due to phase mixing. Hence it is flat within the island, and joins smoothly the unperturbed distribution function outside the separatrix. This prescription is consistent with the detailed calculation of the mode evolution done by O'Neil in the same work. Berk's prescription rather corresponds to a discontinuous distribution function that is flat in the island, and joins the unperturbed distribution right at the separatrix [5]. A third option is to choose a flat distribution function within an annulus, whose width is the island width [2]. Of course other choices are possible as to the contribution of particles near the island separatrix.

## 2.1. Summary of the O'Neil calculation

Let us start from an arbitrary integrable Hamiltonian system that can be described by a set of angle/action variables $(\Theta, \mathbf{J})$ such that

$$ \begin{aligned} \frac{d\Theta}{dt} &= \frac{\partial H_{eq}}{\partial \mathbf{J}} = \mathbf{\Omega}_{eq}(\mathbf{J}) \\ \frac{d\mathbf{J}}{dt} &= -\frac{\partial H_{eq}}{\partial \Theta} = 0 \end{aligned} $$

where $H_{eq}(\mathbf{J})$ is the unperturbed Hamiltonian. This system is perturbed by an oscillation, i.e.

$$ H(\Theta, \mathbf{J}, t) = H_{eq}(\mathbf{J}) - h \cos (\mathbf{n} \cdot \Theta - \omega t) \quad (1) $$

We look for an asymptotic steady state, i.e. the amplitude $h$ of the perturbed Hamiltonian is a constant (in the original O'Neil paper, the evolution of $h(t)$ was computed). We also introduce the resonance frequency $\Omega(\mathbf{J}) = \mathbf{n} \cdot \mathbf{\Omega}_{eq}(\mathbf{J}) - \omega$. The resonant surface is defined by the condition $\Omega(\mathbf{J}) = 0$ in the action space. Near the resonant surface, the actions are of the form $\mathbf{J} = \mathbf{J}_R + \mathbf{n}I$, where $\mathbf{J}_R$ spans the resonant surface, while the action $I$ measures the distance to the resonant surface. The resonant frequency reads $\Omega = CI$, where $C$ is related to the Hessian matrix of the Hamiltonian calculated at the resonant surface

$$ C = n_i n_k \left. \frac{\partial^2 H}{\partial J_i \partial J_k} \right|_{\Omega=0} \quad (2) $$

We assume that $C$ is positive without any loss of generality. The trajectory equations read

$$ \begin{aligned} \frac{d\xi}{dt} &= \Omega \\ \frac{d\Omega}{dt} &= -\omega_b^2 \sin \xi \end{aligned} $$

where $\xi = \mathbf{n} \cdot \Theta - \omega t$ is an angle, and $\omega_b = \sqrt{Ch}$ is the bounce frequency of deeply trapped particles. The variable $\Omega$ is quite convenient since it scales as a pulsation. Nevertheless the variable $I$ presents also some advantages as it is an action conjugate to $\xi$. It is reminded here that the Jacobian is equal to 1 when moving from one to another set of conjugate variables. The new Hamiltonian $\mathcal{H}$ is defined as

$$ \mathcal{H} = \frac{1}{2} \Omega^2 - \omega_b^2 \cos \xi \quad (3) $$
$\mathcal{H}$ is a motion invariant, so that the system is still integrable. Its contour lines exhibit the shape of an island (Fig.1). Eqs.(3) show that $(\xi, \Omega)$ are conjugate variables with respect to the Hamiltonian $\mathcal{H}$. The set $(\xi, \Omega)$ can be completed by a set of conjugate variables $(\mathbf{\Theta}_R, \mathbf{J}_R)$, which span the resonant surface. The phase space volume element reads

$$d^3\mathbf{\Theta}d^3\mathbf{J} = d^2\mathbf{J}_R d^2\mathbf{\Theta}_R dId\xi \quad (4)$$

To simplify the notations, the explicit dependence on $(\mathbf{\Theta}_R, \mathbf{J}_R)$ will be omitted in the following. In other words, the problem becomes 2D, and decoupled from the dynamics along the resonant surface. In the original O'Neil calculation, which addresses a 1D electrostatic wave, the following correspondences hold

$$h = e\phi \quad ; \quad \Omega = kv - \omega \quad ; \quad \xi = kx - \omega t \quad ; \quad C = \frac{k^2}{m} \quad ; \quad \omega_b = k \left(\frac{e\phi}{m}\right)^{1/2} \quad (5)$$

where $k$ is the wavenumber, $x$ the position, $v$ the velocity, and $\phi$ the electric potential. We introduce an average ("bounce" average) over time for any function $f(\mathcal{H}, \xi, \epsilon_\Omega)$ where $\epsilon_\Omega$ is the sign of $\Omega$, which can be written as

$$\langle f \rangle = \frac{1}{\tau} \int_{-\pi}^{\pi} \frac{d\xi}{2\pi} \frac{1}{|\Omega|} f(\mathcal{H}, \xi, \epsilon_{\Omega}) \quad (6)$$

for passing particles $(\omega_b^2 \le \mathcal{H} \le \infty)$, while

$$\langle f \rangle = \frac{1}{2\tau} \int_{-\xi_0}^{\xi_0} \frac{d\xi}{2\pi} \frac{1}{|\Omega|} (f(\mathcal{H}, \xi, 1) + f(\mathcal{H}, \xi, -1)) \quad (7)$$

for trapped particles $(-\omega_b^2 \le \mathcal{H} \le \omega_b^2)$. Here $\Omega$ is to be understood as a function of $(\mathcal{H}, \xi, \epsilon_\Omega)$

$$\Omega = \epsilon_{\Omega} \sqrt{2} (\mathcal{H} + \omega_b^2 \cos \xi)^{1/2} \quad (8)$$

The bounce angle $\xi_0$ is the positive root of $\omega_b^2 \cos \xi_0 = \mathcal{H}$. The normalizing time $\tau$ is defined as

$$\tau = \int_{-\xi_0}^{\xi_0} \frac{d\xi}{2\pi} \frac{1}{|\Omega|} \quad (9)$$

with $\xi_0 = \pi$ for passing particles. Note that $\tau$ is a continuous function of $\kappa$, and differs from the bounce period $T_b$. More precisely $\tau = \frac{T_b}{2\pi}$ for passing particles and $\tau = \frac{T_b}{4\pi}$ for trapped particles.

## 2.2. Saturation level

A subtle point is that the kinetic energy should be calculated in the laboratory frame, while the calculation of the distribution function is done in the wave frame. The Hamiltonian in the laboratory frame $\mathcal{H}_{lab}$ is related to the one in the wave frame $\mathcal{H}$ via the relation $\mathcal{H} = \mathcal{H}_{lab} - \omega I$, or equivalently $\mathcal{H}_{lab} = \mathcal{H} + \Omega \frac{\omega}{C}$. The Hamiltonian $\mathcal{H}$ does not participate in the resonant transfer of energy due to parity reason. Indeed $f = F - F_{eq}$, where $F_{eq}$ is the unperturbed distribution function, is an odd function of $\Omega$. Therefore the resonant transfer of kinetic energy is

$$\Delta E = \frac{\omega}{C} \int_{-\infty}^{+\infty} dI \int_{-\pi}^{+\pi} \frac{d\xi}{2\pi} f\Omega \quad (10)$$
<figure><img src="image_1.png" /><figcaption>Figure 1. Island shape ($\omega_b = 1$).</figcaption></figure>

<figure><img src="image_2.png" /><figcaption>Figure 2. Schematic view of the distribution function at the island O point</figcaption></figure>

Moreover $F_{eq}$ which can be expanded near the resonant surface $\Omega = 0$, i.e. $F_{eq} = F_0 + F'_0\Omega$. The phase mixing prescription states that the relaxed distribution function is the average of the initial one over trajectories, i.e.

$$f = F - F_{eq} = F'_0 [\langle\Omega\rangle - \Omega] \quad (11)$$

Note that $\langle\Omega\rangle$ is zero for trapped particles, since $\Omega$ changes sign when the particle moves back and forth due to trapping. A schematic view of the distribution function at the O point of the island is shown in Fig.(2). A change of variables yields the following result

$$\Delta E = 2\frac{\omega}{C^2}F'_0\int_{-\omega_b^2}^{+\infty} d\mathcal{H}\tau \{\langle\Omega\rangle^2 - \langle\Omega^2\rangle\} \quad (12)$$

It appears that most of the energy given by particles to the mode comes from trapped particles, and to a lesser extent from barely passing particles. When balancing the energy transfer with the mode energy density $\Lambda_\omega h^2$, one gets a prediction for the field amplitude $\Delta E + \Lambda_\omega h^2 = 0$. The final result can be formulated as

$$\omega_b = \mathcal{K}_{BT}\gamma_L \quad (13)$$

where $\gamma_L = \frac{\pi}{4}\omega\frac{F'_0}{\Lambda_\omega}$ is the linear growth rate. The constant $\mathcal{K}_{BT}$ is given by the following relation

$$\mathcal{K}_{BT} = \frac{32}{\pi} \int_{0}^{\infty} \frac{d\kappa}{\kappa^3} \frac{\tau}{\omega_b} \{\langle\Omega^2\rangle - \langle\Omega\rangle^2\} \quad (14)$$

where

$$\Omega = \frac{2\omega_b}{\kappa}\epsilon\Omega \left[1 - \kappa^2 \sin^2\left(\frac{\xi}{2}\right)\right]^{1/2} \quad (15)$$

and $\kappa$ is the trapping parameter defined as $\kappa^2 = \frac{2\omega_b^2}{\mathcal{H}+\omega_b^2}$. Trapping corresponds to $1 \le \kappa^2 \le \infty$ while passing particles are characterized by $0 \le \kappa^2 \le 1$.
### 3. Zonal flow residual

#### 3.1. Formal expression

In 1998 Rosenbluth and Hinton calculated the time evolution of a zonal flow and predicted the asymptotic value, called zonal flow residual [1]. The calculation consists in initializing the electric potential with a perturbation of the form $\phi(\mathbf{x}, t = 0) = \phi_0 e^{ik\psi}$, where $k$ is the radial wavenumber, and $\psi$ is minus the poloidal magnetic flux normalized to $2\pi$, also equal to $-A_\varphi R$, where $A_\varphi$ is the toroidal component of the vector potential, and $R$ the major radius. It is a function of the radial coordinate, labeled $r$. The evolution of $\phi(\mathbf{x}, t)$ is calculated by solving the gyrokinetic Vlasov equation coupled to the electroneutrality equation. The ratio of the potential residual to its initial value is

$$ \frac{\phi_{\infty}}{\phi_0} = \frac{1 - \int d^3\mathbf{p} F_{eq} J_0^2}{1 - \int d^3\mathbf{p} F_{eq} J_0^2 J_b^2} \qquad (16) $$

where $J_0$ and $J_b$ are related to averages over the cyclotron and drift motions, $F_{eq}$ is the unperturbed distribution function. These are scalars for Fourier modes.

#### 3.2. Field shielding

The gyroaverage function is the usual one, and is approximated by the following form in the large scale limit $J_0^2 \simeq 1 - \frac{1}{2}k_\perp^2\rho_c^2$, where $\rho_c$ is the kinetic ion Larmor radius. The guiding center radial coordinate reads $\psi = \bar{\psi} + \hat{\psi}$, where $\bar{\psi}$ is the time average position and $\hat{\psi}$ describes the shift with respect to the reference magnetic surface $\psi = \bar{\psi}$. For a periodic potential, the transit average operator is defined as

$$ J_b = \langle e^{ik\hat{\psi}} \rangle \simeq 1 - k^2 \langle \hat{\psi}^2 \rangle \qquad (17) $$

where the bracket is an average over the drift motion. At first order in $k^2\rho_c^2$, the shielding term reads

$$ 1 - \int d^3\mathbf{p} F_{eq} J_0^2 J_b^2 \simeq k^2 \left\langle \frac{1}{2}\rho_c^2 + \int d^3\mathbf{p} F_{eq} \hat{\psi}^2 \right\rangle = k^2 (\rho_i^2 + \delta^2) \qquad (18) $$

where $\rho_i = \frac{mT}{eB}$ is the thermal ion Larmor radius, and

$$ \left(\frac{d\psi}{dr}\right)^2 \delta^2 = \int d^3\mathbf{p} F_{eq} \langle \hat{\psi}^2 \rangle \qquad (19) $$

is the radial size of the drift orbit squared. The expression of the radial displacement $\hat{\psi}$ is obtained by using the conservation of the canonical toroidal kinetic momentum $P_\varphi = m\frac{I}{B}v_\parallel - e\psi$, namely

$$ \hat{\psi} = \frac{m}{e} I \left( \frac{v_\parallel}{B} - \left\langle \frac{v_\parallel}{B} \right\rangle_b \right) \qquad (20) $$

where the index 'b' indicates a bounce average for particles trapped in the magnetic field due to the mirror force, while it labels a transit time average for passing particles. Here $I = B_T R$ is the product of the toroidal field by the major radius, $B$ is the modulus of the magnetic field, $v_\parallel$ is the particle velocity along the magnetic field, $m$ and $e$ are the ion mass and charge. This gives an expression of the width $\delta$

$$ \left(\frac{d\psi}{dr}\right)^2 \delta^2 = \frac{1}{2} \left(\frac{mI}{e}\right)^2 \int d^3\mathbf{p} F_{eq} \left[ \left\langle \left(\frac{v_\parallel}{B}\right)^2 \right\rangle_b - \left\langle \frac{v_\parallel}{B} \right\rangle_b^2 \right] \qquad (21) $$

This is the main result derived by Rosenbluth and Hinton in arbitrary geometry. We now restrict the analysis to circular concentric magnetic surfaces with large aspect ratio, so that $r$ is now the minor radius of a magnetic surface, and $\theta$ is the poloidal angle. The hamiltonian reads

$$ H = \frac{1}{2}mv_\parallel^2 + \mu B_0 (1 - \epsilon \cos \theta) \qquad (22) $$
where $\epsilon = r/R_0$ is the inverse aspect ratio. This Hamiltonian is obviously identical to the one addressed in the previous section, with the following correspondences

$$h = \mu B_0 \epsilon ; \quad \Omega = \frac{v_{\parallel}}{q R_0} ; \quad \xi = \theta ; \quad C = \frac{1}{mq^2 R_0^2} ; \quad \omega_b = \frac{1}{q R_0} \left( \frac{\mu B_0 \epsilon}{m} \right)^{1/2} \quad (23)$$

and $\mathcal{H} = C(H - \mu B_0)$. Hence the problem is related again to particle trapping. However the reader should be warned that trapping here is due to a magnetic mirror force, not to the mode itself. The amplitude of the latter is very small, so that the problem can be treated linearly with respect to $\phi$. Moreover, it can be easily verified that for a Maxwellian distribution function

$$F_{eq} d^3 \mathbf{p} = (2\epsilon)^{1/2} \frac{2}{\sqrt{\pi}} dE E^{1/2} e^{-E} \frac{2d\kappa}{[(1-\epsilon)\kappa^2 + 2\epsilon]^{3/2}} \omega_b \tau \quad (24)$$

where $E$ is the energy normalized to the temperature. In the limit of large aspect ratio $\epsilon \to 0$ and precluding that the largest contribution comes from the separatrix neighborhood $\kappa \simeq 1$, a compact expression of the zonal flow residual Eq. (16) is found

$$\frac{\phi_{\infty}}{\phi_0} = \frac{1}{1 + C_{RH} \frac{q^2}{\sqrt{\epsilon}}} \quad (25)$$

where

$$C_{RH} = 3\sqrt{2} \int_{0}^{\infty} \frac{d\kappa}{\kappa^3} \frac{\tau}{\omega_b} \left\{ \langle \Omega^2 \rangle - \langle \Omega \rangle^2 \right\} \quad (26)$$

## 4. Relationship between residual zonal flows and bump-on tail saturated instabilities

### 4.1. Numerical values and metrology

The two constants $C_{RH}$ and $K_{BT}$ are related to the quantity $\mathcal{I}$

$$\mathcal{I} = \int_{0}^{\infty} \frac{d\kappa}{\kappa^3} \frac{\tau}{\omega_b} \left\{ \langle \Omega^2 \rangle - \langle \Omega \rangle^2 \right\} \quad (27)$$

through the relations $K_{BT} = \frac{32}{\pi}\mathcal{I}$ and $C_{RH} = 3\sqrt{2}\mathcal{I}$. A straightforward calculation shows that

$$\begin{aligned} \mathcal{I} = & 2 \int_{0}^{1} \frac{d\kappa}{\kappa^4} \left[ \frac{2}{\pi} E(\kappa^2) - \frac{\pi}{2} K(\kappa^2) \right] \\ & + 2 \int_{1}^{+\infty} \frac{d\kappa}{\kappa^3} \left[ \left( \frac{1}{\kappa^2} - 1 \right) \frac{2}{\pi} K\left( \frac{1}{\kappa^2} \right) + \frac{2}{\pi} E\left( \frac{1}{\kappa^2} \right) \right] \end{aligned} \quad (28)$$

where $K$ and $E$ are complete elliptical functions of first and second kinds. A numerical calculation yields $\mathcal{I} \simeq 0.38$ [12], which implies $C_{RH} \simeq 1.61$ and $K_{BT} \simeq 3.84$. The value $C_{RH} \simeq 1.6$ has been recovered in a number of gyrokinetic codes, and is now considered as a "must do" verification test in gyrokinetic computation. However less success was met for the bump-on tail instability since numerical simulations rather find $K_{BT} \simeq 3.2$ [9, 11]. An alternative recipe consists in postulating that trapped particles only contribute to the mode saturation. It turns out that $\mathcal{I}$ is analytic in that specific case, namely $\mathcal{I} = \frac{8}{9\pi}$, and $K_{BT} = \frac{256}{9\pi^2} \simeq 2.9$ [5]. Another calculation based on O'Neil prescription yields $K_{BT} \simeq 3.5$ [10]. Finally a radical option is to choose a flat distribution function within an annulus, whose width is the island width. This procedure is fully analytic and leads to $K_{BT} = \frac{32\sqrt{2}}{3\pi} \simeq 4.8$ [2]. All these results differ from the numerical findings. A value that is too low is found when minimizing the role of passing particles, while an overweight of passing particles leads to an estimate that is too large. So it appears that
only a fraction of particles near the separatrix participate to the energy exchange. In particular it looks like barely passing particles do not contribute as expected. This is not surprising since the bounce time is very large near the separatrix (and infinite right on it), so that a mixing argument is questionable for particles in this region of the phase space. This is consistent with the analysis of the contributions of particles to the Rosenbluth-Hinton constant $C_{RH}$ (and therefore of $K_{BT}$). It appears that barely trapped particles provide the largest contribution (see Fig. 3). This structure in the velocity space was verified numerically by gyrokinetic codes [13]. Another important test is the dynamics of the mode evolution, which is also analytical for both problems, and quite different. Time evolution has been verified numerically with success in both cases. As a final note, we stress that a steady solution is not unique by far for bump-on tail instabilities, in particular when dissipation is accounted for. Quite often the system evolves to unsteady solutions, which can be periodic, or chaotic, thus leading to frequency chirping [14, 15, 11, 6]. These important developments are beyond the scope of the present paper.

<figure><img src="image_3.png" /><figcaption>Figure 3. Contour lines of the integrand of $\mathcal{I}$ in the $v_\parallel, v_\perp$ space at $\theta = 0$.</figcaption></figure>

## 4.2. Pure coincidence or some common physics?

The previous sections show that a relationship connects the bump-on tail constant with the Rosenbluth-Hinton constant, i.e. $K_{BT} = \frac{16\sqrt{2}}{3\pi}C_{RH}$. This is rather unexpected as the two results correspond to different physics. On the one hand, the bump-on tail saturation comes from a flattening of the distribution function near the resonance surface, leading to a weakening of the drive until saturation. On the other hand, the zonal flow residual comes from a shielding effect of the initial perturbation due to the particle response to the field. So the fast conclusion is that this relationship is fortuitous. There are nevertheless common features : particle trapping, absence of collisions, and more importantly momentum conservation. The relation Eq.(10) shows that the energy transfer in the bump-on-tail instability per unit time can be related to an exchange of momentum between the wave and particles. The power is the drag force $F = \frac{dP}{dt}$ times the wave phase velocity $\frac{\omega}{k}$. Here $P$ is the momentum exchange between the wave and resonant particles

$$P = k \int_{-\infty}^{+\infty} dI \int_{-\pi}^{+\pi} \frac{d\xi}{2\pi} If = m^2 \int_{0}^{\lambda} \frac{dx}{\lambda} \int_{-\infty}^{+\infty} dv \left(v - \frac{\omega}{k}\right) f \quad (29)$$

where $\lambda = \frac{2\pi}{k}$ is the wavelength. The zonal flow residual is also related to momentum conservation, though in a different way. As mentioned above, the shielding effect comes from
the particle displacement due to drift motion. The poloidal variation of the magnetic field in a tokamak can be seen as a static wave, that exchanges momentum with particles. The change of parallel velocity induces a radial displacement, which results from the conservation of toroidal kinetic momentum. This establishes a link between the particle radial shift and the parallel velocity Eq.(20). The displacement is directly proportional to the variance of velocity, hence the *I* integral. In the case of the bump-on tail instability, it comes from the structure of the distribution function. This is probably the reason why the analytic calculation works better for the Rosenbluth-Hinton shielding effect than for the bump-on tail mode saturation: the distribution function results from a complex dynamics, and the assumption of phase mixing is easily broken, in particular near the separatrix.

**5. Conclusion**

In summary it is found that the saturation of a bump-on tail instability in collisionless regime, and the zonal-flow residual are connected. More precisely, the two constants that appear in these classical results, namely the ratio of the trapping frequency to the linear growth rate, and the zonal flow residual normalized to its initial value, are connected via a simple formula. From the numerical point of view, the prediction for the zonal flow residual works better than the one for the bump-on tail instability. This indicates that for zonal flow decay, all particles play their expected role, while for the bump-on tail problem, particles near the separatrix do not behave as foreseen. This means that the distribution function is flat within the island, as predicted, but exhibits a transition near the separatrix sharper than expected on the basis of phase mixing argument. The relationship between the two problems seems at first sight fortuitous, since the physics that is involved is different. However particle trapping and momentum conservation appear as common denominators, thus explaining partially the link between these expressions.

**Acknowledgements**

Drs H. Berk and A. Smolyakov are warmly acknowledged for fruitful discussions. This work has also taken benefit from discussions within the NLED Eurofusion Enabling Research project.

**References**

[1] Rosenbluth M N and Hinton F L 1998 *Phys. Rev. Lett.* **80** 724

[2] Galeev A A and Sagdeev R S 1973 in Review of Plasma Physics vol.7, p.1 "Nonlinear Plasma Theory" [english version: 1979 Consultant Bureau, New York]

[3] O'Neil T 1965 *Phys Fluids* **8** 2255

[4] Mazitov R K 1965 *Zhurnal prikladnoi mekhlaniki i tekhnicheskoi fiziki* **1** 27

[5] Berk H L, Breizman B N, Fitzpatrick J, Pekker M S, Wong H V, and Wong K L, Physics of Plasmas 1996 **3** 1827

[6] Liu Chen and Zonca F 2016 Rev. Mod. Phys. **88** 015008

[7] Isichenko M B *Phys. Rev. Lett.* **78** 2369 (1997)

[8] Manfredi G *Phys. Rev. Lett.* **79** 2815 (1997)

[9] Fried B D, Liu C S, Means R W, and Sagdeev R Z, National Technical Information Service Document No. AD730123 PGG-93, University of California, Los Angeles 1971. Copies may be ordered from the National Technical Information Service, Springfield, Virginia 22161.

[10] Levin M B, Lyubarskii M G, Onishchenko I N, Shapiro V D, and Shevchenko V I 1972 *Zh. Eksp. Teor. Fiz.* **62** 1725 [1972 *Sov. Phys. JETP* **35** 898]

[11] Lesur M, Idomura Y, and Garbet X 2009 *Phys. Plasmas* **16** 092305

[12] Rosenbluth M N et al 1996 *Nucl. Fusion* **36** 55. See Eq.(50) and paragraph thereafter.

[13] Watanabe T-H and Sugama H Nucl. *Fusion* 2006 **46** 24

[14] Berk H L, Breizman B N and Petviashvili N V 1997 *Phys. Lett. A* **234** 213

[15] Vann R G L et al. 2003 *Phys. Plasmas* **10** 623