Generation of neutral and high-density electron–positron pair plasmas in the laboratory

Electron–positron pair plasmas represent a unique state of matter, whereby there exists an intrinsic and complete symmetry between negatively charged (matter) and positively charged (antimatter) particles. These plasmas play a fundamental role in the dynamics of ultra-massive astrophysical objects and are believed to be associated with the emission of ultra-bright gamma-ray bursts. Despite extensive theoretical modelling, our knowledge of this state of matter is still speculative, owing to the extreme difficulty in recreating neutral matter–antimatter plasmas in the laboratory. Here we show that, by using a compact laser-driven setup, ion-free electron–positron plasmas with unique characteristics can be produced. Their charge neutrality (same amount of matter and antimatter), high-density and small divergence finally open up the possibility of studying electron–positron plasmas in controlled laboratory experiments.

E lectron-positron (e À /e þ ) plasmas are emitted, in the form of ultra-relativistic winds or collimated jets, by some of the most energetic or powerful objects in the Universe, such as black holes 1,2 , pulsars 3 and quasars 4 . These plasmas are associated with violent emission of gamma-rays in the form of short-lived (milliseconds up to a few minutes) bursts, which are among the most luminous events ever observed in the Universe. These phenomena represent an unmatched astrophysical laboratory to test physics at its limit and, given their immense distance from Earth (some more distant than several billion light years), they also provide a unique window on the very early stages of our Universe [5][6][7] . Arguably, one of the most intriguing questions is how these gamma-ray bursts are produced. It is generally accepted that gamma-ray bursts should arise from synchrotron emission of relativistic shocks generated within an electron-positron beam 8,9 . This radiative mechanism requires a strong and long-lived (t ! 1; 000o À 1 p , with o p being the electron-positron plasma frequency) magnetic field; however, Weibel-mediated shocks generate magnetic fields that should decay on a fast timescale ðt ' o À 1 p Þ due to phase-space mixing 9 . Also, diffusive Fermi acceleration, a proposed candidate for the acceleration of cosmic rays 9 , requires magnetic field strengths that are much higher than the average intergalactic magnetic field (CnT) 10 . These and other questions could be addressed by ad hoc laboratory experiments; however, the extreme difficulty in generating e À /e þ populations that are dense enough to permit collective behaviour 11,12 is still preventing laboratory studies and the properties of this peculiar state of matter are only inferred from the indirect interpretation of its radiative signatures and from matching numerical models. The intrinsic symmetry between negatively charged (e À ) and positively charged (e þ ) particles within the plasma makes their dynamics significantly different from that of an electron-ion plasma or from a purely electronic beam. In the first case, the mass symmetry of the oppositely charged species induces different growth rates for a series of kinetic and fluid instabilities 13 , and significantly affects the possibility of generating acoustic or drift waves. In the second case, the overall beam neutrality forbids the generation of current-driven magnetic fields that would hamper the onset of transverse instabilities. Different schemes have been proposed for the laboratory generation of e À /e þ plasmas: in large-scale conventional accelerators, the possibility of recombining high-quality electron and positron beams via magnetic chicanes 14 is envisaged and a different approach is foreseen in confining low-energy positrons using radioactive sources with Penning traps 11,15 . The proposed APEX experiment 12 builds on this idea, accumulating a large number of positrons in a multicell Penning trap, before injection into a stellarator plasma confinement device. The major challenge of these schemes is the recombination of these separate electron and positron populations. Alternative schemes have been proposed in which electrons and positrons are generated in situ [16][17][18][19][20][21] , thus avoiding the aforementioned recombination issues. Despite the intrinsic interest of these results, the low percentage of positrons in the electron-positron beam (of the order, if not o10%) and the low-density reported (collision-less skin depth much greater than the beam size, forbidding plasmalike behaviour) prevent their application to the laboratory study of e À /e þ plasmas. All these previous experimental attempts have thus not been able to generate e À /e þ beams that present charge neutrality and a plasma-like behaviour, both fundamental prerequisites for the laboratory study of this state of matter 14 .
We report here on the first experimental evidence of the generation of a high-density and neutral electron-positron plasma in the laboratory. Its high density n e À =e þ ' 10 16 cm À 3 À Á implies that the collision-less skin depth in the plasma is smaller than the plasma transverse size effectively allowing for collective effects to occur. These characteristics, together with the charge neutrality, small divergence y e À =e þ % 10 À 20 mrad À Á , and high average Lorentz factor (g AV E15 with a power-law spectral distribution, comparable to what observed in astrophysical jets 22 ) finally open up the possibility of studying the dynamics of e À /e þ plasmas in a controlled laboratory environment.

