Modeling of three-dimensional betatron oscillation and radiation reaction in plasma accelerators

The new generation of accelerators, using plasma as the acceleration media, offer high acceleration gradient in the order of 10–100 $\rm GV/m$ and strong transverse focusing field [1, 2, 3]. Depends on the driver type, the plasma accelerators are named laser wakefield accelerators (LWFAs), which are driven by laser pulses, and plasma wakefield accelerators (PWFAs), which are driven by charged particle beams. When a high intensity laser pulse ($\gtrsim 10^{18}\ \rm W/cm^{2}$) or a high current particle beam ($\gtrsim 1\ \rm kA$) propagates through an underdense plasma, the radiation pressure of the laser or the space charge of the beam expels all plasma electrons away from axis radially, leaves behind a nearly uniform ion channel. This high-intensity three-dimensional (3D) regime has been referred to as the blowout regime [4, 5]. In this regime, the expelled electrons are pulled back by the ion channel and thereby bubble-like plasma wake wave is created. The wake consists of a longitudinal electric field that is a function of distance behind the driver, and transversal electromagnetic fields that are proportional to the off-axis distance. Consequently, in addition to the longitudinal acceleration / deceleration, the electrons reside in the wake also perform radial oscillation, called betatron oscillation (BO), under the action of transverse focusing field, with the frequency $\omega_{\beta}=\omega_{p}\kappa/\sqrt{\gamma}$, where $\omega_{p}$ is the plasma frequency, $\gamma$ is the relativistic factor of the electron, and $\kappa$ is the focusing constant which takes $1/\sqrt{2}$ for the blow-out regime [6, 7].

Electrons emit synchrotron radiation when performing BO [8, 9, 10], which affects the electrons in return. Such effect is called the radiation reaction (RR) and its classical expression is the Lorentz–Abraham–Dirac (LAD) equation or the Landau-Lifshitz equation [11, 12]. Because the RR force is proportional to the classical electron radius $r_{e}\approx 2.81\times 10^{-15}\ \rm m$, it is generally negligible unless under extreme conditions [13, 14] or for sufficiently long interaction time [15]. The BO in a plasma accelerator is another good case for such long interaction time. The radiation leads to the energy loss of electrons and in return affects the energy-dependent betatron frequency, as well as the other beam properties, such as the energy spread and emittance [16, 17, 18, 19, 20, 21, 22, 23].

Although there are many established theories on the long-term RR damping effect of BO, their models assert the electron only moves in one plane, spanned by the longitudinal direction and one transverse direction, thus only linear polarization is considered. Moreover, these models usually neglect the longitudinal and energy oscillations during one betatron period. In this paper, we establish a 3D BO model with RR effect, which generalizes the betatron polarization from linear to elliptical, and considers both the longitudinal and energy oscillations. Long-term equations (LTEs), without resolving the betatron period, are derived and verified by numerical methods. The LTEs reproduce the previous results, such as longitudinal and transverse damping due to RR, as the special cases, and meanwhile reveal new phenomena such as betatron phase shift and polarization change.

The rest of this paper is organized as follows. Sec. II gives the original form of the force, and shows the equations of motion expressed by the transverse motion only. Sec. III derives the LTEs by averaging the equations of motion through one betatron period. Sec. IV discusses the different phenomena in the RR dominant regime and the betatron phase shift dominant regime. Sec. V numerically verifies the LTEs by comparing with the code PTracker which solves the equations with the original form of force. Before start, it is worth noting that we use plasma normalization units described in Appx. A, and some symbols and calculation rules used often during the derivation are described in Appx. B.

Consider an electron with $\gamma\gg 1$ is trapped in a plasma wakefield with the longitudinal co-moving coordinate $\zeta=z-\beta_{w}t$, where the wake is moving in the $z+$ direction with the phase velocity $\beta_{w}$. Neglect the interaction between the beam particles, the electromagnetic field provided by the wake can be modeled as [7]

where $E_{z0}=\left.E_{z}\right|_{\zeta=\left<\zeta\right>}$, $\lambda=\left.dE_{z}/d\zeta\right|_{\zeta=\left<\zeta\right>}$, and $\vec{r}=\left(x,y\right)$ is the transverse offset. The force can be expressed as

where $\beta_{x}=\dot{x}$, $\beta_{y}=\dot{y}$, $\beta_{z}=\dot{z}=\beta_{z0}+\dot{\zeta_{1}}$, $\beta_{z0}=\beta_{w}+\dot{\left<\zeta\right>}=\left<\beta_{z}\right>$, and $\vec{f}^{\rm rad}$ is the RR force. The formulas of 3D BO with RR can be written in the form of transverse terms only (see Appx. C)

