Electron heating by intense short-pulse lasers propagating through near-critical plasmas

We investigate the electron heating induced by a relativistic-intensity laser pulse propagating through a near-critical plasma. Using particle-in-cell simulations, we show that a specific interaction regime sets in when, due to the energy depletion caused by the plasma wakefield, the laser front profile has steepened to the point of having a length scale close to the laser wavelength. Wave breaking and phase mixing have then occurred, giving rise to a relativistically hot electron population following the laser pulse. This hot electron flow is dense enough to neutralize the cold bulk electrons during their backward acceleration by the wakefield. This neutralization mechanism delays, but does not prevent the breaking of the wakefield: the resulting phase mixing converts the large kinetic energy of the backward-flowing electrons into thermal energy greatly exceeding the conventional ponderomotive scaling at laser intensities > 10 21 W cm − 2 and gas densities around 10% of the critical density. We develop a semi-numerical model, based on the Akhiezer–Polovin equations, which correctly reproduces the particle-in-cell-predicted electron thermal energies over a broad parameter range. Given this good agreement, we propose a criterion for full laser absorption that includes field-induced ionization. Finally, we show that our predictions still hold in a two-dimensional geometry using a realistic gas profile.

Although the mean energy ( h ) of the hot electron distribution is a complex function of the laser-plasma parameters, it is commonly assumed to obey the ponderomotive scaling,  , where w º a eE m c e 0 0 0 is the normalized laser amplitude [59,60]. Fewer studies have addressed the case of intense short-pulse lasers interacting with near-critical (ñ n e c ) plasmas [61][62][63]. This relative lack of interest is probably due to the experimental difficulty of achieving, in a controlled way, the required high gas densities (~-n 10 cm e 20 21 3 ) [64] and, also, to the fact that such systems are a priori ill-suited to efficient wakefield electron acceleration. In recent years, however, high-density gases have attracted increasing attention as promising high-energy ion sources, based on a variety of mechanisms. First, as in solid foils, the so-called target normal sheath acceleration (TNSA) [65,66] can arise from the spacecharge field set up at the plasma boundaries by the laser-driven hot electrons [2]. Second, in a dense enough plasma, the front-side ions are pushed by the 'laser piston', i.e., the space-charge field resulting from the laser ponderomotive force on the electrons. Third, under specific conditions, the laser piston or the electron pressure gradients created in an inhomogeneous plasma can launch a collisionless electrostatic shock [67][68][69][70], accelerating a fraction of the background ions to energies possibly larger than through TNSA [60,[71][72][73][74]. Further, for tightly focused lasers, ion acceleration can proceed via the electric field induced by magnetic vortices moving down density gradients [75,76].
The objective of this work is to illuminate the processes of electron energization during the interaction of high-intensity, ultrashort (100 fs) laser pulses with near-critical plasmas. Developing such an understanding is a prerequisite for modeling the above ion acceleration processes that closely depend on electron pressure gradients and the level of electron heating. Such is the case of TNSA, collisionless shock acceleration (CSA) and magnetic vortex acceleration, all necessitating high mean electron energies so as to optimize the energy of the accelerated ions. Regarding particularly CSA, a relatively large plasma density is also needed to ensure a short enough shock formation time, of the order of the inverse ion plasma frequency, wpi 1 [71]. Such conditions imply strong absorption of the laser through the plasma, and therefore entail a laser depletion length shorter than caused by wakefield excitation in dilute plasmas [19,20,77]. For the above reasons, an in-depth insight into the electron heating mechanisms at play in near-critical plasmas is desirable. In the following, through extensive 1D particle-in-cell (PIC) simulations performed with the CALDER code [78] and a reduced model, we clarify the dependence of the final mean electron energy with the laser intensity and plasma density. Most notably, we find that very high electron energization levels (i.e. well above the ponderomotive scaling) can be achieved provided (i) the laser has propagated a distance large enough that it has developed a steep rising front; (ii) breaking of the wakefield and phase mixing between the trapped electrons and return electrons have occurred. Moreover, as this interaction regime involves efficient propagation of the laser pulse through an initially neutral gas jet, we derive a simple criterion for laser propagation that includes field ionization. Finally, we show that our 1D simulation and model results are also relevant to laser pulses focused to moderate spot sizes into unionized gases. This paper is organized as follows. In section 2.1 is presented a tutorial 1D PIC simulation serving to illustrate the dominant interaction processes leading to efficient electron heating in near-critical plasmas. In section 2.2, using a set of reduced fluid equations and several assumptions supported by PIC simulations, we estimate the final mean electron energy as a function of the laser-plasma parameters. In section 2.3, the predictions of this model are compared with a 1D-simulation scan over a broad parameter space. In section 3, we provide a criterion for full laser absorption, taking account of field-induced ionization. In section 4, 2D PIC simulations are carried out to investigate the effects of a finite laser width. For completeness, we consider the case of a laser interacting with a realistic high-density gas jet produced at the Laboratoire d'Optique Appliquée (LOA). Finally, our main results are summarized in section 5.

Collisionless electron heating
While multidimensional kinetic simulations are essential to reproduce experimental conditions, a 1D description is most convenient to highlight the main physical processes and to extract some tendencies and scalings. In this section, a typical 1D PIC simulation of a laser pulse interacting with a near-critical gas is presented to reveal the basic electron heating mechanisms. Ion acceleration-via TNSA or an electrostatic shock wave-is also shown to be a promising application of such interaction regimes. We then propose a simple model that reproduces the salient features of electron energization. Finally, extensive 1D simulations are performed to support the model results and, in particular, the scaling of the electron energy with laser intensity and plasma density.  ). This cold-fluid structure is the signature of the nonlinear electrostatic wakefield induced by the laser pulse. The associated normalized electrostatic potential ( f -e m c e 2 ) is overlaid as a white solid curve in figure 1(b). In the absence of wave breaking, and assuming stationary propagation in the laser's coordinate system ( -

