Ultrafast Electron Holes in Plasma Phase Space Dynamics
UUltrafast Electron Holes in Plasma Phase Space Dynamics
S. M. Hosseini Jenab , I. Kourakis , G. Brodin , J. Juno Department of Physics, Chalmers University of Technology , G¨oteborg, Sweden Department of Mathematics, Khalifa University of Science and Technology, Abu Dhabi, UAE Department of Physics, Ume˚a University , Ume˚a, Sweden Institute for Research in Electronics and Applied Physics, University of Maryland, College Park, MD 20742, USA (Dated: September 29, 2020)Electron holes (EH) are localized modes in plasma kinetic theory which appear as vortices inphase space. Earlier research on EH is based on the Schamel distribution function (df). A noveldistribution function is proposed here, generalizing the original Schamel df in a recursive manner.Nonlinear solutions obtained by kinetic simulations are presented, with velocities twice the electronthermal speed. Using 1D-1V kinetic simulations, their propagation characteristics are traced andtheir stability is established by studying their long-time evolution and their behavior through mutualcollisions.
Plasma phase-space dynamics is tacitly characterizedby the occurrence of electron holes, a term describinga localized plasma region where electrons are trappedby the electric potential stemming from their own self-generated density variation, as a localized electron de-pletion region occurs in a self-consistent manner. Anelectron hole is thus manifested as a localized “trapped”electron population traveling alongside an electrostaticpotential disturbance . Electron-holes present two maincharacteristics : a localized positive potential structurewhich traps electrons, and a symmetry in the electricpotential profile around the peak. In addition, electronholes are a type of Bernstein, Greene, and Kruskal (BGK)mode . Electron holes have been observed and studied inlaboratory experiments , in space measurements andin kinetic simulations .In order to construct electron holes in a self-consistentmanner within a kinetic model, one may either start withan arbitrary potential profile and then proceed by de-riving the distribution function (df) of an electron hole,or, inversely, start with a predefined df for the trappedelectrons and thus derive the associated potential pro-file. The former (integral equation) method, due toBernstein, Greene and Kruskal leads to an infinity ofsolutions whose dynamical stability is not prescribed.The latter (differential equation) method, suggested bySchamel , is based on a parametrized df (henceforthreferred to as “the Schamel df”) allowing one to prescribethe shape of the trapped population (i.e. by assigning avalue to parameter β associated with the inverse temper-ature of the trapped population). The Schamel methodcombined with the Sagdeev pseudopotential approach may provide initial conditions for a controlled numericalinvestigation of EH dynamics . Recent studies haveshown that the Schamel-Sagdeev approach can producenonlinear solutions with Mach numbers 1 . < M < . . < M < . and to survive mutualcollisions . In other words, structures are destabilized asthe Mach number increases. This has been suggested inother kinetic simulations . For very high Mach number ( M > β (values) .Despite these theoretical challenges, the existence ofhigh-speed electron holes is a topic of intense study,first getting attention due to observations by the FASTsatellite . Saeki et al studied electron holes ex-perimentally using a Q-plasma machine and also via ki-netic simulations; they reported structures moving at theelectron thermal speed, which they identified as solitons.Solitons are nonlinear structures that can survive mutualcollisions and are characterized by a phase shift duringa collision . We note however, Saeki et al did notconsider the phase shift separating the hole trajectoriesbefore and after collisions.The aim of this study is to characterize high-speed elec-tron holes by establishing their occurrence in a kineticframework, and by investigating their stability profileand probing their soliton-like features. For this purpose,a novel distribution function ( df ), the ‘ ELIN df’, is intro-duced as a generalization of the Schamel df. The
ELIN df adjusts the distribution function of the trapped pop-ulation of electrons by relying on a dynamically varyingparameter β so that its moments can fit a predeterminedcurve and all of the desired featured of the Schamel dfare retained, such as consistency and smoothness in bothspatial and velocity spaces inside the trapped region.To show the stability of our nonlinear solutions, threeseries of simulations are reported. Firstly, by consider-ing the long-time evolution of an initial condition we willconfirm the stability of the solution’s profile during prop-agation, thus establishing them as solitary waves. Then,two types of mutual collisions are reported, i.e. head-on collisions (with no overlapping in velocity space) andovertaking collisions (moving in parallel and with over-lapping). The aforementioned phase shift through colli-sions has also been investigated, to corroborate the factthat electron holes behave as solitons.The scaled Vlasov-Amp`ere system of equations form- a r X i v : . [ phy s i c s . p l a s m - ph ] S e p FIG. 1: Different distribution functions for the trapped elec-tron population are presented. The Maxwellian df (in theabsence of trapped particles) is shown for sake of compari-son (blue, thin dotted line). Three shapes of the Schameldf are displayed, namely flat (brown dashed, β = 0), hollow(red, dashed-dotted, β = −
2) and a bump (green, dashed-dotted, β = 2), for φ = 25. The ELIN df (black thick line)is shown when ten carving ( φ = 2 . φ = 5, φ = 7 . φ = 25) is carried out with their corresponding β ( β = − , β = − . , β = − . , ..., β = 0). ing the basis of our simulation reads: ∂f s ( x, v, t ) ∂t + v ∂f s ( x, v, t ) ∂x + Υ s E ( x, t ) ∂f s ( x, v, t ) ∂v = 0 , (1) ∂E ( x, t ) ∂t = (cid:88) q s J s ( x, t ) (2) where s = i, e represents the corresponding species, i.e.ions and electrons respectively. The factor Υ s takesthe values Υ e = − i = 1. The normalizedcharges are q e = − q i = 1. The above equa-tions are coupled by integrations for each species, viz. J s ( x, t ) = (cid:82) f s ( x, v, t ) vdv in order to form a closed setof equations for J , denoting the current (contribution)generated by by species s . To derive the above (dimen-sionless) equations, all physical quantities were normal-ized to suitable scales related with ionic parameters, i.e.mass ( m s ) was divided by the ion mass ( m i ), temper-ature ( T s ) by ion temperature ( T i ), charge ( q s ) by theelementary charge ( e ), time ( τ ) by the ion plasma period( ω / pi = (cid:0) n i e m i (cid:15) (cid:1) − ), and length ( L ) by the ion Debyelength ( λ Di = (cid:113) (cid:15) K B T i n i e ). Here, K B is Boltzmann’s con-stant and (cid:15) is the permittivity of free space.We have employed the Gkeyll simulation framework to solve the Vlasov-Ampere system of equations . Gkeyll discretizes the equations using the discontinu-ous Galerkin finite element method in space, with astrong stability-preserving Runge-Kutta method in time.We have adopted a piecewise cubic Serendipity Elementspace for the basis expansion (further details can befound in Refs. Juno et al. and Hakim and Juno ).The Gkeyll method has been compared to standard PICmethods and showed to yield better results .In our study, the temperature and mass ratio are T e T i = 100 and m i m e = 1836. The initial distribution func-tion f is considered to be the Maxwellian df (= D m ). The size (length) of the simulation box is l = 1000in the x -direction. In the v direction for each species,we have different limits: for the electrons we have v =( − , v th e = ( − , v = ( − , v th e = (cid:113) T e T i m i m e ≈ . n X = 2000, n V = 1000 for both electronsand ions. The time step d τ ≈ − is chosen in order tofulfill Courant-Friedrichs-Lewy (CFL) condition .The electron hole speed ( v EH ) is expressed by the“Mach number”, which is defined as the ratio M = v EH c s ,where c s = (cid:113) γ e T e + γ i T i m i is the ion sound speed. As-suming γ e = γ i = 3 , T e = 100 T i and m i = 1, the ionsound speed in our simulations is c s = √ ≈ . (cid:0) φ ( x ) (cid:1) and by choosing the value of the electronhole speed ( v EH ). We then use the ELIN df to producethe electron distribution function. Given that the poten-tial profile provides the charge density ( ρ ( x )), and usingthe Schamel df for the ions to obtain n i ( x ), we then usethe total charge density (profile) n e ( x ) = ρ ( x ) − n i ( x )as a guiding equation for the ELIN df and thus constructthe electron hole. We have adopted, to start with, thesimplest form of potential profile suggested for electronholes i.e. φ = A sech p ( x/L ) in which p = 2 and A and L are the EH amplitude and length, respectively. The am-plitude and length (values) are chosen randomly; how-ever the system will damp/break the forced profile if it isnot close-enough to a self- consistent nonlinear solution.The resulting electron hole may have different size andvelocity, but with an iterative process, one can find thecombination of { A, L } for which the solution will be sta-ble enough for a specific (chosen) velocity value. Since weare not aware of the nonlinear dispersion relation, i.e. arelationship between { A, L, M, p } for the exact nonlinearsolution(s), a sequence of trials is performed to iterate tothe correct combination of { A, L, p } for a given M . Inthe simulations presented here three electron holes werestudied, e.g. EH M = 45 , A = 19 , L = 22 . , P = 2 EH M = − , A = 9 . , L = 22 . , P = 2 EH M = 30 , A = 19 , L = 22 . , P = 2The Schamel distribution function can be written interms of the energy as: f = af ( ε K ) in which a is nor-malization constant and f ( ε K ) = D g ( ε K − ε φ ) if ε K > ε S + ε φ D g ( ε S ) if ε K = ε S + ε φ D g ( ε S ) D m (cid:16) β (cid:0) ( ε K − ε φ ) − ε S (cid:1)(cid:17) if ε K < ε S + ε φ (3) where we have defined the kinetic energy of particles ε K = mT v and the Maxwellian distribution function D m ( ε K ) = exp( − ε K ). Here, the potential energy associ-ated with the potential profile at a given point in the x di-rection is ε φ = qφ ( x ), and the kinetic energy of particlestraveling at the electron hole velocity is ε S = mT v EH .The first and the third lines respectively represent the FIG. 2: The electron phase space is shown for for case EH τ = 12.FIG. 3: The electrostatic potential/E-field profile of EH τ = 0 (red dotted curve) is compared with τ = 12 . free and the trapped population(s). The second line isthe separatrix between the free and trapped population.By D g we denote a general distribution function satis-fying the Vlasov equation, i.e. in principle any functiondepending on the constant(s) of motion. Here, the energyis used to construct a valid function, hence D g = f ( ε K ).Well-known examples of D g are the Maxwell-Boltzmanndf, the κ df and the Cairns distribution func-tion(s).The df of trapped population can be broken down totwo parts: f trapped = f base × f shape in which f base = D g ( ε S ) and f shape = D m ( β ( ε K − ε φ − ε S )). The basedistribution function f base is just the distribution func-tion at the separatrix where ε K = ε S + ε φ (and hence D g ( ε K − ε φ ) = D g ( ε S )). The second component, f shape (controlled by β ), provides the shape of the trapped dfaround the f base within trapped region. It may appearin three qualitative shapes, i.e. flat , a bump or a hollow FIG. 4: The electrostatic potential profile of EH EH τ = 1 . τ = 1 .
6) and c) af-ter ( τ = 1 .
9) the collision. Dotted lines represent the initialcondition for each of the solitary wave as if they are propagat-ing without any numerical noise or collisions. Red/blue is forEH1/EH2 which is propagating to the right/left. After thecollision, the overall shape and velocity of the solitary waveremains intact. curve, if β = 0, β > β <
0, respectively (see Fig. 1).One can understand the Schamel df as carving up agiven general distribution function ( D g ) around a partic-ular velocity (hole velocity) and inserting a Maxwelliandf with arbitrary temperature inside the hole ( f trapped ).Now, in order to generate the new distribution functionfor another φ ( x ), one only needs to update the value of f base in our approach and repeat the process. Using thismethod of multiple carvings, a specific φ can be createdby n recursive iterations, φ = φ + φ + φ + ... + φ n = (cid:80) n φ i . In other words, one can use the previous distribu-tion function f i − as the D g for the current distributionfunction ( f i = af ( φ i )) and thus arrive at the ELIN dis-tribution function: f ( φ i ) = f i − ( ε K − ε φ i ) if ( ε K − ε S ) /q > φ i f i − ( ε S ) if ( ε K − ε S ) /q = φ i f i − ( ε S ) D m ( β i ( ε K − ε φ i − ε S )) if ( ε K − ε S ) /q < φ i (4) in which f is the initial unperturbed df (here assumingMaxwellian df, i.e. f = D m ). β i can change arbitrarilyin order for moments of df to fit a “guiding equation”(here, the equation for the electron density). To obtain asmooth distribution function in the x direction, one can FIG. 5: An overtaking collision between EH EH EH τ = 2 . τ = 3 .
20) and c) after ( τ = 4 . rmsech before (blue) and after (red) the collision, for EH
1. A shiftin the position of the first EH can be witnessed (note the difference between the red and the blue curves) manifesting a phaseshift, as intuitively expected.FIG. 6: The electron phase space is presented for an overtak-ing collision between EH EH EH
1. There is a substantial overlap in velocity direction(a). During collision, the interaction is strong (b). After thecollision, EH increase n until the numerically-desired level of smooth-ness is achieved: φ i = δφ . An example of the ELIN dfprofile is presented at Fig. 1 which shows 10 successive(carving) iterations with β approaching zero from below(negative side).Figs. 2 and 3 display the temporal evolution of EH
1. The initial condition and the last step of temporal evolu-tion can be compared and show that the overall shape ofthe electron hole (fig.2) and the corresponding potentialor field profile (fig.3)stay unperturbed.Fig. 4 depicts a head-on collision between EH EH
2. After the collision (0 < τ < EH EH
3. Fig. 5 presents the temporal evolution of elec-tric field/potential around the collision time τ = 3 . M = 45. Both EHs survive the col-lision, and their respective velocity stays the same. Fo-cusing on EH
1, displacement can be witnessed after thecollision. A phase shift can be measured by comparingEH profile with the red line, which is an extrapolationof an unperturbed path of this EH. This displacement issimilar to the well-known effect of “phase shift” whichobserved to happen in mutual collisions of solitons .We show in Fig.6 the electron df during the overtakingcollision, which demonstrates the considerable interac-tion between the EHs during the collision and their over-lapping on velocity direction. Yet after the collision the EH curve form approximates the numerical data better thanany other exponent, including the (expected, arguably)sech form (see Eq. 39 in Schamel ).In summary, we have provided a method to producehigh-speed nonlinear solutions which move at a speed be-yond the electron thermal speed. We showed that theseelectron holes are stable, retain their profile through col-lisions and remain so in the entire duration of the sim-ulation. For mutual collisions with considerable overlapin the velocity direction, the EHs display a “phase shift”This phase shift represents a signature of soliton behav-ior and hence suggests that these EHs can be consideredas solitons (at least approximately). This has been sug-gested for much lower-speed EHs before but without theobserved “phase shift” reported here. Acknowledgments
This project has received funding from the EuropeanResearch Council (ERC) under the European Union’s Horizon 2020 research and innovation programme undergrant agreement No 647121. Simulations were performedon resources of High Performance Computing CenterNorth (HPC2N) provided by the Swedish National In-frastructure for Computing (SNIC). J. Juno was sup-ported by a NASA Earth and Space Science Fellowship(Grant No. 80NSSC17K0428). One of us (IK) grate-fully acknowledges support from ADEK (Abu Dhabi De-partment of Education and Knowledge) in the form ofan AARE-2018 (ADEK Award for Research Excellence2018) grant. H. Schamel, Phys. Reports , 161 (1986). I. H. Hutchinson, Physics of Plasmas , 055601 (2017). I. B. Bernstein, J. M. Greene, and M. D. Kruskal, PhysicalReview , 546 (1957). K. Saeki, P. Michelsen, H. P´ecseli, and J. J. Rasmussen,Physical Review Letters , 501 (1979). R. Ergun, C. Carlson, J. McFadden, F. Mozer, G. Delory,W. Peria, C. Chaston, M. Temerin, I. Roth, L. Muschietti,et al., Geophysical Research Letters , 2041 (1998). J. R. Franz, P. M. Kintner, and J. S. Pickett, Geophysicalresearch letters , 1277 (1998). H. Matsumoto, H. Kojima, T. Miyatake, Y. Omura,M. Okada, I. Nagano, and M. Tsutsui, Geophysical Re-search Letters , 2915 (1994). H. Kojima, H. Matsumoto, S. Chikuba, S. Horiyama,M. Ashour-Abdalla, and R. Anderson, Journal of Geophys-ical Research: Space Physics , 14439 (1997). X. Deng, R. Tang, H. Matsumoto, J. Pickett, A. Fazaker-ley, H. Kojima, W. Baumjohann, A. Coates, R. Nakamura,D. Gurnett, et al., Advances in Space Research , 1373(2006). B. Eliasson and P. K. Shukla, Phys. Reports , 225(2006). H. Schamel, Plasma Physics , 491 (1971). H. Schamel, Journal of Plasma Physics , 1 (1972). H. Schamel, Plasma Physics , 905 (1972). H. Schamel, Physica Scripta , 306 (1979). R. Sagdeev, Reviews of Plasma Physics , 23 (1966). S. Hosseini Jenab, F. Spanier, and G. Brodin, Physics ofPlasmas , 072304 (2018). S. Hosseini Jenab and G. Brodin, Physics of Plasmas ,022303 (2019). V. Turikov, Physica Scripta , 73 (1984). L. Muschietti, I. Roth, R. Ergun, and C. Carlson, Physicsof Fluids nonlinear processes in geophysics (1999). L. Muschietti, R. Ergun, I. Roth, and C. Carlson, Geo-physical research letters , 1093 (1999). Y. Nishida and T. Nagasawa, Physical Review Letters , 1626 (1980). F. Verheest, M. A. Hellberg, and W. A. Hereman, PhysicalReview E , 036402 (2012). Y. Nakamura, H. Bailung, and K. Lonngren, Physics ofPlasmas , 3466 (1999). T. R. Marchant, Physical Review E , 046623 (2002). H. Demiray, J. Computational & Applied Math , 826(2007). A. H. Hakim, N. R. Mandell, T. Bernard, M. Francisquez,G. Hammett, and E. Shi, Physics of Plasmas , 042304(2020). J. Juno, A. Hakim, J. TenBarge, E. Shi, and W. Dorland,Journal of Computational Physics , 110 (2018). A. Hakim, M. Francisquez, J. Juno, and G. W. Hammett,Journal of Plasma Physics , 905860403 (2020). A. Hakim and J. Juno, arXiv preprint arXiv:2004.09019(2020). D. N. Arnold and G. Awanou, Foundations of Computa-tional Mathematics , 337 (2011). J. Juno, M. M. Swisdak, J. M. Tenbarge, V. Skoutnev,and A. Hakim, Journal of Plasma Physics , 175860401(2020). R. Courant, K. Friedrichs, and H. Lewy, Mathematischeannalen , 32 (1928). R. Courant, K. Friedrichs, and H. Lewy, IBM journal ofResearch and Development , 215 (1967). V. M. Vasyliunas, Journal of Geophysical Research ,2839 (1968). V. Pierrard and M. Lazar, Solar Physics , 153 (2010). D. Summers and R. M. Thorne, Physics of Fluids B:Plasma Physics , 1835 (1991). M. Hellberg, R. Mace, T. Baluku, I. Kourakis, andN. Saini, Physics of Plasmas , 094701 (2009). R. Cairns, A. Mamum, R. Bingham, R. Bostr¨om,R. Dendy, C. Nairn, and P. Shukla, Geophysical ResearchLetters22