Long-term dynamics driven by resonant wave-particle interactions: from Hamiltonian resonance theory to phase space mapping

The resonant wave-particle interaction is known to be one of the main drivers of dynamics of such space plasma systems as planetary radiation belts (e.g., Thorne, 2010; Menietti et al., 2012), collisionless shock waves (e.g., Balikhin et al., 1997; Wilson et al., 2007, 2012; Wang et al., 2020), auroral acceleration region (e.g., Chaston et al., 2008; Watt & Rankin, 2009; Mauk et al., 2017), and solar wind (e.g., Krafft & Volokitin, 2016; Kuzichev et al., 2019; Tong et al., 2019; Yoon et al., 2019; Roberg-Clark et al., 2019). The classical quasi-linear theory (Vedenov et al., 1962; Drummond & Pines, 1962) and its generalizations for inhomogeneous plasma systems (Ryutov, 1969; Lyons & Williams, 1984) describe well charged particle resonant interaction with low-amplitude broadband waves (Karpman, 1974; Shapiro & Sagdeev, 1997; Tao et al., 2012a; Camporeale & Zimbardo, 2015; Allanson et al., 2020).

One of the important examples of application of the quasi-linear theory is the Earth radiation belt models that describe energetic electron acceleration and losses due to resonances with electromagnetic whistler-mode waves and electromagnetic ion cyclotron (EMIC) waves (see reviews Thorne et al., 2010; Shprits et al., 2008; Ni et al., 2016; Nishimura et al., 2010; Millan & Thorne, 2007, and references therein). Moreover, the natural inhomogeneity of the background magnetic field and plasma density in the radiation belts can significantly weaken the conditions of applicability of the quasi-linear theory (Solovev & Shkliar, 1986; Albert, 2001, 2010). However, this theory meets difficulties in describing resonances with sufficiently intense waves (Shapiro & Sagdeev, 1997), when the nonlinear effects of phase trapping and phase bunching become important (Omura et al., 1991; Shklyar & Matsumoto, 2009; Albert et al., 2013; Artemyev et al., 2018a). Indeed, sufficiently intense whistler-mode waves are frequently observed in the radiation belts (Cattell et al., 2008; Wilson et al., 2011; Agapitov et al., 2014) and contribute significantly to wave statistics (Zhang et al., 2018, 2019; Tyler et al., 2019). Theoretically, phase trapping and bunching (also called nonlinear scattering) effects are responsible for fast acceleration (e.g., Demekhov et al., 2006, 2009; Omura et al., 2007; Hsieh & Omura, 2017; Hsieh et al., 2020) and losses (e.g., Kubota et al., 2015; Kubota & Omura, 2017; Grach & Demekhov, 2020) of energetic electrons and for the generation of coherent whistler-mode waves (Demekhov, 2011; Katoh, 2014; Katoh & Omura, 2016; Tao, 2014; Omura et al., 2008, 2013; Nunn & Omura, 2012). There are many observational evidences of such nonlinear resonant wave generation (Titova et al., 2003; Cully et al., 2011; Tao et al., 2012b; Mourenas et al., 2015) and of the related electron acceleration/losses (e.g., Foster et al., 2014; Agapitov et al., 2015b; Mourenas et al., 2016a; Chen et al., 2020).

The quasi-linear diffusion theory describes a sufficiently weak scattering in energy/pitch-angle space and operates with a Fokker-Planck diffusion equation for the charged particle distribution function (Andronov & Trakhtengerts, 1964; Kennel & Engelmann, 1966; Lerche, 1968). In contrast to this description, the nonlinear phase trapping assumes a fast transport in energy/pitch-angle space (e.g., Artemyev et al., 2014a; Furuya et al., 2008), when even a single resonant interaction changes significantly the electron’s energy/pitch-angle (e.g., Albert et al., 2013; Artemyev et al., 2018a). This essentially non-diffusive process cannot be directly included into the Fokker-Planck equation. One possible approach is the construction of an operator that would describe fast charged particle jumps in the energy/pitch-angle; this operator can be constructed with the numerical (test-particle) approach (e.g., Hsieh & Omura, 2017; Zheng et al., 2019) or with the analytical calculation of jumps’ probabilities (e.g., Vainchtein et al., 2018). The main advantage of this approach is the inclusion of almost arbitrary (as realistic as needed) wave spectrum and characteristics (e.g., wave modulation and frequency drifts, see Kubota & Omura, 2018; Artemyev et al., 2019b; Hiraga & Omura, 2020). The main disadvantages are an accumulation of numerical errors with running time, and the almost intractable fine details of the energy/pitch-angle space binning needed to simultaneously resolve large jumps due to trapping and small changes due to drift/diffusion.