where

and $r_{e}$ is also normalized to $k_{p}^{-1}$. One may note the second terms in Eqs. (7) - (10), which come from the oscillation of $\beta_{z}$ and the modulation of $\gamma$ due to transverse oscillation, were neglected in previous works. In the following sections we show these terms lead to new phenomena.

In this section we use the same averaging method as Ref. [20]. We firstly introduce two complex variables

where

is the betatron phase. Obviously $\left|U_{1}\right|\ll\left|\left<U\right>\right|$ and $\left|V_{1}\right|\ll\left|\left<V\right>\right|$ are satisfied, and we apply the rules in Appx. B often in the following. Because the equations for $x$ and $y$ directions are symmetric, we may derive for $x$ direction only, then exchange $x$ and $y$, $U$ and $V$ for the $y$ direction. With the help of Eqs. (7), (9) and (11), we may write the time derivative of Eq. (13) as

In the following, we omit $\left<\right>$ on $U$ and $V$ for convenience, so that all $U$ and $V$ actually mean $\left<U\right>$ and $\left<V\right>$. We perform average on Eq. (16) to obtain (note only the terms with $e^{i0\varphi}$ survive after averaging)

By asserting $V=0$ and omitting the last three terms in Eq. (17), which comes from the second terms in Eqs. (7) - (10), we can reproduce Eq. (19) in Ref. [20].

The average of Eq. (7) leads to

with the second term reproduces Eq. (B2) in Ref. [18]. $E_{z0}$ is a function of $\left<\zeta\right>$, which obeys

where $\gamma_{w}=\left(1-\beta_{w}^{2}\right)^{-1/2}$ and we have used Eq. (46).

The above averaged equations Eqs. (17), (18) and (19) are already enough to predict the long-term behavior of BO. However, the equations for the complex variables are not explicit. To make them more physically meaningful, we introduce

$\left|U\right|$ has the meaning of the BO amplitude in the $x$ direction, and $\Phi_{x}$ the phase shift. For the $y$ direction they are similar. Thus $\Delta\Phi$ is the phase difference of the two directions. By Applying $d\left|U\right|/dt=\left(\dot{U}U^{*}+U\dot{U}^{*}\right)/2\left|U\right|$ and $\dot{\Phi}_{x}=\left(\dot{U}U^{*}-\dot{U}^{*}U\right)/2i\left|U\right|^{2}$ we get

Note when doing exchange of $U$ and $V$ for the $y$ direction, one also has to change the $\pm$ sign of $\Delta\Phi$.

To further simplify we notice Eq. (23) can be rewritten with the help of Eq. (18)

which reproduces Eq. (66) in Ref. [22] if $V=0$. Introduce

which have the physical meaning of the areas (divided by $2\pi$) of the ellipses encircled by the particle trajectory in $x$-$p_{x}$ and $y$-$p_{y}$ phase spaces. Then Eqs. (18), (19), (23), (24) and (25) can be rewritten as

It is generally safe to take $\beta_{z0}=1$ here. But Eq. (12), or $\beta_{z0}=1-\frac{1}{2}\left[\left<\gamma\right>^{-2}+\frac{1}{2}\kappa\left<\gamma\right>^{-3/2}\left(S_{x}+S_{y}\right)\right]$, gives a better accuracy. The above long-term equations, Eqs. (29) - (33), show that the BO experiences acceleration (for $E_{z0}<0$) or deceleration (for $E_{z0}>0$), radiation damping, longitudinal phase drift, and betatron phase shift. These equations may be used for the long-term behavior of BO without resolving the betatron period.

From Eqs. (31) - (33) we note two regimes. One is the RR dominant regime, where $r_{e}\left<\gamma\right>^{5/2}\gg 1$, so that the first terms in Eqs. (31) - (33) dominate. This regime has been discussed before [20], although only for the linearly polarized case $\Delta\Phi=0$ (so that the ratio between $S_{x}$ and $S_{y}$ is a constant). The other is the betatron phase shift dominant regime, where $r_{e}\left<\gamma\right>^{5/2}\ll 1$, so that the remaining terms in Eqs. (31) - (33) dominate. These terms were previously proposed [22], but the betatron phase shift is found for the first time in the present work.

In the RR dominant regime, an interesting phenomenon is that an elliptical polarization (in the $x$-$y$ plane) always approaches linear polarization, because $\Delta\Phi$ always approaches the nearest integer multiple of $\pi$ according to Eq. (33). This phenomenon can also be viewed by rotating the $x$ axis to the major axis of the ellipse, so that $S_{x}>S_{y}$ and $\Delta\Phi=\pi/2$. Define $R=S_{y}/S_{x}$ and perform time derivative with the help of Eq. (31)