Results
Experimental setup. The experiment (shown schematically in Fig. 1a) was carried out using the ASTRA-GEMINI laser system at the Rutherford Appleton Laboratory 23 , which delivered a laser beam with a central wavelength l L ¼ 0.8 mm, energy on target E L E14 J and a duration of t L ¼ 42±4 fs. An f/20 off-axis parabola focussed this laser beam (focal spot with full-width half-maximum (27±3 mm) containing B60% of the laser energy, resulting in a peak intensity of C3 Â 10 19 W cm À 2 ) onto the edge of a 20-mm-wide supersonic He gas jet doped with 3.5% of N 2 . A backing pressure of 45 bar was found to be optimum in terms of maximum electron energy and charge of the accelerated electron beam as resulting from ionization injection 24,25 in the  gas jet. Optical interferometry of the laser-gas jet interaction indicates this gas-pressure to correspond to a plasma density of n pl ¼ (6.0±0.2) Â 10 18 cm À 3 . This interaction produced a reproducible electron beam (shot-to-shot fluctuation in charge and maximum energy below 10%) with a broad spectrum with maximum energy of the order of 600 MeV, full-width halfmaximum divergence of 2 mrad and an overall charge of (0.3 ± 0.1) nC, corresponding to (1.9 ± 0.6) Â 10 9 electrons (see Fig. 1b for typical electron spectra and their average). This electron beam was then directed onto a Pb solid target of different thicknesses covering multiples of the material's radiation length (d ¼ 0.5, 1, 1.5, 2, 2.5, 3 and 4 cm, given that the radiation length for Pb is L rad E0.5 cm (ref. 21)). The electrons and positrons escaping from the rear side of the target were then separated and spectrally resolved by a magnetic spectrometer. The details of this detector can be found in the Methods section.
Experimental results. A scan in target thickness was performed in multiples of its radiation length and the obtained positron spectra, each resulting from an average over five consecutive shots, are depicted in Fig. 2 (see Fig. 1c for the raw signal recorded on the LANEX screen for d ¼ 0.5 cm). All spectra are in good agreement with the ones resulting from matching simulations using the Monte Carlo scattering code FLUKA, which accounts for electromagnetic cascades during the passage of an electron beam through a solid target 26 (see Methods section). A maximum positron energy of E MAX ¼ 600 MeV is obtained for dEL rad (that is, 5 mm; Fig. 2a), whereas a maximum positron yield is obtained for dE2L rad . For thicker targets, the maximum energy gradually decreases as it should be expected due to increased probability of energy loss during the propagation of the generated positrons through the rest of the solid target. For a similar reason, a thicker solid target allows a lower number of electrons and positrons to escape it. This is quantitatively shown in Fig. 3, which depicts the measured number of electrons and positrons (energy exceeding 120 MeV; see Methods section) at the exit of a solid target, as a function of its thickness.
In order to quantitatively explain the observed trends, we have employed a simple analytical model for a quantum electrodynamic cascade that only includes the emission of photons by electrons and positrons via bremsstrahlung 27 and the creation of an electron-positron pair by a photon 28 , both processes occurring in the field of a heavy atom. We thus neglect additional energy losses as resulting, for instance, from Compton scattering with the electrons of the atoms and from the ionization of the atoms themselves (see Methods section). This model is able to  ARTICLE qualitatively reproduce the experimental trends (dashed green curves in Fig. 3), provided that a constant re-scaling factor of about 0.75 is adopted for the absolute yield of both the electrons and positrons. This overestimate is easily understood, as the semianalytical model does not take into account a number of energy loss mechanisms, such as Compton scattering and the ionization of atoms 29 . Once this re-scaling factor is applied, the analytical model reproduces the experimental data within a few per cent, clearly indicating that the only processes of bremsstrahlung and electron/positron pair production in the nuclear field are the dominant mechanisms leading to the generation of the detected electron/positron beam.
Let us now turn our attention to the total positron fraction in the leptonic beam as a function of the target thickness (plotted in Fig. 3c). For dEL rad , the positrons account for B8-10% of the overall beam owing to the fact that most of the primary electrons are able to escape the target (consistently with the results reported in ref. 21). However, as we increase the target thickness, this ratio increases up to a point where the positrons account for almost 50% of the leptonic jet (dZ2.5 cm; Figs 3c and 4 for the overall charge imbalance in the leptonic beam and its simulated spatial distribution, respectively). In this case, not only the integrated number of electrons and positrons is similar but also their spectrum (Fig. 2c), further indication that almost all the electrons and positrons escaping the target arise from pair production. A positron percentage in the beam of the order of 50% is preserved also if the target thickness is increased; however, we will focus our attention only on d ¼ 2.5 cm, since it provides the highest density of the neutral e À /e þ beam. Simulations confirm that the majority of positrons are generated with energies of the order of a few MeV following a Jüttner-Synge distribution, which is commonly assumed for relativistic thermalized plasmas 30 (see, as an example, the inset in Fig. 2a). We thus refer to the experimentally measured number of e À and e þ E e AE 4120 MeV ð Þwith the subscript N EXP , whereas we will refer to their simulated number E e AE 42m e c 2 % 1 MeV ð Þ with the subscript N FLUKA . For d ¼ 2.5 cm, we thus have N e À EXP % N e þ EXP % 3Â10 7 and N e À FLUKA % N e þ FLUKA % 1:2Â10 9 (Fig. 5a). Taking the appropriate moment of the distribution function, the averaged Lorentz factor of the beam is typically of the order of a few tens (g AV E15 for d ¼ 2.5 cm). FLUKA simulations indicate a divergence of the beam to be energy dependent in a range of 5-20 mrad (ref. 31).
It must be pointed out that the propagation of an ultrarelativistic electron beam through a high-Z solid target can only asymptotically give a perfectly neutral e À /e þ beam. Additional scattering mechanisms with the atomic electrons, such as Compton, Moller and Bhabha scattering, will in fact slightly increase the electrons number, especially at low energies. FLUKA simulations take all these processes into account and indeed predict an average percentage of positrons, for d ¼ 2.5 cm, of 46%. The discrepancy between electron and positron number is exclusively at low energies (E e AE ' 5 MeV; Fig. 4d). Most importantly, the electron and positron populations present very similar spatial distributions (Fig. 4a,b) leading to an almost uniform positron percentage in the e À /e þ beam (between 45 and 49%; Fig. 4c). As we shall see later, this slight charge imbalance does not affect the plasma dynamics, which can then be effectively considered to be neutral.
A fundamental requisite for the laboratory study of e À /e þ plasmas is that they must present collective behaviour in their dynamics. Collective (that is, plasma-like) effects are likely to occur in the beam only if its transverse size D B is larger than the collision-less skin depth (l skin Cc/o prop , with o prop being the relativistic plasma frequency). The beam density is determined by the temporal duration of the beam (that relates to its longitudinal extent) and its transverse size. The primary electron beam exits the gas jet with a typical temporal duration comparable to half the plasma period within the gas 32 : t pl C(13.0±0.3) fs. The semianalytical model for the quantum cascade inside the Pb indicates an average temporal spreading across different spectral components of the beam of the order of 1-3 fs, resulting in a beam duration of t e À =e þ ' 15 AE 2 fs. As intuitively expected, the lower energy electrons and positrons will escape the solid target in a wider area if compared with their higher energy counterparts. FLUKA simulations confirm this expectation and indicate, for d ¼ 2.5 cm, a maximum transverse size of the beam of the order of D B C200 ± 30 mm. For these parameters, we thus obtain a particle density in the laboratory reference frame of the order of n e C(1.8 ± 0.7) Â 10 16 cm À 3 , implying a beam proper density of n prop ¼ n e /g AV C(1.5±0.5) Â 10 15 cm À 3 (Fig. 5b). The relativistically corrected collision-less skin depth of the beam is thus l skin Cc/o prop (160±30) mm. This value is smaller than the beam transverse size, indicating that the generated particle beam is a neutral e À /e þ plasma. It is interesting to note that the occurrence of collective behaviour (that is, the situation in which D B /l skin Z1) does not depend on the beam transverse size D B since, based on the considerations presented above, it can be expressed as: D B =l skin % 4:1Â10 À 4 ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi N= g AV t pl ½fs À Á q % 1:4 for our experimental parameters (here N indicates the overall number of leptons in the beam).

Discussion
The presented characteristics of the e À /e þ plasmas generated in our experiment are appealing for the laboratory study of the dynamics of this exotic state of matter. As an example, a particularly active area of research in this direction is the determination of the growth and evolution of kinetic instabilities, which are extensively modelled in order to interpret peculiar astrophysical observations such as the emission of gamma-ray bursts [33][34][35][36] . It is widely accepted that these ultra-bright bursts result from synchrotron radiation generated via relativistic shocks triggered during the propagation of an electron-positron beam through the low-density intergalactic medium 37 . This scenario is now reproducible in a laser-driven experiment in which the photoionized residual low-density gas inside the target chamber 38 can act as the background electron-ion plasma. In this case, the growth rate for transverse instabilities can be estimated as: 13), with b spread and o ei being the velocity spread of the e À /e þ beam and the plasma frequency of the e À -ion plasma, respectively. It is worth noticing that in the ultra-relativistic case, the weak dependence of the growth rate on the beam velocity spread significantly relaxes constraints on the spectral shape of the electron-positron beam. We can assume o ei E1.5 Â 10 12 Hz (n ei E6 Â 10 14 cm À 3 as resulting from full photoionization of the background gas) and b spread E0.1 (b ¼ 0.87 and bE1 for a 1 MeV and a 500 MeV particle, respectively). We thus have G TR ¼ 5 Â 10 11 Hz for g AV ¼ 15 implying a typical time for the instability to grow of the order of 2 ps. Numerical simulations indicate, in the initial instants of the instability, that up to 10% of the average particle energy in the beam can be transformed into electromagnetic fields in the plasma implying fields with an amplitude of the order of the megagauss; once saturation is reached, this value drops to B1% (ref. 13). It is worth noticing that this is similar to what expected for gamma-ray bursts (0.1-1%; ref. 39). This timescale and field amplitude are within reach of plasma radiography techniques such as proton imaging 40 , a highly encouraging factor for the application of these plasmas for laboratory astrophysics.
In order to check the validity of our estimates, we have carried out three-dimensional (3D) particle-in-cell (PIC) simulations using the PIC code OSIRIS 41,42 (see Methods section). Simulation results are illustrated in Fig. 6. During its propagation through a denser e À -ion plasma, the e À /e þ is subject to the Weibel/ current filamentation instability leading to the formation of electron and positron filaments with thicknesses of the order of the beam skin depth. The electron and positron filaments spatially separate from each other leading to net localized currents and the generation of the corresponding azimuthal magnetic field structures with maximum amplitudes of the order of 40 T in the middle of the bunch. At early times, the simulations show that the transverse scale length of the filaments is even shorter than the initial beam skin depth. To further understand the impact of charge neutrality on the instability onset, additional 3D simulations were performed using a purely electronic bunch of same characteristics. In this case, the electron bunch generates plasma wakefields, and neither filamentation of the beam (insets in Fig. 6c) nor the generation of strong magnetic fields (inset Fig. 6d) are observed. These results corroborate the expectation that current filamentation instability growth can be controlled by changing the beam overall total charge and it is maximized for a purely neutral e À /e þ plasma. Finally, we performed an additional 3D PIC simulation devoted at studying whether a slight charge imbalance in the e À /e þ plasma could result in a change in the plasma dynamics if compared with the idealized perfectly neutral plasma scenario. We have thus maintained exactly the same conditions as the other simulation, with the only difference that now the positron account for 45% of the plasma population, in order to match our experimental findings more closely. The obtained spatial distribution of the e À /e þ plasma after propagation through the background electron/ion plasma is shown in Fig. 6e,f, indicating essentially no difference if compared with the purely neutral case. This statement is corroborated by the growth of magnetic fields due to Weibel instability. This is plotted in Fig. 6g that shows virtually the same magnetic field growth for the purely neutral case (blue line) and for the slight charge imbalance (red). For the point of view of studying electron-positron plasma dynamics in the laboratory, the e À /e þ plasma generated in our experiment is virtually indistinguishable from the idealized purely neutral beam.
On the other hand, the beam might also be susceptible to longitudinal instabilities 34,43 , which would induce a broadening of the e À /e þ spectrum and generation of strong fields in the background plasma. For d ¼ 4 cm (neutral beam), the measured electron and positron spectra are indeed flatter than the ones predicted by FLUKA, which does not include collective behaviour of the beam particles during propagation through the background e À -ion plasma (Fig. 2c). For d ¼ 0.5 cm (highly charged beam), simulations and experiments agree much more closely. The spectral flattening may also be produced by kinetic self-focusing of the beam 44,45 .
In conclusion, we have reported on the first creation of a neutral electron-positron plasma in the laboratory. Its overall charge neutrality and plasma-like behaviour are an absolute novelty in the field of experimental physics and, in conjunction with the small divergence and high energy of these plasmas, finally allow for the laboratory study of this unique state of matter.

