# Analytic model for electromagnetic fields in the bubble regime of plasma wakefield in non-uniform plasmas

###### Abstract

Based on a model of plasma wakefield in the strongly nonlinear (bubble) regime, we develop a lowest-order perturbation theory for the components of electromagnetic fields inside and outside the bubble using the assumption of small thickness of the electron sheath on the boundary of the bubble. Unlike previous models, we derive simple explicit expressions for the components of electromagnetic fields not only in the vicinity of the center of the bubble, but in the whole volume of the bubble (including areas of driving or accelerated bunches) as well as outside it. Moreover, we apply the results to the case of radially non-uniform plasma and, in particular, to plasma with a hollow channel. The obtained results are verified with 3D particle-in-cell (PIC) simulations which show good correspondence to our model.

## Introduction

Currently, plasma acceleration methods are considered promising for obtaining high-energy electron bunches. Esarey, Schroeder, and Leemans (2009); Kostyukov and Pukhov (2015) Their main idea is the use of a driver (a relativistic electron bunch Rosenzweig et al. (1988) or an intense laser pulseTajima and Dawson (1979)) for excitation of plasma wakefield whose longitudinal electric field can be used to accelerate charged particles. The achievable accelerating gradients can be of orders of magnitude higher than those in conventional accelerators. In experiments with the laser-wakefield acceleration (LWFA), electrons with the energy of have been obtained at the acceleration length of .Leemans et al. (2014) For plasma-wakefield acceleration (PWFA), in which a relativistic bunch is used as the driver, the possibility of energy doubling from to has been experimentally demonstrated.Blumenfeld et al. (2007)

Particularly interesting is the so-called “bubble” or “blow-out” regime of plasma wakefield in which plasma electrons are almost completely expelled by the driver, leading to the formation of a spherical cavity (or a “bubble”) behind the driver.Pukhov and Meyer-ter-Vehn (2002) On the boundary of this cavity, a thin electron sheath composed of the expelled electrons is formed. The densities and currents in this sheath screen the surrounding plasma from the wakefield, so its structure is important for processes happening on the boundary of the plasma bubble and outside it. However, no self-consistent theory of the bubble regime which can predict the full structure of the wakefield in the bubble regime based on the properties of the driver exists.

Due to its complex nonlinear nature, the bubble regime is commonly studied numerically with the use of the particle-in-cell (PIC) method. Pukhov (2016) However, analytic description of this regime is also of substantial interest. Different theoretical approaches have been used to make significant advancement in the understanding of this regime, which have led to the creation of phenomenological models of the bubble,Kostyukov, Pukhov, and Kiselev (2004); Kostyukov et al. (2009) the development of a similarity theory Gordienko and Pukhov (2005) and an analytic model of the bubble’s boundary.Lu et al. (2006) Following the discovery of the advantages of plasmas with hollow channels in the bubble regime,Pukhov et al. (2014) this analytic model has recently been generalized to describe plasmas with non-uniform transverse density profile.Golovanov et al. (2016a); Thomas et al. (2016)

The aforementioned theories focus mostly on the shape of the bubble and on the dynamics of relativistic accelerated particles, while paying less attention to the structure of electromagnetic field components inside the bubble and at its boundary. Knowledge of this structure can be of significant interest for the processes of particle injection and self-injection into the wakefield,Kostyukov et al. (2009); Yi et al. (2013) when the particles under consideration are not ultra-relativistic. In principle, the approach which can be used to calculate the spatial distributions of the fields is presented by Yi et al. (2013) However, this approach involves extensive numerical calculation of derivatives and integrals. In the current paper, we make an additional step and develop a perturbation theory based on the assumption that the ratio between the thickness of the electron sheath at the boundary of the bubble and the size of the bubble is small. This assumption allows to significantly simplify all expressions in the lowest order with respect to this ratio. The novelty of our work compared to the model by Yi et al. (2013) lies in the following three points: (i) we obtain simple explicit expressions for the components of the electromagnetic field both inside and outside the bubble; (ii) we do not make an a priori assumption about the shape of the electron sheath and consider the most general case; (iii) our results are applied to the case of radially non-uniform plasmas.