Basic mechanisms
where v g is the laser group velocity), this wakefield would be described by the well-known Akhiezer-Polovin (AP) equations [80,81], and hence would make the plasma electrons oscillate around p x = 0 [16,19,20]. In the present case, however, wave breaking and phase mixing have occurred, giving rise to a broad distribution (   p m c 70 150 x e ) of energetic electrons behind the laser head. This hot electron population is dense enough to neutralize the cold-fluid electron flow during its backward acceleration by the electrostatic field. This process, referred to as beam loading in wakefield accelerators [82], is well illustrated by the electrostatic potential, which, after abruptly dropping behind the laser front ). Further behind the laser beam, the nonlinear wakefield initially sustained by the cold electron flow, finds itself in the presence of the counterstreaming cold and hot electron beams. The hot electron flow initially acts to neutralize the wakefield, hence preventing its rapid breaking (i.e., after just one oscillation). The occurrence of wavebreaking, however, is only delayed for a few plasma periods, probably as a consequence of the generic unstable behavior of relativistic plasma waves [81,83]. The perturbations responsible for this wavebreaking result from the hot electron stream but, also, from the non-stationary character of the wave itself (which, in contrast to a pure AP solution, does not propagate at a uniform phase velocity). This phenomenon is further discussed in section 2.2 B.450. Through phase mixing, the drift energy of the cold-fluid electrons ends up being converting into heat, with a mean energy ( g á -ñ 1 , blue solid line in figure 1(b)) in the ultra-relativistic range, g á -ñ1 40. As a consequence of the electron heating, the ions are set into motion through TNSA in the density downramp, as illustrated in figure 2 where the ion - x p x phase space is plotted at successive times. At times w = t 1450 0 and w = t 1450 0 , the ion phase space takes the form of a self-similar profile (i.e., roughly scaling as x/t) interrupted by a velocity spike, increasing in time as i n , (L n is the initial electron density scale length). This structure is reminiscent of the wave steepening occurring in an expanding nonuniform plasma [84]. Further away, a secondary self-similar expansion takes place due to charge separation (here the local Debye length exceeds the ion density scale length) and the rarefaction wave coming from the front of the initially truncated ion profile (at w = x c 1400 0 ). At wt 5000 0 , the velocity spike has evolved into a collisionless shock [68,69]. The associated moving potential barrier reflects part of the TNSA-expanding upstream ions, thereby boosting their non-relativistic velocity by a factor of ∼2 [60]. Since the shock velocity roughly scales as the sound speed of the upstream plasma, the final ion velocity is proportional to Note, however, that in a multidimensional geometry and for a small laser spot size, the formation of a magnetic dipole vortex in the plasma down-ramp [75,85] may change the scaling of the accelerated ions, as observed in previous simulation studies [76,86].