An alternative approach to the construction of such an operator is a generalization of the Fokker-Planck equation to include effects of phase trapping and phase bunching (Solovev & Shkliar, 1986; Artemyev et al., 2016b, 2017). This approach is based on a fine balance of trappings and bunchings for a single wave system (e.g., Shklyar, 2011; Artemyev et al., 2019a). The main advantage of this approach is that the evolution of charged particle distribution function can be investigated in arbitrary details in presence of phase trapping, phase bunching, and diffusion (Artemyev et al., 2018b, 2019a, e.g.,). The main disadvantage is that there is no straightforward generalization of this approach for multi-wave (multi-resonance) systems. A single-wave resonance results in charged particle transport in the energy/pitch-angle space along 1D curves, so-called resonance surface curves (e.g., Lyons & Williams, 1984; Summers et al., 1998), and the Fokker-Planck equation with trapping was derived for such a quasi-1D system (Artemyev et al., 2016b).

Another alternative for the description of charged particle distribution evolution driven by nonlinear wave-particle interaction (phase trapping and bunching) is the mapping technique that describes the characteristics of the Fokker-Planck equation (Van Kampen, 2003). The classical example of this approach is the Chirikov map (Chirikov, 1979), which describes particle diffusion and is widely used for systems with wave-particle resonances (e.g., Vasilev et al., 1988; Zaslavskii et al., 1989; Benkadda et al., 1996; Khazanov et al., 2013, 2014). Such a map has been constructed for a single-wave system with phase-trapping and phase bunching effects (Artemyev et al., 2020b). In this study we show the generalization of this map for a multi-resonance system.

We consider a strong magnetic field system, where charged particle motion is well gyrotropic and magnetic moments are well conserved away from the resonances. Thus, 3D velocity space can be reduced to 2D energy/pitch-angle space. The mapping for this space should describe 2D charged particle motion due to energy/pitch-angle jumps with the time-intervals between jumps equal to the interval between passages through the resonances. Diffusive jumps (with zero mean values) and jumps driven by nonlinear phase bunching and phase trapping depend on the resonant phase $\varphi_{R}$, i.e. a variable proportional to the particle gyrophase, which changes fast. In low wave intensity systems this phase is randomly distributed over entire ($\varphi_{R}\in[0,2\pi]$) range, and the phase dependence $\sim\sin\varphi_{R}$ can be directly included into the map (Vasilev et al., 1988; Zaslavskii et al., 1989; Benkadda et al., 1996; Khazanov et al., 2013, 2014). The phase bunching and phase trapping operate in certain $\varphi_{R}$ ranges (e.g., Albert, 1993; Itin et al., 2000; Grach & Demekhov, 2020), whereas jumps depend on $\varphi_{R}$ quite nonmonotonically (see Artemyev et al. (2014b, 2018a)). However, due to phase randomization between two successive resonances (see Appendix in (Artemyev et al., 2020b)), the phase-dependence can be reduced to a simplified determination of ranges corresponding to phase trapping $\varphi_{R}\in[0,2\pi\Pi]$ and phase bunching $\varphi_{R}\in[\Pi,2\pi]$ where $\Pi<1$ is the probability of trapping (see, e.g., Artemyev et al. (2018a)). The phase gain between two resonances is a large value depending on particle energy and pitch-angle, but this dependence can be omitted in the leading approximation (see discussion in Artemyev et al. (2020b)). Therefore, in this study we consider charged particle transport in the energy/pitch-angle space due to nonlinear resonant interaction under assumption of resonant phase randomization (limitations of this assumption have been studied in Artemyev et al. (2020a)).

The paper structure includes a description of the basic system properties and examples of multi-resonant systems observed in the Earth’s radiation belts (Sect. 1). We present three examples: with two whistler-mode waves providing two cyclotron resonances, with one oblique whistler-mode wave providing cyclotron and Landau resonances, and with one whistler-mode wave and one EMIC wave providing two different cyclotron resonances. Then we focus on the first example and construct the map for this system (Sect. 2). Theoretical results derived from this map are verified with test particle simulations. At the end of the paper we discuss the constructed map and possible extensions of the proposed approach (Sect. 3).