In section I, we write out basic equations for the potentials and electromagnetic fields in any kind of plasma wakefield. In section II, we introduce a model for the bubble regime which assumes general structure of the electron sheath on the bubble’s boundary. In section III, we make an assumption that the electron sheath is thin compared to the size of the bubble. This allows us to obtain a simpler equation for the boundary of the bubble. In section IV, based on the aforementioned assumption, we develop a perturbation theory for the components of the electromagnetic field inside and outside the bubble and find the distributions of all these components. All derivations are accompanied and verified by PIC simulations for plasma with a hollow channel which show excellent agreement with our model.

## I General equations

We consider boundless plasma in which a laser pulse or an electron driver propagates along the axis and excites plasma wakefield in it. The plasma is assumed to be non-uniform perpendicular to the -axis. Under the assumption of axial symmetry, the electron density depends only on the coordinate in cylindrical geometry. For simplicity, we use dimensionless units in which all charges are normalized to , densities to , time to , coordinates to , momenta and energies to and , respectively, and electric and magnetic fields to . Here, is the elementary charge, is the electron mass, is the speed of light in vacuum, is typical plasma density (for example, for the case of a plasma with a hollow channel it may be the density outside the channel), is a typical electron plasma frequency.

Plasma wakefield is propagating with the velocity determined by the velocity of the driver and is close to the speed of light (equal to in dimensionless units). As the shape of the wakefield changes slowly during its propagation through plasma, we use the quasi-stationary approximation and assume that electromagnetic fields in the wakefield depend on the longitudinal coordinate and time through their combination . Plasma fields are described using the vector potential and the wakefield potential , where is the scalar potential. The electromagnetic fields can be easily retrieved from the potentials,

(1) | |||

(2) |

All other components of the electromagnetic field are equal to zero due to axial symmetry.

We use the Lorenz gauge for the potentials

(3) |

which allows us to express through the wakefield potential

(4) |

thus leaving only and as independent potentials. The combined Maxwell’s equations for these potentials under the Lorenz gauge can be separated

(5) | |||

(6) |

where is the longitudinal current density, is the charge density. By integrating these equations we get

(7) | |||

(8) |

Here, we introduce . The equation for is integrated only once because its radial derivative is sufficient to retrieve the fields in Eqs. (1), (2).

So, if we know the spatial distributions of and in the wakefield, we can easily find spatial distributions of the electromagnetic field components. So far, all equations can be applied to any kind of wakefield. Next, we introduce phenomenological models of these two sources for the bubble regime of wakefield.

## Ii Model of the bubble

Based on the properties of the bubble regime observed in the PIC simulations, we introduce a model of a bubble in which we assume that no plasma electrons are present inside the bubble, and there is a thin electron sheath on its boundary determined by a function . Far outside the bubble, the plasma is assumed non-perturbed. In this case, the source is modeled asGolovanov et al. (2016a)

(9) |

Inside the bubble, , where is the ion charge density. For relativistic electron bunches in the bubble, the longitudinal velocity is close to (the dimensionless speed of light), thus , so they do not contribute to . The function describes the shape of the electron sheath at the bubble’s boundary. We assume that it is normalized, so that and , where

(10) |

are generalized moments of the function . The parameter determines the typical width of the electron sheath. In the previous works, exponential Yi et al. (2013) and rectangular Lu et al. (2006); Thomas et al. (2016) profiles have been used.

The function is obtained using the condition and Eq. (8)

(11) |

where

(12) |

The same answer for can be obtained from the continuity equation and the assumption that when .Thomas et al. (2016)

A similar model is used for the longitudinal current . Inside the bubble, is the current created by relativistic electron bunches (either driver or witness bunches). Thus,

(13) |

