The Weibel instability beyond Weibel's Bi-Maxwellian Anisotropy
aa r X i v : . [ phy s i c s . p l a s m - ph ] F e b The Weibel instability beyond Weibel’s Bi-Maxwellian Anisotropy
T. Silva, ∗ B. Afeyan, † and L. O. Silva ‡ GoLP/Instituto de Plasmas e Fus˜ao Nuclear, Instituto Superior T´ecnico,Universidade de Lisboa, 1049-001 Lisbon, Portugal Polymath Research Inc., 94566 Pleasanton, CA, USA (Dated: February 11, 2021)The shape of the anisotropic velocity distribution function, beyond the realm of strict Maxwellianscan play a significant role in determining the evolution of the Weibel instability dictating the dy-namics of self-generated magnetic fields. For non-Maxwellian distribution functions, we show thatthe direction of the maximum growth rate wavevector changes with shape. We investigate differentlaser-plasma interaction model distributions which show that their Weibel generated magnetic fieldsmay require closer scrutiny beyond the second moment (temperature) anisotropy ratio characteri-zation.
The anisotropy of typical plasma velocity distributionfunctions can fuel instabilities. The Weibel instability, inparticular, is associated with some of the most strikingastrophysical phenomena, such as gamma-ray bursts, col-lisionless shocks, and the universe’s magnetogenesis [1–5].It can play a crucial role in laser plasma interactions aswell since it can directly couple non-local heat transport,parametric instabilities and spontaneous magnetic fieldgeneration. Weibel was the first to discover this instabil-ity in Maxwellian plasmas with anisotropic temperatures[6]. While in these same conditions, plasmas are stableto electrostatic perturbations, Weibel demonstrated thatsuch plasmas contain unstable electromagnetic modes forsufficiently small wavevectors. These modes do not prop-agate with time, as they have a pure imaginary frequency(they either grow or damp).The free energy source for the instability is typi-cally attributed to and characterized by an anisotropyof the parallel and perpendicular temperatures assum-ing Maxwellians: current density perturbations in theplasma interact with magnetic fields perturbations; if thecurrent density perturbation is in the hotter directionand the wavevector of the magnetic fields in the colder,these fields act to reinforce the current perturbation, thusgenerating an unstable feedback loop which leads to thespontaneous growth of magnetic fields. Here we challengethis attribution by showing that even for equal temper-atures but different shapes in the parallel and perpen-dicular directions, Weibel unstable modes do arise. Wealso show that for widely different shaped electron dis-tribution functions oblique modes dominate beyond thedictates of conventional theory.Due to its important potential role in high energy den-sity plasmas and its connection to astrophysics, therehas been renewed recent effort to measure it in the lab-oratory. Typically, counter-propagating plasma flows ∗ [email protected] † [email protected] ‡ [email protected] [7, 8] or laser-plasma interactions [9–12] are used to drivethe instability. Frequently, laboratory and astrophysi-cal plasmas operate in collisionless regimes, making sig-nificant deviations from Maxwellian distributions ratherlikely. Many astrophysical plasmas, such as the solarwind, are known for having hot non-thermal tails [13].Laser-plasma interactions are also known for generat-ing hot [14] or depleted tails [15, 16]; the latter had theWeibel instability and its characteristics extensively stud-ied in connection with fusion [17–21]. Non-Maxwelliandistributions are rarely considered in the literature whenperforming kinetic analysis. Nevertheless, distributionfunction details can be essential to estimate growth ratesfor the Weibel [22] and other kinetic instabilities [23–26].In this Letter, we extend the treatment of the Weibelinstability by exploring previously unidentified and un-expected features. Our results enable a deeper funda-mental understanding of the instability mechanism. Togo beyond standard theory, we consider a more genericspectrum of wavevectors and distribution functions. Wefocus our attention on laser-plasma interaction generateddistribution functions and use particle-in-cell simulationsto confirm theoretical predictions and to peer beyond thelinear theory regime.In the kinetic theory for the Weibel instability, one tra-ditionally assumes Maxwellian velocity distribution func-tions (VDFs) with different temperatures in distinct di-rections, i.e., f ( v x , v y , v z ) = f x ( v x ) f y ( v y ) f z ( v z )= f Max x ( v x ) f Max y ( v y ) f Max z ( v z ) , (1)where f Max i ( v i ) = (2 πT i ) − / e − v i / T i and T i is the tem-perature along the i th direction in units of m e c , where m e is the electron mass and c is the speed of light in vac-uum. Considering immobile ions, the dispersion relationfor electromagnetic modes with k = k x ˆ x and k = k y ˆ y are k x − ω + 1 + Z ∂f x /∂v x v x − ω/k x dv x Z v y f y dv y = 0 , (2a) k y − ω + 1 + Z v x f x dv x Z ∂f y /∂v y v y − ω/k y dv y = 0 . (2b)Here, v is given in units of c , ω in units of the electronplasma frequency ω p = (4 πn e e /m e ) / ( e and n e arethe elementary charge and the plasma density), and k in units of ω p /c . Assuming the VDF in Eq. (1) and T x > T y = T z , solutions of Eq. (2a) are always dampingmodes for any k x , while Eq. (2b) has growing modes for k y < ( T x /T y −
1) [27]. Additionally, it is guaranteed thatthe wavevector with highest growth rate has k = k y ˆ y ,i.e., one of the solutions of Eq. (2b).We notice that Eq. (2b) dependence on the hot direc-tion comes only as the temperature (we are generalizingthe concept of temperature as the second moment of anydistribution), so one could imagine that replacing f x byother distribution would lead to identical results as for aMaxwellian of the same temperature. While that is truefor Eq. (2b), we can no longer guarantee that there areno growing solutions of Eq. (2a) or that the maximumgrowth rate is a solution of Eq. (2b). The wavevectorwith highest growth rate does not require to be alignedwith either axis.One can derive a more general dispersion relation bysolving the Vlasov-Maxwell system using the method ofcharacteristics. We consider an initially field free plasma.The solutions have the form D · E = 0 [28], where E isthe electric field, D ij = k i k j − k δ ij + ω ε ij , δ ij is theKronecker delta, and ε ij is the dieletric tensor ω ε ij = (cid:0) ω − (cid:1) δ ij + k x Z d v v i v j ( ∂f /∂v x ) ω − k x v x − k y v y + k y Z d v v i v j ( ∂f /∂v y ) ω − k x v x − k y v y , (3)where we assume that k = k x ˆ x + k y ˆ y . We now fo-cus on a specific class of VDFs, namely f ( v x , v y , v z ) = f x ( v x ) f Max y ( v y ) f Max z ( v z ), with f x ( v x ) = f x ( − v x ) to guar-antee current neutrality, and T y = T z = T ⊥ . Underthese assumptions, we note that the temperature ten-sor T ij ∝ R d vv i v j f is diagonal, which guarantees thatthe phenomena observed henceforth is not due to a par-ticular choice of axes. This type of distribution is ex-tremely relevant for high-energy-density physics, wheresome laser-plasma interaction processes modify the VDFprimarily in one direction. The transverse distribution v y v z isotropy justifies the choice of k , as it is possibleto change coordinates such that k z = 0. More generally,there will be a continuum of transverse wavevectors thatgrow. As a consequence, the results will look more com-plicated, but they should be a superposition of the modesthat will be studied here. Non-trivial solutions of the Vlasov-Maxwell systemlead to the dispersion relation and require det( D ) = 0,i.e., ( D xx D yy − D xy D yx ) D zz = 0 , (4)as one can verify that D xz = D yz = D zx = D zy ≡ f ansatz in Eq. (3) and performing in-tegration over v y and v z , the relevant components of D are D xx = ω − k y + T x T y − − k x k y p T y Z dv x v x f ′ x Z ( ξ )+ 1 T y Z dv x v x f x ξZ ( ξ ) , D xy = D yx = k x k y + k x k y − k x k y Z dv x v x f ′ x ξZ ( ξ )+ s T y Z dv x v x f x ξ [1 + ξZ ( ξ )] , D yy = ω − k x − k x k y p T y Z dv x f ′ x ξ [1 + ξZ ( ξ )]+ 2 Z dv x f x ξ [1 + ξZ ( ξ )] , D zz = ω − k x − k y + T z T y − − k x k y T z p T y Z dv x f ′ x Z ( ξ )+ T z T y Z dv x f x ξZ ( ξ ) , where ξ = ω − k x v x p T y k y and Z ( ξ ) = π − / Z e − v v − ξ dv. We drastically reduce the computational costs and nu-merical errors when solving Eq. (4) by using the ansatz where f y ( v y ) and f z ( v z ) are Maxwellian. Our theory alsoreproduces the results from Ref. [29] for separable VDF’sand pure wavevector. Nevertheless, we will demonstratethat general solutions of the Weibel instability must in-clude mixed modes. Considering these modes not onlyis necessary to estimate maximum growth rates for somenon-Maxwellian distributions, but also sheds light on themechanism behind the instability.Equation (4) has two kinds of solutions D xx D yy − D xy D yx = 0 , (5a) D zz = 0 . (5b)Equation (5b) represents solutions with E = E z ˆ z , and,due to the electromagnetic nature of the Weibel insta-bility, unstable solutions of this equation results in thegrowth of magnetic field components B x and B y . Anal-ogously, unstable solutions of Eq. (5a) results in thegrowth of the B z component. FIG. 1. (a) Growth rate rescaled by the √ T , (b) wavevec-tor magnitude, and (c) wavevector angle for the mode withhighest growth rate as a function of the DLM exponent m for different values of T (symbols). In (a), the solid line isEq. (7). (d) Maximum growth rate for m = 2 and m = 3 asfunction of T x for fixed T ⊥ = 0 . Henceforth, we use two example VDFs to highlight fun-damental news aspects when using our theory. They aresuper-Gaussians (or DLM [23]) and Maxwellians with hottails, i.e., f x ( v x ) = A m exp (cid:20) − B m T m/ x | v x | m (cid:21) , (6a) f x ( v x ) = 1 − δn (2 πT c ) / exp (cid:20) − v x T c (cid:21) + δn (2 πT h ) / exp (cid:20) − v x T h (cid:21) , (6b)where A m and B m are such that R v x f x dv x = T x and R f x dv x = 1, T c and T h are the cold and hot popu-lations temperatures, and δn is the fraction of parti-cles in the tail. These two distributions could resultfrom laser-plasma interaction phenomena such as inverse bremsstrahlung [16] and stimulated Raman scattering[30], respectively.We first explore DLM distributions in cases where T x = T ⊥ ≡ T . Under this assumption, Equation (2b)has only damping solutions. Numerical analysis showsthat this is also true for Eq. (2a) for any exponent m ,i.e., pure k x and k y modes are always stable. Unexpect-edly, allowing mixed modes unlocks growing solutions ofEq. (5a) [Eq. (5b) has only damping modes]. Figures1(a-c) characterize the solution of Eq. (5a) with highestgrowth rate as function of m . Figure 1(a) verifies thatthe growth rate γ is proportional to √ T and increasesmonotonically with m . We find a engineering formulafor the growth rate as function of T and m to be γω p = 4 . × − (cid:18) T (cid:19) tanh [0 . m − / ] , (7) FIG. 2. For a DLM with m = 4 and T = 1 keV. (a) Theo-retical growth rate (solutions of Eq. (5a)) for a wide range of k x and k y . (b) B z magnetic field component in Fourier spacetaken from simulations during the Weibel instability linearstage. (c) Energy in all components of the electromagneticfields taken from simulations. (d) Saturated B z in configura-tion space. which is the line in Fig. 1(a). Figures 1(b) and (c) dis-play the wavevector magnitude | k | = ( k x + k y ) / andangle θ k = tan − ( k y /k x ), respectively. We notice that | k | and θ k are independent of T . For increasing m , thewavenumber and angle vary from 0 and π/ ∼ . ω p /c and ∼ π/
9, respectively.In Fig. 1(d), we fix T ⊥ = 0 . T x . The dashed-gray line rep-resents standard Weibel theory approximate solution forlow anisotropy ( A ≡ T hot /T cold − ≈
0, where T hot/cold is the higher/lower of T x and T ⊥ ), which predicts no in-stability when T x = T ⊥ . This trend is only observed inthe m = 2 (Maxwellian) case. For m = 3, we noticethat the growth rate remains appreciable for all valuesof T x . Thus, the growth rate could differ by several or-ders of magnitude depending on the VDF shape. When A ≈
0, it is fundamental to consider the VDF shape toaccurately predict growth of the Weibel instability. Forsufficiently large anisotropy, the maximum growth ratebecomes identical to the Maxwellian case and, if A ≫ γ/ω p = ( T hot /m e c ) / scaling [27].To confirm our theoretical predictions, we performedtwo-dimensional particle-in-cell simulations using theOSIRIS framework [31]. Simulation details are in Ref.[32]. We compare theory and simulations in a case where m = 4. Figure 2(a) shows the growth rate predicted byEq. (5a) for a wide range of k . Figure 2(a) can be di-rectly compared with Fig. 2(b), the magnetic field inFourier space in the simulation at t = 3000 ω − p , i.e., dur-ing the linear stage of instability, showing an excellentagreement between theory and simulation. Figure 2(c)displays the energy evolution of all components of theelectromagnetic fields. Only the B z component presentsexponential growth, as predicted earlier for solutions ofEq. (5a). Although the modes observed are oblique, thefact that the electric field does not grow [Fig. 2(c)] andthe consistency of the magnetic component with a pureelectromagnetic mode makes us confident that this is amanifestation of the Weibel, rather than the oblique in-stability [33, 34]. Additionally, the VDF [Eq. (6a)] isnot prone to the two-stream instability, hence there is nosource for the electrostatic modes, which are a compo-nent of the oblique instability. For completeness, Figure2(d) shows the saturated ( t = 6000 ω − p ) magnetic field B z in configuration space.The VDF shape plays a dominant role to determinefundamental aspects of the Weibel instability, such as un-stable wavevectors and their growth rates and which partof the VDF population is the most relevant for the in-stability. To demonstrate those points, we compare fourexamples, all in which the effective temperature in the x -direction is T x = R v x f x dv x = 1 . T ⊥ = 1 keV.Figure 3 shows theory and particle-in-cell simulation re-sults, with each row showing the initial VDFs f x and f y (left panel), the theoretical growth rate for a range of k (center panel), and B z in Fourier space at the linear stageof the instability taken from simulations (right panel).Since the effective temperature T x is the same in all ex-amples, theory predicts the same k = k y ˆ y modes. That isobserved in Fig. 3. The remaining unstable solutions aresubstantially different for each example. Figures 3(a-c)show results for a Maxwellian distribution function. Inaddition to the k = k y ˆ y expected modes, we notice thatthe unstable branch extends up to k x = 0 . ω p /c . Fig-ures 3(d-f) present results for a DLM ( m = 4) distribu-tion, where we observe a wide range of modes which areoblique and have higher growth rates than the k = k y ˆ y .In the case of Figs. 3(g-i), we explore a hot tail distribu-tion [Eq. (6b)], where the cold part has the same tem-perature as the transverse distribution T c = T ⊥ = 1 keV,and 1% of the particles are in the hot distribution with T h = 21 keV. We observe a narrow unstable branchwith k ≈ k y ˆ y , confirmed in simulations. Figures 3(j-l)explore a second example of hot tail distribution, with T c = 0 . δn = 0 .
01, and T h = 40 . k ≈ k y ˆ y , but this time thereis a larger branch with modes k ≈ k x ˆ x that extendsup to k y = 0 . ω p /c . The latter branch is commensu-rate with the solution when δn = 0, i.e., a bi-Maxwellianwith T x = 0 . T ⊥ = 1 . k = k x ˆ x ), once again demonstrating thatit does not necessarily lies along the cold direction fornon-Maxwellian VDFs.While the field free Vlasov equation has infinite solu- FIG. 3. Initial VDF (left column), theoretical growth rate(center column), and B z component of the magnetic field fromsimulations (right column) for four different examples. (a-c)Maxwellian with T x = 1 . m = 4 and T x = 1 . T c = 1 keV, T h = 21 keV, and δn = 0 .
01; and (j-l) hot tail with T c = 0 . T h = 40 . δn = 0 . T ⊥ = 1 keV in all the examples. tions of the kind f = f ( v ), anisotropic (but currentneutral) distribution functions evolve in ways that gen-erate magnetic fields in the plasma through the Weibelinstability. The temperature tensor T ij ∝ R d vv i v j f isalready diagonal in the v x v y coordinates for the exam-ples in Fig. 3. Thus the traditional temperature mea-surement cannot explain the presence of oblique modesand multiple branches. We conjecture that the instabil-ity will rise from any anisotropy in the VDF, and thegenerated magnetic field wavevector angle is perpendicu-lar to the direction of maximum velocity spread. To testthe conjecture, we define the quantity ε θ = Z ∞−∞ v f x ( v cos θ ) f y ( v sin θ ) dv Z ∞−∞ f x ( v cos θ ) f y ( v sin θ ) dv , (8)the VDF dispersion along the direction θ with respect to v x . This quantity mainly differs from temperature tensorin that integration is made along one direction and notall velocity space.Figure 4(a-d) shows ε θ as a function of θ for the ex-amples of Fig. 3 in the order that they appear. Wenotice that, for all the examples, ε = 1 . ε π/ = 1 . θ varies significantly for all the examples. The maximaof ǫ θ indicate the presence of unstable branches. Pan-els (a) and (c), for f x Maxwellian and hot tail with T c = T y , the maximum of ε θ for θ = 0 points to un-stable branches at θ = π/
2, as observed in Figs. 3(b)and 3(h). Additionally, the unstable branch size seemsto be related to the excess area above the minimum valueof ε θ = π/ = 1 . θ ≈ . π radimplies the unstable branch at θ ≈ . π rad, which is ingood agreement with the theoretical maximum growthrate wavevector for tan − ( k y /k x ) ≈ . π rad. An anal-ogous calculation for the example of Fig. 2 gives theunstable branch at θ ≈ . π rad, while the theoreticalmaximum growth rate is for tan − ( k y /k x ) ≈ . π rad.Finally, the presence of the different maxima at θ = 0and θ = π/ ε θ is related to the unstable branch size andgrowth rates, we define the quantityΓ = Z γ> γ ( k x , k y ) dk x dk y , (9)i.e., the area of the Fourier space with unstable solutions( γ >
0) weighted by the growth rate of each mode. Wecompare with the quantity∆ ε = Z π (cid:0) ε θ − ε min θ (cid:1) dθ, (10)where ε min θ is the minimum value of ε θ . Notice that ∆ ε =0 if f is isotropic, i.e., a stationary solution of Vlasov’sequation, and there is no instability.Figure 4(e) compares Γ and ∆ ε for several f x and pa-rameters (see Fig. 4 legend and caption for details) with f y being a Maxwellian with T ⊥ = 1 keV. The relationbetween Γ and ∆ ε is well described by a power lawΓ = 390 ω p c × (cid:18) ε min θ (cid:19) − . × (cid:18) ∆ εm e c (cid:19) . , (11)the solid line in Fig. 4. We verified Eq. (11) for other val-ues of ε min θ not shown in Fig. 4(e). While we do not havea proof that this is the most general mechanism for theWeibel instability, it assuredly leads to a better under-standing of the instability. Previous metrics rely on thetemperatures or the anisotropy parameter [35] and canonly reliably determine the maximum growth rate for bi-Maxwellian distribution functions. Even for these VDFs,there is a range of unstable oblique wavevectors that sig-nificantly contribute to the magnetic field generated, as FIG. 4. Panels (a-d): ε θ given in Eq. (8) for the same exam-ples in Fig. 3. Panel (e): relation between ∆ ε and Γ for sev-eral examples. Hot tail (a) : T c = 1 . δn = 0 .
01, with T h varying from 3 keV to 11 keV. Hot tail (b) : T c = 1 . δn = 0 .
05, with T h varying from 3 keV to 11 keV. Hottail (c) : T c = 1 . δn = 0 .
05, with T h varying from2 keV to 20 keV. DLM with m = 4 and T x = 1 keV to 10 keV.Maxwellian with T x = 1 keV to 10 keV. T ⊥ = 1 keV in all theexamples. The solid line represents the scaling given in Eq.(11). demonstrated in Fig. 3(b). In addition, the metric ε θ predicts instability even when f x = f y if f = f ( v ),which was confirmed in simulations with f x and f y beingDLMs ( m = 4) with T x = T y = 1 keV.In this Letter, we have studied the Weibel instabilityfor non-Maxwellian distribution functions and allowingoblique wavevectors. We have shown that the VDF shapeplays a significant role to determine the maximum growthrate and unstable modes. We have derived an empiricalformula for the maximum growth rate when the VDF issuper-Gaussian along one direction and the temperatureanisotropy is small. We have also shown that a bettermeasurement for the Weibel instability is based on pop-ulation spread excess along a particular direction ratherthan the temperature. Such quantity leads to a bettergrasp of the full-range Weibel unstable modes than byonly looking at the temperatures. We have explored ex-amples typical from laser-plasma interactions, thus show-ing that Weibel fields generated in those scenarios mayneed to consider the shape of the distribution function tocorrectly characterize the magnetic fields observed.TS is thankful for discussions with Dr. Kevin Schoef-fler. We gratefully acknowledge computing time onthe GCS Supercomputer SuperMUC (Germany); andPRACE for awarding us access to MareNostrum at BSC(Spain). This work was supported by the European Re-search Council through the InPairs project Grant Agree-ment No. 695088, FCT (Portugal) grant PTDC/FIS-PLA/2940/2014 and UID/FIS/50010/2023. The work ofBA was supported by a grant from the DOE FES-NNSAJoint program in HEDLP DE-SC 0018283. [1] A. Gruzinov and E. Waxman,The Astrophysical Journal , 852 (1999).[2] M. V. Medvedev and A. Loeb,The Astrophysical Journal , 697 (1999).[3] R. Schlickeiser and P. K. Shukla,The Astrophysical Journal , L57 (2003).[4] M. V. Medvedev, L. O. Silva, M. Fiore, R. A. Fonseca,and W. B. Mori, Journal of The Korean AstronomicalSociety , 533 (2004).[5] M. V. Medvedev, L. O. Silva, and M. Kamionkowski,The Astrophysical Journal , L1 (2006).[6] E. S. Weibel, Phys. Rev. Lett. , 83 (1959).[7] C. M. Huntington, F. Fiuza, J. S. Ross, A. B. Zyl-stra, R. P. Drake, D. H. Froula, G. Gregori, N. L. Kug-land, C. C. Kuranz, M. C. Levy, C. K. Li, J. Meinecke,T. Morita, R. Petrasso, C. Plechaty, B. A. Remington,D. D. Ryutov, Y. Sakawa, A. Spitkovsky, H. Takabe, andH.-S. Park, Nature Physics , 173 (2015).[8] E. N. Nerush, D. A. Serebryakov, and I. Y. Kostyukov,The Astrophysical Journal , 129 (2017).[9] C. Zhang, C.-K. Huang, K. A. Marsh,C. E. Clayton, W. B. Mori, and C. Joshi,Science Advances , eaax4545 (2019).[10] T. Silva, K. Schoeffler, J. Vieira,M. Hoshino, R. A. Fonseca, and L. O. Silva,Phys. Rev. Research , 023080 (2020).[11] G. Raj, O. Kononenko, M. F. Gilljohann, A. Doche,X. Davoine, C. Caizergues, Y.-Y. Chang, J. P. Coupe-rus Cabada˘g, A. Debus, H. Ding, M. F¨orster, J.-P.Goddet, T. Heinemann, T. Kluge, T. Kurz, R. Pausch,P. Rousseau, P. San Miguel Claveria, S. Sch¨obel,A. Siciak, K. Steiniger, A. Tafzi, S. Yu, B. Hid-ding, A. Martinez de la Ossa, A. Irman, S. Karsch,A. D¨opp, U. Schramm, L. Gremillet, and S. Corde,Phys. Rev. Research , 023123 (2020).[12] N. Shukla, K. Schoeffler, E. Boella, J. Vieira, R. Fonseca,and L. O. Silva, Phys. Rev. Research , 023129 (2020).[13] V. M. Vasyliunas, Journal of Geophysical Research , 2839 (1968).[14] H. Figueroa, C. Joshi, H. Azechi, N. A. Ebrahim, andK. Estabrook, The Physics of Fluids , 1887 (1984).[15] A. B. Langdon, Phys. Rev. Lett. , 575 (1980).[16] J. P. Matte, M. Lamoureux, C. Moller, R. Y.Yin, J. Delettrez, J. Virmont, and T. W. Johnston,Plasma Physics and Controlled Fusion , 1665 (1988).[17] A. Bendib, K. Bendib, and A. Sid,Phys. Rev. E , 7522 (1997).[18] K. Bendib, A. Bendib, and A. Sid,Laser and Particle Beams , 473–490 (1998).[19] A. Sid, Physics of Plasmas , 214 (2003).[20] S. Belghit, , A. Sid, and B. K.,International Journal of Physics , 226 (2014).[21] S. Belghit, Z. Beddiaf, O. Boultif, and A. Sid,Journal of Photonic Materials and Technology , 15 (2018).[22] M. Lazar, R. Schlickeiser, and S. Poedts,Physics of Plasmas , 062112 (2010). [23] B. B. Afeyan, A. E. Chou, J. P. Matte, R. P. J. Town,and W. J. Kruer, Phys. Rev. Lett. , 2322 (1998).[24] S. M. Shaaban, M. Lazar, S. Poedts, and A. Elhanbaly,Astrophysics and Space Science , 193 (2016).[25] S. M. Shaaban, M. Lazar, S. Poedts, and A. Elhanbaly,Astrophysics and Space Science , 13 (2016).[26] S. M. Shaaban, M. Lazar, R. A. L´opez, and S. Poedts,The Astrophysical Journal , 20 (2020).[27] N. Krall and A. Trivelpiece, Principles of Plasma Physics , International Seriesin Pure and Applied Physics No. v. 0-911351 (McGraw-Hill, 1973).[28] See Supplemental Material at [URL] for a detailed deriva-tion of Eq. (3).[29] M. Tzoufras, C. Ren, F. S. Tsung, J. W. Tonge,W. B. Mori, M. Fiore, R. A. Fonseca, and L. O. Silva,Physics of Plasmas , 062108 (2007).[30] K. Estabrook, W. L. Kruer, and B. F. Lasinski,Phys. Rev. Lett. , 1399 (1980).[31] R. A. Fonseca et al. , OSIRIS: A Three-Dimensional, FullyRelativistic Particle in Cell Code for Modeling PlasmaBased Accelerators, in Computational Science — ICCS2002 , edited by P. M. A. Sloot, A. G. Hoekstra, C. J. K.Tan, and J. J. Dongarra (Springer Berlin Heidelberg,Berlin, Heidelberg, 2002) p. 342.[32] Simulations were performed with periodic boundaries.The box size is 400 c/ω pe × c/ω pe divided into agrid with square cells with edges 0 . c/ω pe long. Eachcell starts with 15625 particles.[33] A. Bret, The Astrophysical Journal , 990 (2009).[34] A. Bret, L. Gremillet, and M. E. Dieckmann,Physics of Plasmas , 120501 (2010).[35] K. M. Schoeffler and L. O. Silva,Plasma Physics and Controlled Fusion60