Methods
The electron-positron spectrometer. The magnetic spectrometer comprised a pin-hole entrance with a diameter of B15 mm through 5 cm of plastic followed by 5 cm of lead. This plastic-lead wall was indeed necessary in order to shield the particle detectors from noise generated during the electron beam impact onto the solid target. After this, a dipole permanent magnet (B ¼ 0.8 T, length of 10 cm) was inserted to spectrally resolve the electrons and the positrons, which were recorded by two LANEX screens 46 . This arrangement allowed us to resolve particle energies from 120 MeV to 1.2 GeV. The LANEX screens were cross-calibrated using absolutely calibrated Imaging Plates 47 . The small difference in stopping power (of B2%; ref. 48) between electrons and positrons was taken into account in calibrating the LANEX screens. Every electron or positron spectrum shown in the manuscript results from an average over five consecutive shots. The energy resolution of the spectrometer can be approximated in the ultra-relativistic limit, as: Where D s is the distance from the source to the magnet entrance, D l is the distance from the entrance of the magnet to the detector (1 m), R L EE/(ecB) is the radius of curvature of the particle with energy E and charge e in the magnetic field B, y s ¼ 15 mrad is the angular acceptance of the detector, and L m (10 cm) is the length of the magnet. For the energies of interest in our experiment (120rE[MeV]r300), the energy resolution is between 10 and 20%.
FLUKA simulations. FLUKA is a nuclear physics Monte Carlo scattering code that accounts for electromagnetic cascades during the passage of an electron beam through a solid target 26 . The numerical model for the quantum electromagnetic cascade is routinely checked and constantly improved to take into account any refinement in cross-section measurements in conventional accelerators. As an input for the simulation, we assume an electron beam with the spectral shape depicted in Fig. 1b (brown solid line), 2 mrad full-width half-maximum divergence and 10 mm radius source size. The electron beam then interacts with a lead target of different thicknesses and 1 cm transverse size, placed 1 cm downstream of the electron beam source. Iterations (10 6 ) were used in order to achieve a good statistical representation in the Monte Carlo method. Every numerical result reported originates from an average over five identical runs in order to minimize any stochastic error arising from the random seed generator of the code. The results of the simulations, obtained in units of particles per initial electron, were then rescaled with the measured number of primary electrons, giving a good quantitative agreement with the experimental data.
Semi-analytical model for the quantum cascade. We assume a quantum electrodynamics cascade shower involving only electrons, positrons and photons at energies much larger than the electron rest energy m (units with h ¼ c ¼ 1 are assumed hereafter). We thus neglect additional electron and positron energy losses as resulting, for instance, from Compton scattering with the electrons of the atoms and from the ionization of the atoms themselves. The only processes to be included in the kinetic equations are thus the emission of photons by electrons and positrons via bremsstrahlung and the creation of an electron-positron pair by a photon, both processes occurring in the field of a heavy atom. By setting the target thickness d in units of the radiation length L rad , that is, c ¼ d/L rad , the electron/positron distribution functions f±(E,c) and the photon distribution function f g (E,c) satisfy the kinetic equations 29 : where the functions with m 0 ¼ 7/9 À b/3 and b ¼ 1/18 log(183/Z 1/3 ), are related to the cross-section of bremsstrahlung and pair photo-production in the field of a heavy atom with atomic number Z (see ref. 29 for details). By numerically solving these equations, we are able to reproduce the experimental trends well (dashed green curves in Fig. 3 of the manuscript), provided that a constant re-scaling factor of 0.75 is adopted for the absolute yield of both the electrons and positrons. The overestimation of the experimental results by this simplified model is easily understood, as the latter does not take into account a number of braking mechanisms, such as Compton scattering for photons and the ionization of atoms for electrons and positrons 29 .
On the one hand, braking mechanisms such as ionization affect essentially relativistic electrons and positrons in the same way 49 . On the other hand, however, our analytical model cannot predict charge asymmetries brought in by injection of atomic electrons in the cascade following these scattering processes. This is the reason why, for a target thickness of 2.5 cm, our semi-analytical model predicts a 50% percentage of positrons in the leptonic beam, whereas our FLUKA simulations indicate a positron percentage of the order of 46%.
Starting from a simple model, where each electron/positron (photon) after a radiation length emits a photon (transforms into an electron-positron pair) with half of the energy of the initial electron (with the electron and positron sharing half of the energy of the initial photon), it can also be shown that the maximum yield of positrons with an energy exceeding E can be estimated to occur for a target thickness d opt BL rad log(hE e i/E)/log(2) (ref. 29), where hE e i is the average energy of the initial electron distribution (Fig. 1b). In our case, it results hE e iE456 MeV and d opt B1.1 cm (C1.96 L rad ) in good agreement with the experimental results.
The PIC simulations. The simulations were performed with the fully relativistic, massively parallel, PIC code OSIRIS 41,42 . OSIRIS has been extensively used to explore relativistic beam plasma interaction scenarios, and has been widely applied to model the Weibel instability in various configurations (see, for instance, refs 13,14,41,42,50). In OSIRIS, the electric and magnetic fields are defined in a grid. The trajectory of each simulation particle is determined through the relativistic equations of motion by interpolating the grid fields to the position of the particle. Current density is deposited onto the grid, and used to advance the electric and magnetic fields through Maxwell's equations discretized using a finite-difference scheme. In this section, we give the numerical parameters for the simulations. Simulations used a moving window with dimensions 1.5 Â 100 Â 100 (c/o p ) 3 divided into 75 Â 1,000 Â 1,000 cells with 2 Â 1 Â 1 particles per cell for plasma electrons and for beam particles. Here o p is the plasma frequency of the background electron-proton plasma, which has a density of n ei ¼ 10 16 cm À 3 . A charge-neutral beam constituted by electrons and positrons was initialized at the entrance of the plasma. The density profile for electrons and positrons is given by n b ¼ n b0 exp À x 2 s 2 x À r 2 s 2 r where n b0 ¼ 10 n ei ¼ 10 17 cm À 3 , s x ¼ 0.22 c/ o p ¼ 11.7 mm and s r ¼ 10 c/o p ¼ 530 mm are the bunch peak density, length and transverse waist, respectively. The particles' Lorentz factor is initialized to be g e À ¼ g e þ ¼ 700.