As shown by Yi et al. (2013), a different typical width is required to describe the longitudinal current distribution. A different from shape function is also introduced for the purpose of generality. This function is normalized in the same way as with , where are introduced similarly to Eq. (10). The value of can be obtained from the condition , which will be demonstrated in Sec. IV.

So, for any specific shape of the bubble , the models (9), (13) allow us to calculate all of the field distributions. As shown by Golovanov et al. (2016a), model (9) and electron motion equations lead to the equation for the bubble’s boundary

(14) |

where , ,

(15) |

is the source term created by relativistic electron bunches inside the bubble, and the coefficients , , are determined by the shape of the electron sheath and its width . Despite having no laser-related terms, this equation can be applied for the case of a laser driver in regions where the laser pulse is not present. It is also possible to include the influence of the laser field on the boundary of the bubble into this equation, which is not considered in the current paper; for uniform plasma it is done in Vieira et al.Vieira, Fonseca, and Silva (2016)

If we take , , our theory up to this point is completely equivalent to the model used by Yi et al. (2013) However, despite the fact that all previous equations give an analytical model for the fields inside and outside the bubble, the expressions for the field distributions are too complex, so they can be obtained only numerically. In the subsequent sections, we show that it is possible to obtain significantly simpler expressions by making an assumption about the width of the electron sheath.

## Iii Thin sheath approximation

As shown in the previous works,Lu et al. (2006); Thomas et al. (2016); Golovanov et al. (2016a) Eq. (14) for the bubble’s boundary can be significantly simplified by assuming that the width of the electron sheath on the boundary is significantly smaller than the size of the bubble, i. e. . If we also assume that , which corresponds to a sufficiently large bubble, Eq. (14) transforms into

(16) |

where is given by Eq. (12). The coefficients of this equation are determined solely by the plasma profile , while the source on the right-hand side reflects the influence of relativistic electron bunches inside the bubble.

The applicability of this model has already been verified by PIC simulations.Thomas et al. (2016); Golovanov et al. (2016b) To demonstrate this again and to have numerical results for further analysis, we have carried out simulations of a wakefield generated by a relativistic electron bunch using the open-source 3D PIC code Smilei. Smi (); Derouillat et al. (2017) This code operates in the same dimensionless values as our theory, so all initial conditions were set dimensionless. In these simulations, we used a driver (an electron bunch) with parabolic density profile, transverse size , longitudinal size , charge of , and electrons’ energy of . High electron energy was chosen to prevent the evolution of the driver. The surrounding plasma had a hollow channel of radius , corresponding to the electron density of , where is the Heaviside step-function. For comparison, similar simulations have been performed using the PIC code Quill developed in our group. Nerush and Kostyukov (2010) The comparison of the resulting bubble with the model given by Eq. (16) is shown in Fig. 1. As expected, this model produces overall good results which slightly diverge from the actual shape of the bubble when becomes close to and thus approximations used in our model are invalid. The results of the simulations are generally similar in the two codes we have used, so all following figures will be made using only the results from the simulations with Smilei.

Interestingly enough, unlike Eq. (14), Eq. (16) holds no dependence on the width of the electron sheath or its shape whatsoever. Encouraged by this fact, we hope that a similar level of simplification can be achieved for the components of the electromagnetic field. So, in all further calculations, we assume that , corresponding to . We also make a similar assumption regarding the electron sheath current profile: . Using these two assumptions and model , we rigorously calculate distributions of fields both inside and outside the bubble in the lowest possible order with respect to and .

## Iv Fields distributions

### iv.1 Wakefield potential

In order to calculate the fields, we first need to obtain the distribution of the wakefield potential given by Eq. (8). Equations (8), (9), (11) together with the condition give us

(17) |

The second integral is equal to in the first order with respect to . Thus, wakefield potential (8) inside the bubble is

(18) |

The second term can be neglected for but is required for to stay continuous at .

### iv.2 Longitudinal electric field