The Hamiltonian of a relativistic charged particle (e.g., an electron with rest mass $m_{e}$ and charge $-e$) moving in the 2D inhomogeneous magnetic field of the Earth dipole and interacting with electromagnetic waves (in the low amplitude limit with the wave energy $U_{w}$ much smaller than electron energy $\sim m_{e}c^{2}$, where $c$ is the speed of light) can be written as (e.g., Albert et al., 2013; Artemyev et al., 2018b):

	$\displaystyle H$	$\displaystyle=$	$\displaystyle m_{e}c^{2}\gamma+U_{w}\left({s,I_{x}}\right)\sin\left({\phi\pm n\psi}\right)$
	$\displaystyle\gamma$	$\displaystyle=$	$\displaystyle\sqrt{1+\frac{{p_{\parallel}^{2}}}{{m_{e}^{2}c^{2}}}+\frac{{2I_{x}\Omega_{ce}}}{{m_{e}c^{2}}}},$		(1)

where two pairs of conjugate variables are $(s,p_{\parallel})$ (the field-aligned coordinate and momentum) and $(\psi,I_{x})$ (gyrophase and momentum $I_{x}=c\mu/e$; $\mu$ is the classical magnetic moment). The electron gyrofrequency $\Omega_{ce}=eB_{0}/m_{e}c$ is determined by the background magnetic field $B_{0}(s)$, given by, e.g., the reduced dipole model (Bell, 1984). The sign $\pm$ in front of $\psi$ is determined by the wave polarization: $+$ for whistler-mode waves interacting with electrons and $-$ for EMIC waves interacting with electrons. The resonance number is $n=0,\pm 1,\pm 2...$. The wave vector ${\bf k}=(k_{\parallel}(\omega,s),k_{\perp}(\omega,s))$ is given by cold plasma dispersion equation (Stix, 1962) for a constant wave frequency $\omega$ (i.e., $\partial\phi/\partial s=k_{\parallel}$, $\partial\phi/\partial t=\omega$). For a finite angle $\theta={\rm arctan}(k_{\perp}/k_{\parallel})$ between the wave vector and the background magnetic field the wave amplitude in Hamiltonian (1) takes the form (Albert, 1993; Tao & Bortnik, 2010; Nunn & Omura, 2015; Artemyev et al., 2018b):

	$\displaystyle U_{w}$	$\displaystyle=$	$\displaystyle\sqrt{\frac{{2I_{x}\Omega_{ce}}}{{m_{e}c^{2}}}}\frac{{eB_{w}}}{k}\sum\limits_{\pm}{\frac{{\cos\theta\pm C_{1}}}{{2\gamma}}J_{n\pm 1}\left({\sqrt{\frac{{2I_{x}k^{2}}}{{m_{e}\Omega_{ce}}}}\sin\theta}\right)}$
		$\displaystyle+$	$\displaystyle\frac{{eB_{w}}}{k}\left({\frac{{p_{\parallel}}}{{\gamma m_{e}c}}+C_{2}}\right)J_{n}\left({\sqrt{\frac{{2I_{x}k^{2}}}{{m_{e}\Omega_{ce}}}}\sin\theta}\right)\sin\theta$

where $B_{w}$ is the wave magnetic field amplitude, $C_{1,2}$ are functions of wave dispersion and $\theta$, and $J_{n}$ are Bessel functions. Equation (LABEL:eq02) shows that for field-aligned waves $\theta=0$ there is only one cyclotron resonance $n=-1$: $U_{w}=\sqrt{2I_{x}\Omega_{ce}/m_{e}c^{2}}eB_{w}/k\gamma$ (with $C_{1}=1$ for $\theta=0$, see (Tao & Bortnik, 2010)). For oblique wave propagation $\theta\neq 0$ the whole set of resonances with different values of $n$ is present.