which suggests that the ellipse monotonically becomes thinner.

The betatron phase shift dominant regime requires a moderate $\gamma$, or straightforwardly $r_{e}\rightarrow 0$, which corresponds to very dilute plasma case, leads to a constant $S\equiv S_{x}+S_{y}$. It can be proved that the time integral of Eq. (30) reproduces Eq. (6) in Ref. [21], which is the $\left<\zeta\right>$ drift, in the case that $\left<\gamma\right>$ linearly depends on $t$. In another case that $\left<\zeta\right>$ drifts around the zero point of $E_{z0}$, the drift frequency can be obtained by using Eqs. (29) and (30), and asserting $E_{z0}=\lambda\left<\zeta\right>$

We also note the angular momentum $L_{z}=\gamma x\beta_{y}-\gamma y\beta_{x}$ and its changing rate

which obeys the law of conservation of angular momentum if $\left<L_{z}\right>=0$ initially, or $r_{e}\rightarrow 0$. Especially for $r_{e}\rightarrow 0$, the particle trajectory in the $x$-$y$ plane generally encircles an ellipse with constant area and shape, which also has precession leading to the rotation of the major and minor axes of the ellipse.

To verify the LTEs, we solve them numerically using the backward-differentiation formulas (BDF) in the SciPy integration package [24]. Meanwhile, the original equations of motion with the force expressions Eqs. (4) - (6) are solved by Runge-Kutta 4th order method using the code PTracker (PT) [25]. We choose four cases with their parameters and initial values listed in Tab. 1, and the comparison results are plotted from Fig. 1 to 4. Note that $\Phi_{x}$ cannot be obtained directly from PT. Thus we perform the following treatment to the PT results

because $x=\left|U\right|\cos\left(\varphi+\Phi_{x}\right)$, where $\varphi$ is obtained by numerical integral based on Eq. (15). Then $\cos\Phi_{x}$ can be obtained by a low-pass filter. Similar treatment is performed to obtain $\cos\Delta\Phi$, according to

A case in the betatron phase shift dominant regime is shown in Fig. 1. We see $S_{x}+S_{y}$ is a constant, although $S_{x}$, $S_{y}$ and $\Delta\Phi$ change gradually. The approximate “phase-locking” is chosen, i.e. $\gamma_{w}\approx\gamma_{z0}$, thus $\left<\zeta\right>$ oscillates near the zero point of $E_{z}$ with the drift frequency $\omega_{\left<\zeta\right>}$ according to Eq. (35). $\left<\gamma\right>$ oscillates with the same drift frequency as shown in Fig. 1 (b).

A second case during the transition of the two regimes is shown in Fig. 2. $S_{x}+S_{y}$ is approximately a constant initially, and starts to decrease near the regime transition $\gamma=r_{e}^{-2/5}$.

A third case in the RR dominant regime is shown in Fig. 3. The initial values are chosen so that the particle trajectory in the $x$-$y$ plane is a ellipse with its major axis laying on the $x$ axis. As shown in Fig. 3 (a), $R=S_{y}/S_{x}$ decreases monotonically, as predicted by Eq. (34).

The last case shown in Fig. 4 is also in the RR dominant regime, but the particle trajectory in the $x$-$y$ plane is a oblique ellipse. As shown in Fig. 4 (c), $\Delta\Phi$ gradually approaches $\pi$, which is in accordance with the discussion in Sec. IV.

In all these plots, the results from PT and LTEs show agreement with high accuracy, demonstrating the correctness of LTEs. Because the BO frequency is the highest frequency in our physical process, the LTEs largely reduce the numerical complexity and meanwhile keep the long-term accuracy.

We have established a three-dimensional betatron oscillation model including radiation reaction to study the long-term behavior of an electron in laser or beam driven plasma wakefield. The original equations of motion have been expressed by the transverse oscillation terms as Eqs. (7) - (10), and then averaged in one betatron period to obtain the long-term equations Eqs. (29) - (33). The conditions of our model are $r^{2}\ll\gamma$, $r^{2}\gamma\gg 1$ and $r\gamma r_{e}/2\alpha\ll 1$, as discussed in Appx. C. Our model, on one hand, reproduces previous results such as longitudinal deceleration and transverse damping, and on the other hand reveals new phenomena such as longitudinal phase drift oscillation, betatron phase shift and betatron polarization change. Two regimes with distinct behaviors, determined by $r_{e}\gamma^{5/2}$, are discussed in Sec. IV, and are demonstrated by numerical methods in Sec. V. The numerical comparisons of the long-term equations and the original equations of motion show the high accuracy of our model. This model can be fundamental for future plasma based high-energy accelerators and colliders [26].