1D model for the bulk electron dynamics
A self-consistent analytic modeling of the collisionless electron heating, which comprehends laser propagation effects, is a challenging task, well beyond the objectives of this study. However, some insight can be gained from a simple model supported by numerical simulations.
As noted above, the bulk electron heating proceeds through a sequence of wavebreaking, two-stream and phase mixing processes. These phenomena are illustrated by a 1D simulation similar to that discussed in the previous section, except that the plasma profile is now chosen to be uniform ( = = n n n 0.032 i e c ). Figure 3 presents the electron - x p x phase space at different times. Early in time (w = t 937.5 0 ), the nonlinear wakefield induced by the laser pulse is already strong enough to have trapped a few electron bunches, and triggered phase mixing near the vacuum/plasma interface (w~x c 300 500 0 ). Despite these growing kinetic processes, the first nonlinear oscillations of the bulk electron flow roughly obey the cold-fluid equations of [16,19,20]. This remains true at w = t 1250 0 , even though phase mixing is now well developed after a few wakefield oscillations. Since the laser pulse is longer than the local plasma wavelength (t w w = > c c c 70 30 p 0 0 ), the wakefield is generated in the laser beam head [19]. The sustained energy transfer from the laser to the wakefield depletes the laser beam head, so that the laser front profile becomes increasingly steep during its propagation [19,23,87]. This self-steepening effect is clearly seen in the insets of figures 3(a)-(c), where the normalized transverse electric field ( w = a eE m c y y e 0 ) is plotted at the same times as the electron phase spaces. The initially Gaussian laser beam progressively develops a short gradient at its front, down to roughly one laser wavelength (inset of figure 3(c)). This 'optical shock' increases the longitudinal ponderomotive force µ ¶ . In this regime, the energy loss rate of the laser beam largely exceeds the depletion rate derived in the 'classical' wakefield regime [19,20,77]. The transition from the classical depletion regime to the optical-shock/phase-mixing-assisted regime is clearly demonstrated by the time evolution of the laser and kinetic particle energies shown in figure 4. At early times (  w t 1800 0 ), the energy loss rate of the laser pulse (and, conversely, the electron energy gain) is small. Once the beam loading and the associated phase mixing become significant (  w t 1800 0 ), the energy loss rate of the laser pulse drastically increases.
Both the mean drift momentum/energy of the cold electron flow and the mean thermal energy of the hot electrons are correlated with the electrostatic potential gap set up behind the laser beam head. The latter can be estimated by a simple extension of the well-known quasi-stationary wakefield model [19]. In our model, we assume that (i) the laser pulse propagates at the speed of light; (ii) the laser front has already been depleted, so that its sharp intensity gradient extends over , the laser field profile is taken to be of the following form . Here, the Lorentz factor writes g = + + p a 1 x y 2 2 due to conservation of the transverse canonical momentum (valid in a 1D geometry). The Heaviside term h µ -+ + ( ) H x t in the right-hand side of equation (4) represents the hot electron flow. The latter is supposed to be a rigid beam sited behind the laser front, at the position, η, where the electric field is neutralized: , 0 x . The (uniform) hot electron beam density is chosen to be equal to half the initial plasma density (n 2 0 ) owing to conservation of the total electron density, as justified below. Looking for quasi-stationary solutions depending on x =x t only, equations (2)-(4) can be recast as .    x p x phase spaces extracted from the corresponding 1D PIC simulations. The simulated electrostatic potential is also plotted as a white solid line. Qualitative agreement is found between the early shapes of the theoretical and PIC curves, i.e., prior to the occurrence of phase mixing effects. The differences between the model predictions and the PIC curves can be ascribed to the non-stationary character of the cold electron response.
The oscillations in the cold electron flow, discernible in the phase spaces of figures 3(c), (d) and figures 5(a), (b), can be attributed to the non-stationarity of the laser field in the coordinate system x =x t. More precisely, owing to the strong discontinuity of the electron bulk density-increasing from n 0 tõ - ( ) n v n 1 x 0 0 over a laser wavelength-at the laser front, the latter propagates at a velocity v L significantly different from the laser phase velocity v f . This difference induces oscillations in the laser ponderomotive force, which causes the bulk electrons to leave the interaction region with an initially oscillating negative momentum, of period p - The wakefield can thus be viewed as an AP wave perturbed by large oscillations in velocity and density, further enhanced by the collective interaction with the counterstreaming hot electrons. Such a modulated AP wave is known to be prone to fast wave breaking, at amplitudes well below the limiting amplitude of a pure AP wave [81,83]. Under the present conditions, we find that phase mixing occurs after several tens of relativistic plasma periods g n f 0 . We point out that the sole two-stream instability between the cold and hot electron flows (both supposed of uniform density and velocity profiles) cannot account for the fast phase mixing observed in the simulations. The maximum two-stream growth rate (G max ) can be estimated assuming a simple waterbag distribution function, , a value well above the observed phase-mixing time (~10 3 ). A tentative description of the observed phase mixing process in terms of the two-stream instability should therefore include the nonlinear wakefield that modulates both electron flows. A similar calculation, yet restricted to the nonrelativistic limit, has demonstrated that G max strongly increases with the amplitude of the preexisting electrostatic wave [90,91]. The extension of this problem to the relativistic regime is left for a future work.
Neglecting the ion plasma expansion, the conversion from electron kinetic to thermal energy can be assumed to be isochoric. Therefore, from energy conservation, the final bulk electron energy should approximately verify , where g á ñ h is the mean energy of the hot electrons As the electron energies are significantly higher than the potential variations df in the phase-mixing region, and the bulk electron energy is expected to converge to the hot electron energy well behind the laser pulse ( g g á ñ á ñ  h ), one obtains for the final bulk electron energy: The mean electron energy predicted by equation (9) is plotted in figure 6(a) over the parameter space   a 5 100 0 and   n 0. 1 2 0 . Some numerical studies [60,92] suggest that the mean electron energy follows a ponderomotive-like scaling (i.e., g á -ñ µ a 1 0 ) in the case of a laser pulse interacting with a marginally overcritical plasma. To compare our results with the ponderomotive scaling, we display in figure 6