To describe the long-term evolution of electron dynamics in the energy/pitch-angle space, we propose to develop a map providing relations for each resonant interaction $\Delta\gamma=\Delta\gamma(\gamma,\alpha_{eq})$, $\Delta\alpha_{eq}=\Delta\alpha(\gamma,\alpha_{eq})$. Changes $\Delta\gamma,\Delta\alpha_{eq}$ are due to phase bunching (nonlinear scattering) and phase trapping. Thus, the first step in the construction of such a map is to derive equations for $\Delta\gamma,\Delta\alpha_{eq}$ driven by both these processes. We start with Hamiltonian (LABEL:eq05) and follow the standard procedure of Hamiltonian expansion around the resonant $I_{0},I_{1}$ values (Neishtadt, 2014; Artemyev et al., 2018a), which are defined by equations $\partial H_{I}/\partial I_{i}=0$:

$$\frac{{k_{i}I_{iR}}}{{m_{e}c}}=-\frac{{P+k_{i^{\prime}}I_{i^{\prime}}}}{{m_{e}c}}-\frac{{\Omega_{ce}}}{{k_{i}}}+\frac{1}{{\sqrt{\left({k_{i}c/\omega_{i}}\right)^{2}-1}}}\sqrt{1-\left({\frac{{\Omega_{ce}}}{{k_{i}c}}}\right)^{2}-2\frac{{\Omega_{ce}}}{{k_{i}c}}\frac{{P+\left({k_{i^{\prime}}-k_{i}}\right)I_{i^{\prime}}}}{{m_{e}c}}}$$

(12)

where $i’=0$ for $i=1$ and $i’=1$ for $i=0$. Expansion of Hamiltonian (LABEL:eq05) around $I_{i}=I_{iR}$ gives

	$\displaystyle H_{Ii}$	$\displaystyle\approx$	$\displaystyle\Lambda_{i}+m_{e}c^{2}\frac{1}{2}g_{i}\left({I_{i}-I_{iR}}\right)^{2}+u_{iR}\sin\varphi_{i}$
	$\displaystyle\Lambda_{i}$	$\displaystyle=$	$\displaystyle m_{e}c^{2}\gamma_{iR}-\left({\omega_{0}I_{0}+\omega_{1}I_{1}}\right)_{I_{i}=I_{iR}}$		(13)
	$\displaystyle\gamma_{iR}$	$\displaystyle=$	$\displaystyle\frac{{\left({k_{i}c/\omega_{i}}\right)}}{{\sqrt{\left({k_{i}c/\omega_{i}}\right)^{2}-1}}}\sqrt{1-\left({\frac{{\Omega_{ce}}}{{k_{i}c}}}\right)^{2}-2\frac{{\Omega_{ce}}}{{k_{i}c}}\frac{{P+\left({k_{i^{\prime}}-k_{i}}\right)I_{i^{\prime}}}}{{m_{e}c}}}$
	$\displaystyle u_{iR}$	$\displaystyle=$	$\displaystyle\sqrt{\frac{{2\Omega_{ce}\left({I_{0}+I_{1}}\right)_{I_{i}=I_{iR}}}}{{m_{e}c^{2}}}}\frac{e}{{\gamma_{iR}}}\frac{{B_{w,i}}}{{k_{i}}}$
	$\displaystyle g_{i}$	$\displaystyle=$	$\displaystyle\left.{\frac{{\partial^{2}\gamma}}{{\partial I_{i}^{2}}}}\right|_{I_{i}=I_{iR}}=\frac{{k_{i}^{2}}}{{m_{e}^{2}c^{2}}}$

where $\varphi_{i}$ are fast variables and $I_{i}-I_{iR})$ and $(s,P)$ are slow variables (note that $\Lambda_{i}$ does not depend on fast variables). Next, we introduce new variables $P_{\varphi i}=I_{i}-I_{iR}$ with the generating function $Q_{i}=(I_{i}-I_{iR})\varphi_{i}+s\tilde{P}_{i}$. New Hamiltonians are

	$\displaystyle F_{i}$	$\displaystyle=$	$\displaystyle\Lambda_{i}\left({\tilde{s},\tilde{P}}\right)+m_{e}c^{2}\frac{1}{2}g_{i}P_{\varphi i}^{2}+u_{iR}\sin\varphi_{i}$
		$\displaystyle\approx$	$\displaystyle\Lambda_{i}\left({s,P}\right)+\left\{{\Lambda_{i},I_{iR}}\right\}_{s,P}\varphi_{i}+m_{e}c^{2}\frac{1}{2}g_{i}P_{\varphi i}^{2}+u_{iR}\sin\varphi_{i}$