Expressions (18), (19) for the potential are sufficient to obtain the longitudinal field by using Eq. (2). For the field inside the bubble, we differentiate only the first term in Eq. (18) and get:

(20) |

which is a well-known formulaGolovanov et al. (2016a); Thomas et al. (2016) for the field inside the bubble. This field does not depend on the radial coordinate because depends only on .

Differentiating Eq. (19) and leaving only the zeroth-order term in leads to

(21) |

for the longitudinal field component outside the bubble . It vanishes to zero according to . For example, if , then , so the electric field shows exponential decay with the typical width .

In order to compare the analytical results to the simulations, we need to determine the function and the thickness of the electron sheath. Based on the results of the simulations corresponding to Fig. 1, we assume that . In this case, gives the best fit to the simulations. The comparison of the simulations results to analytical expressions (20), (21) is shown in Fig. 2. Figure 3 also shows longitudinal and transverse slices of the field distribution. This comparison shows that our model adequately describes both the distribution of the electric field inside the bubble—which is indeed uniform in the transverse direction—and the way the electric field vanishes to zero outside the bubble. As already mentioned, difference is seen only at the front edge of the bubble where our model is not applicable.

### iv.3 Azimuthal magnetic field

The magnetic field in our model has only the azimuthal component . According to Eqs. (1), (2), (4), (7), it can be obtained by solving an equation

(22) |

By setting , we can get

(23) |

which can be easily integrated. Inside the bubble (for ), , thus integrating the first term gives us by its definition (15). As the longitudinal electric field inside the bubble, according to Eq. (20), does not depend on , integration of the second term is reduced to multiplication by . Thus,

(24) |

To simplify the results, we have also used Eq. (16).

Using a similar approach, it is straightforward to calculate the magnetic field for . The result is

(25) |

where

(26) |

Obviously, . So, for the longitudinal coordinates not located outside electron bunches (), the azimuthal component of the magnetic fields grows linearly with . This linear behavior does not depend on the radial plasma profile. Also, the sign of is always negative.

In order to calculate outside the bubble (), we need to use model (13) for and Eq. (21) for and substitute them into Eq. (22). The integration results in

(27) | |||

(28) |

where . The value of can be found from the condition and is

(29) |

Finally,

(30) |

This equation describes how converges to zero for , and together with Eq. (25) they give the complete spatial distribution of the azimuthal magnetic field. One more interesting result following from Eq. (29) is that the longitudinal current changes its sign at the points where , which has been observed in the simulations by Yi et al. (2013)

In order to perform comparison to the simulations, we have assumed and set to provide the best fit to the results of the simulations. The comparison is shown in Figs. 4 and 5. It again shows the validity of our approximate theory for the description of the fields both inside and outside the bubble. The longitudinal perturbations of the magnetic field in the simulations are not physical but numerical. There are two numerical sources of these perturbations: numerical Cherenkov radiation from the relativistic driverGodfrey (1974) and reflections from the simulation box boundaries. The box used in simulations is large enough ( in all dimensions) to mitigate the reflections. However, numerical Cherenkov radiation of ultrarelativistic electrons cannot be mitigated by changing the simulation parameters; methods for its suppression are not available yet in the Smilei 3D PIC code.

### iv.4 Radial electric field

After we have calculated the azimuthal component of the magnetic field, it is a very simple task to calculate the radial electric field as, according to Eq. (2), it is obtained from by adding , and has already been calculated in Sec. IV.1. Thus, for , we obtain

(31) |

where is given by Eq. (25), and for , we have

(32) |

Inside the bubble, and have different signs. For uniform plasma, grows linearly with and is typically larger than , thus has a different sign than and also shows linear dependence on , which is well-known from previous models. Kostyukov, Pukhov, and Kiselev (2004) However, for plasmas with a hollow channel, inside the channel, and therefore can have different signs inside and outside the channel. This property of the radial electric field in a hollow channel has been observed by Pukhov et al. (2014) Also, the radial component of the electric field, like , turns to zero at the boundary of the bubble when .