1D PIC simulations: electron energy scaling
To test the accuracy of the above scaling for the electron heating, we performed an extensive PIC simulation scan over the parameter space of figures 6(a), (b). In all (3600) 1D simulations, the ions were kept fixed to avoid electron cooling from plasma expansion. The laser temporal profile was taken to be a l 10 0 -long plateau surrounded by two l 1 0 -long ramps so as to mimic optical shock formation. The plasma length was adjusted as a function of the (uniform) plasma density in order to allow the laser pulse to propagate through the whole gas length. Each cell, of fixed size D = x 0.075, was filled by 100 electrons and 100 protons in the plasma region. Absorbing boundary conditions were applied to both particles and fields. To reduce the statistical fluctuation level, the mean electron energy, g á -ñ 1 , was computed over the whole simulated gas length, at the time when the laser pulse reached the end of the plasma. Figures 6(c), (d) display the resulting values of g á -ñ 1 and g á -ñ + -( ) a 1 1 2 1 0 2 in the same parameter space ( ) a n , 0 0 as in figures 6(a), (b). The PIC results are observed to corroborate the model predictions. The quantitative agreement between figures 6(a)-(d) as well as the scalings that these maps evidence constitute the main results of this paper.
The electron heating process under consideration relies on enhanced phase mixing in the plasma wakefield as a result of optical shock formation and beam loading. The two latter phenomena should then occur fast enough that the electron heating takes place over most of the gas length and the laser pulse is efficiently absorbed. The minimum gas length (or areal density) for complete laser absorption can be determined from the above scalings as done in the next section.

A criterion for complete laser absorption
Given the final electron energy behind the laser beam, one can infer the minimum plasma areal density, * s = Z n L a abs abs (where n a is the atomic density, Z * is the final ionization degree and L abs is the absorption gas length), that yields complete laser absorption. In principle, the ionization dynamics of the gas should be included in the wakefield model of section 2. Since, for simplicity, we have assumed a plasma of fixed ionization degree, a knowledge of the dependence of Z * with the laser-plasma parameters is required. Here we propose a simple model that predicts, for a few-cycle laser pulse, the final ionization state as a function of the laser intensity I L , the atomic number Z and the atomic density n a .
The coupled dynamics of the variously charged ion species is ruled by the following set of continuity equations: where n i j is the density of the jth charged ions and n j is the ionization rate by an oscillating laser field [93][94][95], given by (in atomic units): . Following [97], the final ionization state can be deduced as follows. Starting from a weak laser intensity, the first atomic shell can be considered to be fully ionized over one laser cycle if the ionization rate is comparable with the laser frequency, i.e., n » 1 (in physical units), which is plotted in figure 7 as a function of Z * for various gases. Note that the ionization state reached after N laser cycles can also be estimated by solving n = ( ) E N 1 j . Yet, the threshold intensity, varying as ( ) N log , is weakly sensitive to this correction. From the knowledge of the electron density and mean energy following the laser interaction and phase mixing in the plasma wakefield, we can infer the minimum electron areal density, s abs , that yields full absorption of the laser pulse. The ion kinetic and electrostatic field energies being negligible in this regime, the laser energy is essentially converted into electron internal energy. The minimum areal density required for complete laser absorption can therefore be approximated as where t L is the laser pulse duration. Substituting g á ñ from equation (9), one can calculate s abs in terms of a 0 and * º n Z n a 0 , assuming a constant laser amplitude durant 10 laser cycles (i.e., t p = 20 L in normalized units). Owing to the strong electron heating occurring at low plasma densities and high laser intensities, s abs decreases with decreasing plasma density. This somewhat counterintuitive result holds insofar as the lengths required for optical shock formation and phase mixing are much smaller than the gas length. Obviously, the minimum gas length for full laser absorption is always a decreasing function of the electron density, as illustrated in figure 8(b).