where $\tilde{s}=s-(\partial I_{iR}/\partial P)\varphi_{i}$, $\tilde{P}=p+(\partial I_{iR}/\partial s)\varphi_{i}$, $\{…\}$ are Poisson brackets, and we expand $\Lambda(\tilde{s},\tilde{P})$ over small $\partial I_{iR}/\partial s$, $\partial I_{iR}/\partial P$ terms. Hamiltonian $F_{i}$ is the sum of $\Lambda_{i}(s,P)$ describing slow variable dynamics and pendulum Hamiltonian describing fast variable dynamics:

$$F_{\varphi i}=m_{e}c^{2}\frac{1}{2}g_{i}P_{\varphi i}^{2}+\left\{{\Lambda_{i},I_{iR}}\right\}_{s,P}\varphi_{i}+u_{iR}\sin\varphi_{i}$$

(15)

where the coefficients depend on the slow variables. Figure 7 shows phase portraits of $F_{\varphi i}$ for systems with $a_{i}=|u_{iR}/\{\Lambda_{i},I_{iR}\}|<1$ (panel a) and with $a_{i}=|u_{iR}/\{\Lambda_{i},I_{iR}\}|>1$ (panel b). For low wave amplitude $a_{i}<1$ the phase portrait does not contain closed orbits, i.e., all particles cross the resonance $\dot{\varphi}_{i}=m_{e}c^{2}g_{i}P_{\varphi,i}=0$ within an interval of about one period of $\varphi_{i}$. There are only weak scatterings in this regime with zero mean changes of $I_{i}$, and such scatterings can be described by the quasi-linear diffusion model for inhomogeneous plasma (e.g., Karpman, 1974; Albert, 2010; Grach & Demekhov, 2020). For sufficiently high wave amplitude $a_{i}>1$, however, the phase portrait contains both closed and open orbits, i.e., there are now phase trapped particles oscillating around the resonance $\dot{\varphi}_{i}=m_{e}c^{2}g_{i}P_{\varphi,i}=0$ for a long time. Scattering (crossing of the resonance along the open orbits) would result in phase bunching with a small, yet nonzero mean change of $I_{i}$ (see reviews by Shklyar & Matsumoto, 2009; Albert et al., 2013, and references therein), whereas phase trapping would significantly change $I_{i}$. We would like to include this nonlinear regime of wave-particle interaction into the map in energy/pitch-angle space. For this reason, we derive expressions for changes of $I_{i}$ due to phase bunching, $\Delta_{scat}I_{i}$, and due to phase trapping $\Delta_{trap}I_{i}$. As $\Delta_{scat}I_{i}$ is local, i.e. depends on particle and system characteristics at the resonance, we can keep slow variables unchanged for $\Delta_{scat}I_{i}$ evaluations:

	$\displaystyle\Delta_{scat}I_{i}$	$\displaystyle=$	$\displaystyle 2\int\limits_{-\infty}^{t_{iR}}{\frac{{\partial H_{Ii}}}{{\partial I}}dt}=\frac{{2u_{iR}}}{{m_{e}c^{2}g_{i}}}\int\limits_{-\infty}^{\varphi_{iR}}{\frac{{\cos\varphi_{i}d\varphi_{i}}}{{P_{\varphi i}}}}$
		$\displaystyle=$	$\displaystyle\sqrt{\frac{{2u_{iR}}}{{m_{e}c^{2}g_{i}}}}\int\limits_{-\infty}^{\varphi_{iR}}{\frac{{\sqrt{u_{iR}}\cos\varphi_{i}d\varphi_{i}}}{{\sqrt{F_{\varphi i}-\left\{{\Lambda_{i},I_{iR}}\right\}_{s,P}\varphi_{i}-u_{iR}\sin\varphi_{i}}}}}$
		$\displaystyle=$	$\displaystyle\sqrt{\frac{{2u_{iR}}}{{m_{e}c^{2}g_{i}}}}\int\limits_{-\infty}^{\varphi_{iR}}{\frac{{\sqrt{u_{iR}}\cos\varphi_{i}d\varphi_{i}}}{{\sqrt{\left\{{\Lambda_{i},I_{iR}}\right\}_{s,P}\left({\varphi_{Ri}-\varphi_{i}}\right)+u_{iR}\left({\sin\varphi_{Ri}-\sin\varphi_{i}}\right)}}}}$
		$\displaystyle=$	$\displaystyle\sqrt{\frac{{2u_{iR}}}{{m_{e}c^{2}g_{i}}}}\int\limits_{-\infty}^{\varphi_{iR}}{\frac{{\sqrt{a_{i}}\cos\varphi_{i}d\varphi_{i}}}{{\sqrt{\left({\varphi_{Ri}-\varphi_{i}}\right)+a_{i}\left({\sin\varphi_{Ri}-\sin\varphi_{i}}\right)}}}}=\sqrt{\frac{{2u_{iR}}}{{m_{e}c^{2}g_{i}}}}f_{i}\left({\varphi_{Ri},a_{i}}\right)$