### iv.5 Radial force

The Lorentz force acting on accelerated particles in the bubble is important to understand their dynamics. General distributions of , , and are sufficient to find the force acting on any particle. However, there is a special case of relativistic particles which move alongside the wakefield with . This special case is typical for particles accelerated in the bubble and for the driver; for this case, the Lorentz force is significantly simplified. Kostyukov, Pukhov, and Kiselev (2004) The longitudinal component of the Lorentz force acting on such an electron is determined only by the longitudinal field distribution, , while the transverse force

(33) |

These forces depend only on the distribution of the wakefield potential .

As we have already calculated the distributions of and , the distribution of is calculated simply by taking the difference between them. The force is

(34) |

This force is always focusing for electrons (and defocusing for positrons) and turns to zero inside a hollow channel where .

The comparison of analytical and numerical distributions is shown in Figs. 8 and 9. Unlike distributions of and , the distribution of obtained from PIC simulations shows almost perfect correspondence to the one predicted by our model, which further confirms that the longitudinal perturbations visible in Figs. 4 and 6 are transverse electromagnetic waves radiated by a driver via the mechanism of numerical Cherenkov radiation.

## Conclusions

In this paper, we have considered an analytical model of a bubble excited by a relativistic electron bunch or a laser pulse in transversely non-uniform plasma. Special attention in our theory is paid to the electron sheath formed at the boundary of the bubble. The distributions of densities and currents in this sheath is important for the properties of the electromagnetic field outside the bubble. Unlike most previous theories, we use a more general model for this electron sheath. Using the assumption that the electron sheath is thin compared to the size of the bubble, we have developed a lowest-order perturbation theory and obtained explicit expressions for all components of the electromagnetic field both inside and outside the bubble. 3D particle-in-cell simulations with two PIC codes (Smilei and Quill) confirm the validity of our model. Noticeable difference between our model and the simulations is observed only in the small areas near the front and the back edges of the bubble, where our approximate theory is not applicable. Hence, the obtained simple analytic expressions can be used to describe the Lorentz force acting on charged particles in the wakefield and therefore their motion. The dynamics of charged particles can be of considerable interest for the processes of self-injection and external injection which play important role in determining the properties of the accelerated bunches.

Despite the fact that our theory is able to correctly predict the shape of the bubble and fields both inside and outside it, it still requires knowledge of several parameters from simulations, including the profile and the width of the electron sheath and the size of the bubble. The possibility of creation of a self-consistent theory which is able to predict the full structure of the wakefield based on the properties of the driver is still an important task for further investigation. One more important phenomenon our model does not consider are the so called bow waves Esirkepov, Kato, and Bulanov (2008); Luo et al. (2016) which are streams of electrons detaching from the bubble. This phenomenon is visible in our simulations in Fig. 1, but it does not have a significant influence on the fields spatial distribution for the parameters used in the simulations. However, the possibility of this phenomenon being more important under different set of parameters is a topic for future consideration.

###### Acknowledgements.

This work has been supported by the Russian Science Foundation through Grant No. 16-12-10383.## References

## References