2D simulations with realistic gas profiles
Multidimensional processes stemming from the finite size of the laser focal spot are likely to affect the electron energization mechanisms previously discussed. In the range of laser intensities under consideration, and for a focal spot larger than the relativistic plasma skin depth, p g  w n 2 L 0 , the pulse may break up into smallscale filaments [98,99]. In the small focal spot limit,  p g w n 2 L 0 , the plasma electrons can be transversely expelled by the laser ponderomotive force, leading to ion-filled cavities in the wake of the beam [21,23].
To explore these effects, we performed 2D Cartesian PIC simulations using a typical gas profile obtained at LOA and presented in figure 9(a). In figure 9(b), the electron areal densities for fully ionized Helium ( * = Z 2, blue solid line) and partially ionized Argon ( * = Z 16, black solid line) gases are plotted against the distance (y) from the nozzle. Note that a fully ionized Ar gas would require laser intensities   x . This decrease probably originates from the transverse breakup of the laser into three subbeams. The ponderomotive expulsion of the plasma electrons that results from the filamentary laser field is clearly seen in figure 10(d). Interestingly, such a laser breakup is not observed in the weaker focusing conditions of figures 10(a), (b). Figures 10(b), (d) show the development of electron cavitation, which tends to modify the scalings (of the electrostatic potential and of the electron energy gain) derived above in a 1D geometry. By contrast, in a 2D/3D geometry and in case of laser filamentation or cavitation, the electrons are transversely ejected while being pushed forward by the ponderomotive force and subsequently drawn back by the ions [23].As a result, they  reach weaker longitudinal momenta than predicted in the 1D framework of section 2. Note that their backward motion may be affected by the magnetic field induced by the forward-moving hot electrons, which is predicted to scale as g á ñ  B n z 0 [76]. While this field may help confine the hot electrons in the laser wake, we do not expect that the final mean electron energy will be significantly modified compared to the 1D picture. ), the mean electron energy is 4 times higher than the ponderomotive scaling.
One should note that the relativistic nonlinearities responsible for laser focusing and filamentation effects in our 2D simulations are probably underestimated in our 2D simulations compared to a realistic 3D configuration [102,103]. The 3D radial focusing or filamentation, depending on the initial beam waist, can substantially increase the laser intensity, and hence modify the nonlinear plasma wakefield. A more refined model for the electron energization should also take account of the 3D laser beam deformation and the quasistatic fields induced by partial electron depletion in the laser channel. Working out such a sophisticated model is well beyond the scope of the present paper.
Finally, to assess the robustness of our criterion for full laser absorption, we also conducted a simulation in an Ar gas with a small laser spot size, p = w 12 L ( l 6 0 ). The other parameters are identical to those described above. Figure 12(a) shows that the laser-driven electron density perturbations extend only over the first half of the plasma length. This means that the laser has been fully absorbed, as predicted from figure 9(b). By contrast, and, and as expected, figure 12(b), which is recorded at time t=6750 (i.e., much later than figure 10(d)), evidences partial laser transmission through the He gas. Note that the thin electron density shell that extends up to  x 3000 is the signature of an electrostatic shock accompanied by a magnetic dipole vortex [76].  The above condition for efficient laser absorption may be used in the future to determine the occurrence and/or efficiency of the ion acceleration mechanisms (TNSA or CSA) resulting from the electron energization.

Conclusions and perspectives
We have revealed an efficient, as yet overlooked, electron heating mechanism induced during the propagation of an ultra-intense laser pulse through a near-critical plasma. This mechanism arises once the laser pulse has developed a steep ( l 0 ) front profile, and results from the interaction between the bulk electrons and a population of hot electrons originating from the breaking of the laser-driven wakefield. The former experience a strong backward acceleration before being thermalized. Mean thermal energies several times higher than the ponderomotive scaling can be reached-of the order of~100 MeV and beyond-at laser field strengths a 50 100 0 and plasma densities~n 0.1 0.2 c . We have worked out a simple model, based on modified AP equations, which predicts the final mean electron energy in the ( ) a n n , e c 0 parameter space, in good agreement with extensive PIC simulations. From these results, we have derived a practical criterion for efficient laser absorption through an initially neutral gas jet, taking account of field-induced ionization. Our study is of particular interest for applications necessitating strong electron energization. Such is the case of laser-induced electrostatic shocks, which are believed to form when driving strong electron pressures within inhomogeneous plasma profiles.