where $t_{iR}$ is the time of passage through the resonance, $\varphi_{iR}$ is the wave phase at this time, and we use $\dot{\varphi}_{i}=m_{e}c^{2}g_{i}P_{\varphi i}=2^{1/2}\sqrt{F_{\varphi i}-\{\Lambda_{i},I_{iR}\}\varphi_{i}-u_{Ri}\sin\varphi_{i}}$, $F_{\varphi i}=\{\Lambda_{i},I_{iR}\}\varphi_{Ri}+u_{Ri}\sin\varphi_{Ri}$ (resonant energy $F_{\varphi i}$ value evaluated at $P_{\varphi i}=0$). Note Eq. (3) describes $\Delta_{scat}I_{i}$ change for the particle motion thought the resonance from $-\infty$ to resonant $\varphi_{iR}$, whereas the motion in opposite direction would result in change of sign of $\Delta_{scat}I_{i}$. Function $f_{i}(a_{i},\varphi_{iR})$ is periodic for $\varphi_{iR}$, see Fig. 7(c). Although the sign of $f_{i}$ changes within one $\varphi_{iR}$ period, the mean value of this function for $a_{i}>1$ is not zero, providing the effect of phase bunching. To consider the precise $\Delta_{scat}I_{i}$ dependence on $\varphi_{iR}$ in the mapping, one would need to keep information about resonant phase $\varphi_{i}$ and calculate the phase gain between resonances. However, the phase is fast rotating, and even a small change of $\varphi_{i}$ at the resonance would result in a significant change of phase gain between resonances. Therefore, we can assume that $\varphi_{iR}$ is a random variable with a uniform distribution of the resonant energy $F_{\varphi i}(\varphi_{iR})$ at $P_{\varphi i}=0$ axis (see justification of this assumption in Itin et al. (2000); Artemyev et al. (2020a, b)), and all resonant particles with the same slow variables (same energy and pitch-angle) at the resonance would experience the same $\Delta I_{i}$ change equal to $\langle\Delta I_{i}\rangle$ averaged over the resonant energy (Artemyev et al., 2020b).

An important property of $f$ function from Eq. (3) is that being averaged over energies in resonance, $F_{\varphi i}=\{\Lambda_{i},I_{iR}\}\varphi_{Ri}+u_{Ri}\sin\varphi_{Ri}$, this function gives

$$\left\langle{f_{i}}\right\rangle=-\sqrt{\frac{m_{e}c^{2}g_{i}}{2u_{iR}}}\frac{S}{{2\pi}}=-\sqrt{\frac{8|u_{iR}|}{a_{i}m_{e}c^{2}g_{i}}}\int\limits_{\varphi_{i-}}^{\varphi_{iR}}{\sqrt{\left({\varphi_{Ri}-\varphi_{i}}\right)+a_{i}\left({\sin\varphi_{Ri}-\sin\varphi_{i}}\right)}d\varphi_{i}}$$

(17)

where $S$ is the area surrounded by the separatrix in the phase portrait in Fig. 7(b) (see details of Eq. (17) derivations in Neishtadt (1999) and Artemyev et al. (2018a)). Therefore, the $\Delta_{scat}I_{i}$ change due to phase bunching is equal to $-S/2\pi$ and for $a_{i}\gg 1$ (i.e. for very weak magnetic field inhomogeneity; note $\{\Lambda_{i},I_{iR}\}\sim\partial/\partial s$) we have $\Delta_{scat}I_{i}=-8\sqrt{2u_{iR}/m_{e}c^{2}g_{i}}$ where $S=16\sqrt{2u_{iR}/m_{e}c^{2}g_{i}}$ is the width of the resonance for large amplitude waves (Palmadesso, 1972; Karimabadi et al., 1990).