- Esarey, Schroeder, and Leemans (2009) E. Esarey, C. B. Schroeder, and W. P. Leemans, Rev. Mod. Phys. 81, 1229 (2009).
- Kostyukov and Pukhov (2015) I. Yu. Kostyukov and A. M. Pukhov, Phys.-Usp. 58, 81 (2015).
- Rosenzweig et al. (1988) J. B. Rosenzweig, D. B. Cline, B. Cole, H. Figueroa, W. Gai, R. Konecny, J. Norem, P. Schoessow, and J. Simpson, Phys. Rev. Lett. 61, 98 (1988).
- Tajima and Dawson (1979) T. Tajima and J. M. Dawson, Phys. Rev. Lett. 43, 267 (1979).
- Leemans et al. (2014) W. P. Leemans, A. J. Gonsalves, H.-S. Mao, K. Nakamura, C. Benedetti, C. B. Schroeder, Cs. Tóth, J. Daniels, D. E. Mittelberger, S. S. Bulanov, J.-L. Vay, C. G. R. Geddes, and E. Esarey, Phys. Rev. Lett. 113, 245002 (2014).
- Blumenfeld et al. (2007) I. Blumenfeld, C. E. Clayton, F.-J. Decker, M. J. Hogan, C. Huang, R. Ischebeck, R. Iverson, C. Joshi, T. Katsouleas, N. Kirby, W. Lu, K. A. Marsh, W. B. Mori, P. Muggli, E. Oz, R. H. Siemann, D. Walz, and M. Zhou, Nature 445, 741 (2007).
- Pukhov and Meyer-ter-Vehn (2002) A. Pukhov and J. Meyer-ter-Vehn, Appl. Phys. B 74, 355 (2002).
- Pukhov (2016) A. Pukhov, CERN Yellow Rep. 1, 181 (2016).
- Kostyukov, Pukhov, and Kiselev (2004) I. Kostyukov, A. Pukhov, and S. Kiselev, Phys. Plasmas 11, 5256 (2004).
- Kostyukov et al. (2009) I. Kostyukov, E. Nerush, A. Pukhov, and V. Seredov, Phys. Rev. Lett. 103, 175003 (2009).
- Gordienko and Pukhov (2005) S. Gordienko and A. Pukhov, Phys. Plasmas 12, 043109 (2005).
- Lu et al. (2006) W. Lu, C. Huang, M. Zhou, M. Tzoufras, F. S. Tsung, W. B. Mori, and T. Katsouleas, Phys. Plasmas 13, 056709 (2006).
- Pukhov et al. (2014) A. Pukhov, O. Jansen, T. Tueckmantel, J. Thomas, and I. Yu. Kostyukov, Phys. Rev. Lett. 113, 245003 (2014).
- Golovanov et al. (2016a) A. A. Golovanov, I. Yu. Kostyukov, A. M. Pukhov, and J. Thomas, Quantum Electron. 46, 295 (2016a).
- Thomas et al. (2016) J. Thomas, I. Yu. Kostyukov, J. Pronold, A. Golovanov, and A. Pukhov, Phys. Plasmas 23, 053108 (2016).
- Yi et al. (2013) S. A. Yi, V. Khudik, C. Siemon, and G. Shvets, Phys. Plasmas 20, 013108 (2013).
- Vieira, Fonseca, and Silva (2016) J. Vieira, R. A. Fonseca, and L. O. Silva, CERN Yellow Rep. 1, 79 (2016).
- Golovanov et al. (2016b) A. A. Golovanov, I. Yu. Kostyukov, J. Thomas, and A. Pukhov, Phys. Plasmas 23, 093114 (2016b).
- (19) https://smileipic.github.io/Smilei/.
- Derouillat et al. (2017) J. Derouillat, A. Beck, F. Pérez, T. Vinci, M. Chiaramello, A. Grassi, M. Flé, G. Bouchard, I. Plotnikov, N. Aunai, J. Dargent, C. Riconda, and M. Grech, (2017), arXiv:1702.05128 [physics.plasm-ph] .
- Nerush and Kostyukov (2010) E. N. Nerush and I. Yu. Kostyukov, Probl. Atom. Sci. Tech. 4, 3 (2010).
- Godfrey (1974) B. B. Godfrey, J. Comp. Phys. 15, 504 (1974).
- Esirkepov, Kato, and Bulanov (2008) T. Zh. Esirkepov, Y. Kato, and S. V. Bulanov, Phys. Rev. Lett. 101, 265001 (2008).
- Luo et al. (2016) J. Luo, M. Chen, G.-B. Zhang, T. Yuan, J.-Y. Yu, Z.-C. Shen, L.-L. Yu, S.-M. Weng, C. B. Schroeder, and E. Esarey, Phys. Plasmas 23, 103112 (2016).