The change of $I_{i}$ due to phase bunching (nonlinear scattering) is sufficiently small to consider this process locally in energy/pitch-angle space, i.e., $\Delta_{scat}I_{i}\ll I_{i}$ (see discussion of exceptions for $\Delta_{scat}I_{i}\sim I_{i}$ in Appendix A), whereas the change of $I_{i}$ due to phase trapping is essentially non-local. To evaluate $\Delta_{trap}I_{i}$, we take into account that $I_{i}=I_{iR}$ in the resonance (during the trapping), the trapping time is defined as $2\pi I_{\phi}=\int{P_{\varphi i}d\varphi_{iR}}=S$, and $I_{\phi}$ is conserved during the trapping (because trapped particles oscillate in the $(\varphi_{iR},P_{\phi i})$ plane much faster than the system evolves (much faster than variations of slow variables $s,P$). Thus, the trapping time is defined as the time of arrival to the resonance $t_{R}$ with $\dot{S}(t_{R})>0$ (the growth of the area surrounded by the separatrix allows trapping of particles moving along open trajectories into closed trajectories), whereas the time $t_{R}^{*}$ of escape from the trapping is defined by $S(t_{R}^{*})=S(t_{R})$ and $\dot{S}(t_{R}^{*})<0$ (see scheme in Fig. 7(d)):

$$\Delta_{trap}I_{i}=I_{iR}\left({t_{iR}^{*}}\right)-I_{iR}\left({t_{iR}}\right),\quad S\left({t_{iR}^{*}}\right)=S\left({t_{iR}}\right),\quad\dot{S}\left({t_{iR}}\right)>0,\quad\dot{S}\left({t_{iR}^{*}}\right)<0$$

(18)

At the resonance, an electron can be scattered (i.e., experience the phase bunching) or trapped, and this depends on the $\varphi_{iR}$ value (e.g., Albert, 1993; Itin et al., 2000; Grach & Demekhov, 2018). However, as $\varphi_{iR}$ is a fast oscillating variable, we can consider the so-called probability of trapping instead of tracing the precise $\varphi_{iR}$ value: the range of $\varphi_{iR}$ of trapped particles, i.e., the ratio of trapped particles to the total number of resonant particles for a single resonance, is the probability of trapping, $\Pi_{i}$ (e.g., Arnold et al., 2006, and references therein). For small $\Pi_{i}$, this probability is defined as the ratio of the change of the area under the separatrix, $\dot{S}$, and the total resonant flux $\int_{0}^{2\pi}\dot{P}_{\varphi,i}d\phi=2\pi\{\Lambda_{i},I_{iR}\}$: $\Pi_{i}=\dot{S}/2\pi\{\Lambda_{i},I_{iR}\}=\{S,F_{i}\}/2\pi\{\Lambda_{i},I_{iR}\}$. This definition of the trapping probability has been verified for various plasma systems (e.g., Artemyev et al., 2014b; Leoncini et al., 2018; Vainchtein et al., 2018). Therefore, the resonant interaction can be characterized by $\Pi$, $\Delta_{trap}I_{i}$, and $\Delta_{scat}I_{i}$.

Due to conservation of $H_{i}=m_{e}c^{2}\gamma-\omega_{0}I_{0}-\omega_{1}I_{1}$, changes of $I_{i}$ are directly related to $\gamma$ changes, whereas the $I_{x}=I_{0}+I_{1}$ relation gives the pitch-angle change:

			$\displaystyle\frac{{\omega_{i}\sin^{2}\alpha_{eq}}}{{\gamma^{2}-1}}\left({\Delta_{i}\gamma}\right)^{2}-2\Delta_{i}\gamma\frac{{\Omega_{eq}-\gamma\omega_{i}\sin^{2}\alpha_{eq}}}{{\gamma^{2}-1}}$
		$\displaystyle+$	$\displaystyle\omega_{i}\left({\sin^{2}\left({\alpha_{eq}+\Delta_{i}\alpha_{eq}}\right)-\sin^{2}\alpha_{eq}}\right)=0$

that for small changes (phase bunching) can be rewritten as

$$\Delta_{i}\alpha_{eq}=\Delta_{i}\gamma\frac{{\Omega_{eq}-\gamma\omega_{i}\sin^{2}\alpha_{eq}}}{{\omega_{i}\sin\alpha_{eq}\cos\alpha_{eq}\left({\gamma^{2}-1}\right)}}$$

(19)