Detailed characterization of a laboratory magnetized supercritical collisionless shock and of the associated proton energization

The acceleration of energetic charged particles by a collisionless magnetized shock is a ubiquitous phenomenon in astrophysical environments, in which the most energetic particles are the ultra-high-energy cosmic rays (UHECRs) accelerated in the interstellar medium (ISM).^1,2 In these environments, the sources of collisionless dissipation are self-generated electromagnetic fields, resulting from kinetic instabilities such as the Weibel instability. Particles are also accelerated in the Solar System owing to collisionless magnetized shocks generated by the interaction of the solar wind with planetary magnetospheres^3,4 and, at larger distances, with the ISM.⁵ Here, the source of collisionless dissipation is the pre-existing global electromagnetic structure. This will also be the case for the experiment detailed here, where we apply a global strong magnetic field to a laser-ablated fast plasma. Since these shocks usually have a magnetosonic Mach number M_ms = v_sh/v_ms ≳ 2.7 (where v_sh is the shock velocity, $vms=Cs2+vA2$ is the magnetosonic velocity, and C_s and v_A are the ion sound velocity and Alfvénic velocity, respectively), they belong to the so-called supercritical regime,^6,7 which means that these shocks are not maintained by classical dissipation alone. To help maintain a shock, the additional channel for the expulsion of energy is via reflection of particles back upstream.⁸ 

A variety of acceleration mechanisms have been invoked as a means to transfer energy from shock waves to particles, including shock surfing acceleration (SSA), shock drift acceleration (SDA), and diffusive shock acceleration (DSA). DSA requires a high initial energy before further acceleration,⁹ thus raising the so-called “injection problem,”¹⁰ while SSA and SDA are believed to be responsible for generating the pre-accelerated seed particles, i.e., for the initial process of acceleration from thermal energies. Although it is still under debate whether SSA or SDA dominates the pre-acceleration process in various collisionless shock environments,^11–13 we can distinguish them by the following two aspects. On the one hand, in SSA, charged particles are first reflected at the shock front (owing to the cross-shock potential electric field), and they then surf along the shock front against the convective electric field (E = −v × B), thus gaining energy, while in SDA, charged particles drift (owing to the magnetic field gradient at the shock front) along the convective electric field and then gain energy.¹⁴ On the other hand, SSA requires a shock width that is thin compared with the Larmor radius of the charged particles, while SDA needs the opposite (so that the charged particles can gyrate and drift within the shock layer).^9,15

However, because of the immense spatial scales involved in collisionless phenomenon (e.g., the mean free path is λ_mfp ∼ 1 AU in the Solar System), only a very small sampling of the shock formation and dissipation mechanisms can be realized. As a result, we still do not have a full understanding of the formation and evolution of collisionless shocks, and the question of the effectiveness and relative importance of SDA and SSA is still widely debated in the literature.¹⁶ To further our understanding, laboratory experiments (and their simulations) have proved to be an effective tool, providing highly resolved, reproducible, and controllable multidimensional datasets that can complement astrophysical observations.^17,18 Below, we will briefly review the investigation of collisionless shocks via laboratory experiments.

The route that up until now has been most explored in the laboratory is to produce a shock (mediated by the Weibel filamentation instability) by colliding two ablative unmagnetized flows driven by high-energy nanosecond lasers. This setup has yielded promising results at the Omega Laser Facility^19–21 and the National Ignition Facility (NIF),^22,23 as well as at many other laser facilities all over the world.^24–26 Recently, experiments on collisionless shocks in plasma flows in which there was a significant self-generated magnetic field showed, for the first time, the formation of magnetized collisionless shocks, with the generation of Weibel instability and observation of electron acceleration in the turbulent structure.²⁷ Most recently, the dynamics of the ion Weibel instability have been characterized by local quantitative measurements of ion current filamentation and magnetic field amplification in interpenetrating plasmas via optical Thomson scattering (TS).²⁸ What is more, the generation of subrelativistic shocks, together with relativistic electron acceleration, has been demonstrated to be within the reach of larger-scale, NIF-class laser systems.²⁹ 

Another setup relies on a plasma expanding into a preformed ambient magnetized secondary plasma. Thanks to the magnetization, the target ions create a collisionless magnetic piston that accelerates the ambient plasma to super-Alfvénic velocity, thus creating a high-Mach-number shock with velocity of the order of 1000 km/s.^30–33 Recently, Schaeffer et al.³⁴ have been able to make significant progress in characterizing the formation of collisionless shocks in terms of ion and electron density and temperature, as well as electric and magnetic field strengths, as functions of time at OMEGA.

At the LULI laser facility at the École Polytechnique (France), collisionless shock waves and ion-acoustic solitons have been investigated by proton radiography.³⁵ Moreover, significant electron preheating via lower-hybrid waves has also been achieved in laboratory laser-produced shock experiments with strong magnetic fields, providing a potential mechanism for the famous “injection” problem.³⁶ Additionally, at the VULCAN laser facility at the Rutherford Appleton Laboratory, temporally and spatially resolved detection of the formation of a collisionless shock has been achieved.³⁷ 

In contrast to the above schemes, novel setups have been used with ultra-high-intensity lasers. For example, at the XingGuang III laser facility at the Laser Fusion Research Center in China, using a short (2 ps) intense (10¹⁷ W/cm²) laser pulse, an electrostatic (ES) collisionless shock, together with the filaments induced by ion–ion acoustic instability, could be observed via proton radiography.³⁸ 

In our experimental campaigns at JLF/Titan and LULI2000,³⁹ we investigated shock formation combining laser-produced plasmas, a background medium, and a strong ambient magnetic field (as detailed below). We chose to have an expanding plasma drive a shock into an ambient gas in the presence of a strong external magnetic field. In contrast to Schaeffer et al.,³² in our setup, the expanding plasma and the magnetic field were decoupled as the higher-Z piston evacuated the magnetic field and was thus unmagnetized. This also allowed us to simultaneously have a highly magnetized ambient plasma (with homogeneous and steady magnetic field) and a high-β piston (as can be seen in Table I, the plasma thermal β of the piston is β_ther ≡ P_ther/P_mag ∼ 14.0). Moreover, since our magnetic field strength was more than two times higher,⁴⁰ reaching 20 T compared with the 8 T in Schaeffer et al.,³² we were able to decouple the electrons more strongly from the ions,⁴¹ and the shock was able to fully separate from the piston, which is crucial for its characterization.⁴² As a result, we have been able to characterize the plasma density and temperature, as well as the electric field developed at the shock front, and, more importantly, observe strong nonthermal accelerated proton populations for the first time.

In this paper, we will first show, in Sec. II, that laboratory experiments can be performed to generate and characterize globally mildly supercritical, quasi-perpendicular magnetized collisionless shocks, and detail their characteristics. Then, in Sec. III, we will detail three-dimensional (3D) magnetohydrodynamic (MHD) simulations reproducing the laser-driven piston generation and the subsequent shock formation process. In Sec. IV, with the parameters characterized in the experiment, we will report the results of kinetic particle-in-cell (PIC) simulations, which pinpoint that shock surfing acceleration (SSA) can be effective in energizing protons from the background plasma to 100 keV-level energies.

The experiments were performed at the JLF/Titan (LLNL, USA) and LULI2000 (France) laser facilities with similar laser conditions but using complementary diagnostics, with the latter being mostly a consequence of the availability of different auxiliary laser beams at each facility.

In the experiment at JLF/Titan, as shown in Figs. 1(a)–1(c), the collisionless shock was generated by sending a plasma, generated by having a high-power laser (1 µm wavelength, 1 ns duration, 70 J energy, and 1.6 × 10¹³ W/cm² on-target intensity) irradiating a solid target (Teflon, CF₂), into a low-density (10¹⁸ cm⁻³) H₂ ambient gas pulsed from a nozzle prior to the shot, and in the presence of a 20 T magnetic field that was homogeneous and steady state on the time scale of the experiment. As shown in Figs. 1(b) and 1(c), the magnetic field, created by a Helmholtz coil system,^40,43 was oriented along the y or z axis.

In the subsequent experimental campaign at LULI2000, we set up a similar system, shown in Fig. 1(d). The laser parameters, the ambient H₂ pressure, and the magnetic field strength were the same as in our previous experiment on Titan and detailed above. In this case, owing to geometrical constraints, the CF₂ target was tilted by 60° around the z axis and by 45° around the y axis, in order to allow the main laser to reach it, as well as to leave access to the Thomson scattering probe beam and collection path.

Using an interferometry setup,⁴⁴ the plasma electron density was recorded by optically probing the plasma (with a mJ, 1 ps auxiliary laser pulse). In Fig. 2, we present the overall electron density recorded in three different cases.

For the case with both ambient gas and magnetic field shown in Figs. 2(a) and 2(b), the laser irradiation induced the expansion of a hot plasma (the piston) that propagated along the x axis, and a collisionless shock was formed as a consequence of the plasma piston propagating in the magnetized ambient gas.³² We can clearly see both the piston front and the shock front (indicated by the orange and green arrows, respectively), and indeed they are well detached from each other, enabling us to characterize them separately.

A lineout of the plasma density is shown in Fig. 2(f), where the piston and shock fronts are also well identified by the abrupt density changes. The piston front is steepened by the compression of the magnetic field (see also below). Besides, we can clearly see a “foot” structure ahead of the shock front in the upstream (US) region for the cases with both ambient gas and magnetic field, indicating the formation of the magnetized shock.⁴⁷ 

By contrast, for the case with only a magnetic field but no ambient gas⁴⁵ shown in Figs. 2(c) and 2(d), owing to the lack of ambient gas, no collisionless shock is formed ahead of the piston. For the case with only ambient gas but no magnetic field in Fig. 2(e), again no shock is formed in the ambient gas. From the corresponding lineout in Fig. 2(f), it is clear that only a smooth plasma expansion into the ambient gas (the green dashed line) occurs.

To further characterize the piston, the x-ray ion emission of fluorine compressed within the expanding piston was measured by a focusing spectrometer with high spatial resolution (FSSR)⁴⁸ at both laser facilities. The FSSR was based on a spherically bent mica (2d = 19.9376 Å) crystal with a radius of curvature R = 150 mm. A spatial resolution of 100 µm/pixel was achieved along the plasma expansion. Image plates (Fujifilm TR BAS) were used as fluorescent detectors. The implemented scheme resulted in a 13–16 Å spectral range with a high resolution (λ/dλ > 1000). This range covers the spectral lines of fluorine: resonance H-like (2p–1 s transition) and He-like (3p–1 s, 4p–1 s, 5p–1 s, etc.) transitions, as well as dielectronic satellites to Ly_α. The diagnostic allowed us to measure electron density and temperature profiles of the piston expansion using a quasistationary approach.⁴⁹ The method is based on analyzing the relative intensities of spectral lines of the same charge state and also takes into account the recombining plasma with a “frozen” ion charge.

Figure 3(a) shows that the piston obviously encounters stronger hindrance in the case with both ambient gas (H₂) and magnetic field (B_z) (see the green diamonds), compared with the other cases, i.e., the case with only B_z (red dots) and the case with only H₂ (blue triangles). Being time-integrated, the FSSR provides an average on-axis volumetric density over the few tens of nanoseconds when the plasma density and temperature are high, rather than the time-resolved and line-integrated density measured by optical interferometry. These two diagnostics provide complementary views of the piston. We also see in Fig. 3(b) that the electron temperature in the case of B_z + H₂ becomes the highest at the piston front (between 4 and 7 mm), compared with the other cases. In addition, at the position of 4.5 mm, the evaluated electron density for the case of B_z + H₂ is around 2–3 × 10¹⁸ cm⁻³ and the electron temperature is about 65 eV, which are reproduced well by our FLASH simulations; see Figs. 8(a) and 8(b).

The single shock front was also probed with protons by measuring the local electric field. The probing protons (accelerated by the target normal sheath acceleration process⁵¹ from an auxiliary target and using the short-pulse arm of Titan) were sent parallel to the magnetic field, i.e., along the z axis, as shown in Fig. 1(a).

As shown in Fig. 4(a), we can clearly observe the same structures of the piston front and the shock front, consistent with those observed via the optical probe, as shown in Fig. 2.

The synthetic proton radiography shown in Fig. 4(c) is calculated from a particle tracing code, ILZ.⁵² ILZ is a test-particle code that uses a given 3D distribution of electric and magnetic fields to simulate the trajectories of the protons as they pass through and then ballistically propagate to the detector.

As an input for the tracing code, we used an electric field E_x with a bipolar spatial profile along the x axis with a hemispheric geometry in the xz plane (where z is the axis along which the proton beam propagates) in order to imitate the experimentally observed curvature of the shock front. The same dependence on the coordinates x and z is reproduced along y over a small thickness. The electric field has only a component along x, given by

where E₀ is the maximum field amplitude and L represents the width of the region affected by the electric field. $g(z)=−R/2+R2−z2$ represents the shift of E_x along x in order to take into account the hemispherical geometry in the xz plane and derives directly from the equation for a circumference of radius R centered at (−R/2, 0): (x + R/2)² + z² = R². The hemispherical structure has a radius R = 5.1 mm, which is estimated from the experimental proton radiography relative to that time (t = 5.0 ns after the laser impact), as shown in Fig. 4(a).

We adjust E₀ and L of the ILZ input in order to match the simulated proton dose with the experimental one. A comparison between these two is shown in Fig. 4(c), where the ILZ-simulated modulation was obtained with the field in Eq. (1) with E₀ = 4.33 MV/m and L = 0.12 mm and is shown in Fig. 4(b).

Note that since the protons are sent along the external magnetic field (i.e., in the z direction) and since the strength of the magnetic field components in other directions are orders of magnitude lower than the externally applied one, the Lorentz forces associated with these magnetic field components are much smaller than that induced by the E_x electric field associated with the shock front. As a result, we interpret the proton dose modulation as being caused by the electric field only. This assumption is supported by the fact that it yields a simulated proton deflection that is highly consistent with the one recorded in the experiment, as can be seen in Fig. 4(c). Also note that the integration time in a given film is too short for motion blurring of the moving shock front to be seen. Additionally, the amplitude of the electric field at the shock front inferred from the synthetic proton radiography is of the order of MV/m. We will compare it with the PIC simulation results and discuss these in detail in Sec. IV.

Furthermore, we compare the position of the shock structure seen in the electron density (via interferometry) with that in the electric field (via proton radiography) for the case with both an external magnetic field and ambient gas. For the former, we consider the point where the electron density has a sharp jump, as shown in Fig. 2(f); for the latter, we take into account the external edges of the proton dose accumulation. Figure 5(a) shows the evolution of the piston front and that of the shock front through both diagnostics (see the caption for details), and they clearly reveal the slowing down of the piston and shock fronts over the first few nanoseconds after the laser pulse. Note that when the target was not clearly visible in the radiography, i.e., for the series of points around 5 ns, we have made use of the interferometry results to shift all the points by the right amount, while keeping the distances between the piston and shock fronts constant. The original radiochromic films (RCFs) for the data points at various times are also shown in Figs. 5(b)–5(e).

With a second high-energy auxiliary laser (526.5 nm wavelength, 1 ns, 15 J) available at LULI2000, we were able to perform Thomson scattering (TS) off the electron and ion waves in the plasma (used in a collective mode⁵³ and analyzed by different spectrometers). As shown in Fig. 1(d), the collection of the scattered light was performed at 90° (along the z axis) with respect to the incident direction of the laser probe (the x axis). The light scattered off the ion (TSi) and electron (TSe) waves in the plasma was analyzed by means of two different spectrometers, set to different dispersions (3.1 mm/nm for TSi and 7.5 × 10⁻² mm/nm for TSe), which were coupled to two streak cameras (Hamamatsu for TSe, and TitanLabs for TSi, both equipped with a S-20 photocathode to be sensitive in the visible part of the spectrum, and both with typical 30 ps temporal resolution), allowing us to analyze the evolution of the TS emission in time. The scattering volumes sampled by the instruments were 120 µm along the x and y axes and 40 µm along the z axis for TSi, and 100 µm along the x and y axes and 40 µm along the z axis for TSe. The analysis of the Thomson-scattered light was performed by comparison of the experimental images (recorded by the streak cameras) with the theoretical curves of the scattered spectrum for coherent TS in noncollisional plasmas, with the instrumental function widths of 5.9 nm for the electron spectrometer and 0.12 nm for the ion spectrometer taken into account. Note that the TS laser probe induces some heating in the hydrogen ambient gas. With the estimate of the electron temperature heated by the TS laser obtained through inverse Bremsstrahlung absorption, we confirm that the upper limit of the TS-induced heating is only around 60 eV, which is significantly smaller than the level of temperatures we observe in the shock.³⁹ With the above TS diagnostics, we could access spatially and temporally resolved measurements of the plasma density and temperatures (electron and ion) in the upstream (US), as well as in the downstream (DS) region.

Figure 6 shows the TS measurements in the DS region compared with the shock front for cases with and without an external magnetic field. By comparing the experimental data profiles with the theoretical equation for the scattered spectrum for coherent TS in unmagnetized and noncollisional plasmas, with the instrumental function taken into account, we are able to retrieve the local electron number density, as well as the electron and ion temperatures.⁵⁴ For the case without a magnetic field (i.e., with only ambient gas), both TSe and TSi give n_e ∼ 2.0 × 10¹⁸ cm⁻³ and T_e ∼ 80 eV, and TSi also gives T_i ∼ 50 eV in the DS region, as can be seen in Figs. 6(a) and 6(c). However, for the case with B = 20 T, we see strong heating in the DS region, indicated by the higher temperatures, i.e., T_e ∼ 230 eV and T_i ∼ 250 eV, as can be seen in Figs. 6(b) and 6(d).

Note that with a multi-ionic approach (i.e., including C, F, and H ions), the theoretical scattering curves are calculated from the distribution functions of each species, weighted by their presence fraction. Owing to the natural damping of the other ion species modes, we observe only one mode in the ion-acoustic wave (IAW) signal.

With all the above diagnostics, we summarize the characterization of the piston and shock conditions in Table I. Note that for the collisionality, the local magnetic field strength is used. Specifically, the local magnetic field strength of the piston is taken from our FLASH simulation results (see Sec. III) and that of the shock is from our PIC simulations, while for the calculation of the Mach numbers for the shock, the parameters of the US region are used.

For observation of the nonthermal proton spectrum, we used a standard magnetic spectrometer, equipped with permanent magnets of 0.5 T strength and a pinhole. It was located close to the target (17.5 cm away) in order to maximize its collection efficiency, and it had its main axis along z, the main axis of the external magnetic field (in an alternate mode to performing TS). Note that having the spectrometer collection axis aligned with that of the magnetic field allows measurement of the ions energized out of the plasma,⁴³ which otherwise could not have been recorded, since they would have been deflected away by the 20 T large-scale magnetic field. We also used filters in order to eliminate the possibility that the signal observed in the dispersion plane of the spectrometer was originating from heavy ions rather than protons from the ambient gas.

That spectrometer had been calibrated precisely with a Hall probe and on many previous campaigns using filters to verify its energy dispersion. The protons were deflected by the magnetic field inside the spectrometer and landed after a short drift space onto imaging plates (of TR type), the detector used here. These detectors were absolutely calibrated.⁵⁵ 

The recorded proton spectrum is shown in Fig. 7 by red dots and blue error bars (corresponding to a one-sigma deviation from the average of five shots). On comparing it with the analytical thermal proton spectra (200 eV shown by red dash-dotted lines, as observed in Ref. 39 through TS), it is clear that the proton energization is nonthermal. The cutoff energy reaches about 80 keV.

Note that there is no signal recorded above the noise baseline for cases with only the magnetic field or the ambient gas, indicating that the nonthermal particle populations are indeed coming from the shock. We will pinpoint the underlying proton acceleration mechanism with dedicated PIC simulations in Sec. IV.

We use the 3D MHD code FLASH⁵⁶ to study the dynamics of the plasma plume expansion and piston formation in the ambient gas with a strong magnetic field, using the same parameters as the JLF/Titan experiment. Note that the piston is collisional and is modeled with a MHD code here; whereas the shock is collisionless and is modeled with a kinetic PIC code in Sec. IV. The simulations are initialized in a 3D geometry, using three temperatures (two for the plasma and one for the radiation) with the equation of state of Kemp and Meyer-ter-Vehn⁵⁷ and radiative transport, in the framework of ideal MHD and including the Biermann battery mechanism of magnetic field self-generation in plasmas.⁵⁸ Specifically, the laser beam is normal to a Teflon target foil and has an on-target intensity of 10¹³ W/cm², and the generated plasma plume expands in a hydrogen gas jet having a uniform density of 10¹⁸ cm⁻³. Moreover, the plume expands in a uniform external magnetic field of 20 T (aligned along the z axis, as in the experiment).

Figure 8 shows the FLASH simulation results, i.e., the electron density, electron temperature, and ion temperature from FLASH at t = 2 ns (after the laser irradiation), in two different cases (the upper row is for the case with ambient gas but no magnetic field, while the lower row is for the case with both ambient gas and a magnetic field). As FLASH cannot handle a vacuum, we do not have a FLASH simulation for the case with a magnetic field but no ambient gas). We can observe that the structures of both the hydrodynamic piston and the induced shock, which propagates inside the ambient gas, are qualitatively reproduced compared with the experiment, specifically:

The Teflon expanding piston produces a forward shock in the ambient gas (around x = 1.4 mm), as well as a reverse shock inside the Teflon piston (around x = 0.8 mm). The electron density is ∼1.6 × 10¹⁸ cm⁻³ in the forward shock in the gas and increases to ∼5 × 10¹⁹ cm⁻³ in the reverse shock.

The electron temperatures are between 60 and 70 eV in the forward and reverse shocks. Both correspond quite well to what is measured in the experiment (see the FSSR measurements in Fig. 3 and the TS measurements in Ref. 39). The ion temperature is 15 eV in the forward shock and between 80 and 180 eV inside the reverse shock.

With regard to the electron temperature, the FLASH simulation results are two times lower compared with the TS measurements in the DS region shown in Fig. 6, while for the ion temperature, the situation is worse, since this is ten times less in the forward shock compared with the TS measurements. Also note that we have not seen the foot structure ahead of the shock in the FLASH simulations.

Such discrepancies between the MHD simulations and the experiments show the difficulties in reproducing the shock condition in our case. This points to the fact that the shock evolution is dominated by kinetic effects. This is why we have resorted to using PIC simulations, the initial conditions for which are taken from the experimental measurements. Nevertheless, we can still observe that the FLASH simulations reproduce well the dynamics of the piston that induces the shock.

Since FLASH has the ability to model magnetic field generation through the Biermann battery effect, it allows us to assess the importance of this effect in the present configuration. Biermann battery generation of magnetic fields is typically important only close to the target surface (the order of 1 mm), and it is localized over the steep temperature gradients generated by the laser beam and rapidly decays once the laser beam is off (see, e.g., Refs. 59–61). As the shock is induced by the piston in the ambient gas 1 mm away from the target surface after the laser is off (∼2 ns), as shown in Fig. 8, the Biermann battery effect is negligible compared with that of the strong externally applied magnetic field.

The proton energization via the collisionless shock is modeled with the kinetic PIC code SMILEI.⁶² During the interaction between the shock front and the ambient plasma, as the width of the shock structure (of order millimeters) is much greater than the thickness of the shock (of order micrometers), we can treat this quasi–one-dimensional interaction, or a thin 1D lineout perpendicular to the curved shock (as shown in Fig. 2), via the 1D3V version of the code.

As shown in Fig. 9, the ambient plasma lies in the right half of the simulation box, while the left half is for the shocked plasma, flowing toward the right with an initial velocity v₁ = 1500 km/s. Both plasmas consist of electrons and protons, with the real mass ratio m_p/m_e = 1836. The simulation box size is L_x = 2048d_e = 11 mm, and the spatial resolution is d_x = 0.2d_e = 1.1 µm, in which d_e = c/ω_pe = 5.3 µm is the electron inertial length, with $ωpe=(ne0qe2/meϵ0)1/2=5.6×1013$ s⁻¹ being the electron plasma frequency. Here, c is the speed of light, n_e0 = 1.0 × 10¹⁸ cm⁻³ is the electron number density of the ambient plasma, and m_e, q_e, and ϵ₀ are the electron mass, elementary charge, and permittivity of free space, respectively. Note that the shock width is initialized to be equal to the ion inertial length d_i = 200 µm. The magnetic field is homogeneously applied in the z direction with B_z = 20 T (ω_ce/ω_pe = 0.06, where ω_ce = q_eB/m_e). The simulation lasts for $1.5×105ωpe−1∼2.5$ ns. Inside each cell, we put 1024 particles for each species. From the perspective of the ion Larmor motion, the simulation size is more than 10r_Li, in which r_Li = v₁/ω_ci = m_iv₁/q_eB ∼ 0.8 mm.

For the shocked plasma, the electron number density is n_e1 = 2n_e0 = 2.0 × 10¹⁸ cm⁻³, and the temperatures are T_e1 = 100 eV and T_i1 = 200 eV, all inferred from the TS characterization.³⁹ The boundary conditions for both particles and fields are open, and sufficient room is left between the boundary and the shock that the boundary conditions do not affect the physics of interest. Given the initial low temperature of the ambient plasma in the simulation (T_e0 = 50 eV), the Debye length is small compared with the grid resolution d_x, i.e., $λDe=(ϵ0kTe0/ne0qe2)1/2≈0.01de=0.05dx$⁠. However, we do run a series of simulations with different initial temperatures, showing that the energy conservation for those cases is limited around 0.05% and the physical results are almost the same. The mean free path of the presented case is λ_mfp ≈ 1800d_e, which is larger than the interaction scale, further confirming that the shock is collisionless.

We report in Fig. 10 the results of two PIC simulations, namely, with and without a magnetic field. For the case with an applied magnetic field (the left column), typical structures of a supercritical quasi-perpendicular collisionless shock can be seen.⁸ These include, for example, the overshoot in the DS region (on the left of the red dashed line), the ramps in the shock fronts (both the red dashed line and the cyan dotted line), and the foot in the upstream (US) region, as can be seen in Fig. 10(a). This foot region is formed by the reflected protons at a distance within r_L,i and modulated by the modified two-stream instability.⁶³ The proton density n_i in Fig. 10(e) indicates a compression ratio of n_i,DS/n_i,US ≈ 4, which agrees with the theoretical jump condition prediction.⁶⁴ This density profile, together with the transverse electric field E_y (not shown here), also follows the distribution of the external applied magnetic field B_z. The longitudinal electric field E_x in Fig. 10(c) peaks right at the ramps, providing the electrostatic cross-shock potential to trap and reflect the protons, as can be seen in the phase-space distribution in Fig. 10(g). Because the proton reflection is clearly due to E_x in our case, not to the DS compressed magnetic field,⁹ together with the fact that the ion Larmor radius (about 0.8 mm) is larger than the shock width (around 200 µm), the dominant particle acceleration mechanism is SSA, not SDA. Note that the cyan dotted line indicates one of the periodic shock reformations.⁸ On the contrary, for the case without a magnetic field (the right column), the drifting plasma just penetrates through the ambient gas and no shock is formed; thus, no proton energization can take place, which is in accordance with our experimental observation.

We note that the bipolar structure of the E_x electric field at the front displayed in Fig. 10(c) matches the one that was retrieved from the proton radiography data, as shown in Fig. 4. Such structures seen in the 1D PIC simulations are also verified in the complementary 2D simulations that will be detailed below. Now, we also note that, quantitatively, in Fig. 10(c), the PIC simulation gives an amplitude of the longitudinal electric field E_x ∼ 5 × 10⁸ V/m in the shock layer, which is two orders of magnitude higher than the fitting of the proton radiography in Fig. 4(c). This discrepancy regarding the amplitude of the field between the simulation and the experiment may be due to several reasons: first, the bipolar electric field structure fitted in the proton radiography has a size of 0.4 mm, while the E_x peaks in the PIC simulations are very sharp, with width smaller than 0.02 mm. With a time average of the PIC simulation over 0.2 ns, the E_x profile around the shock front reaches a size of 0.4 mm, and its value drops to 2 × 10⁷ V/m. Second, the PIC simulation represents the tip of the hemispherical expanding shock front at a single slice of the z direction, where the magnetic field is strictly perpendicular to the plasma flow and the shock is the strongest [as sketched in Fig. 4(b)]. However, the proton radiography covers the whole shock front with an integration along the z direction. It includes all other plasma flow directions in the xz plane, which are not perpendicular to the magnetic field and the corresponding shocks are weaker. This is what was already considered in the analysis of the proton radiography data, as shown in Fig. 4. Last, but not least, the amplitude of the electric field decreases with time. In our 1D PIC simulation, the results are at t = 2.7 ns (from the shock formation); while in the experimental case, the proton radiography results are at t = 5.0 ns (from the laser impact). Hence, we can expect the early-time PIC field to be higher than that derived from radiographs taken at a later time.

The particle dynamics of a high-velocity shock (as well as a comparison with the low-velocity case) and of the subsequent shock surfing proton energization are detailed in our previous paper,³⁹ while here we focus on demonstrating the robustness of the SSA mechanism that is at play in our experiment via 2D simulations, taking the nonstationarity⁶⁵ into consideration. Owing to limitations of computational resources, we reduce the 2D simulation scale to an acceptable level: the simulation box sizes are L_x = 8 mm and L_y = 0.8 mm, the simulation time t_end = 0.7 ns, and the resolution d_x = 0.4d_e.

From Fig. 11(a), we can clearly see that the transverse nonstationarity has already occurred, with 2D stripes mostly positioned at or behind the shock layer, while the trajectories of protons with kinetic energy above 30 keV show that they mainly appear at the shock front, traveling down the negative y direction. Note that the convective electric field E = −v × B is toward the positive y direction, i.e., E_y = v_xB_z, and the drifting of the protons against the convective electric field serves as a distinctive feature indicating that the dominant proton acceleration mechanism is SSA, not SDA.¹⁴ Figure 11(b) shows the proton energy spectra at 0.7 ns in both the 1D and 2D cases, which are close to each other, and there is only a 2 keV difference in the highest-energy cut, which could be caused by the numerical heating in the 2D case (with lower spatial resolution). Moreover, from the energy evolution of the protons in the x–t diagram, overlaid on the transversely averaged magnetic field map in the reference frame of the contact discontinuity (CD), it can clearly be seen that the accelerated protons are first reflected at (or picked up by) the shock front in Fig. 11(c) and then surf along the shock front while continuing to gain energy in Fig. 11(d). This is exactly the same picture as we have seen in the 1D simulations,³⁹ proving that SSA is the dominant proton acceleration mechanism at play (even in the multidimensional case).

Nevertheless, the nonstationarity of the shock might further accelerate the protons at a later time, especially after they have passed through the shock front and gyrate in the DS region. However, unfortunately, right now we do not have the computational resources to reveal that scenario. In short, our 2D simulation shows that the nonstationarity does not prevent the protons from being accelerated by SSA (reflecting and surfing), at least not at early times.

In conclusion, we have shown that laboratory experiments can be performed to generate and characterize globally mildly supercritical, quasi-perpendicular magnetized collisionless shocks. More importantly, nonthermal proton spectra have been observed for the first time, and the underlying acceleration mechanism has been pinpointed as being SSA via kinetic simulations, which can remarkably reproduce the experimental proton spectra. Such laboratory studies of proton acceleration, as well as those for electrons reviewed above, can not only further our understanding of shock formation and evolution by complementing spacecraft and remote sensing observations, but also help shed new light on solving the fundamental issue of injection for UHECR production.

Our platform can be tuned in the future to perform a systematic study of collisionless shocks with different magnetic field strengths and orientations, enabling us to capture the transition of magnetized collisionless shocks from the subcritical regime to the supercritical one, allowing us to explore the triggering of the other acceleration scenarios (i.e., SDA and DSA).

The authors would like to thank the teams of the LULI (France) and JLF Laser (USA) facilities for their expert support, as well the Dresden High Magnetic Field Laboratory at Helmholtz-Zentrum Dresden-Rossendorf for the development of the pulsed power generator used at LULI. We thank the Smilei dev-team for technical support; we thank P. Loiseau (CEA-France) for the Thomson scattering analysis code and R. Riquier (CEA, France) for the ILZ code. We also thank Ph. Savoini (Sorbonne University, France), L. Gremillet and C. Ruyer (CEA, France), M. Manuel (General Atomics, USA), and Ph. Korneev (MEPHI, Russia) for discussions. W.Y. would like to thank R. Li (SZTU, China) for discussions. This work was supported by funding from the European Research Council (ERC) under the European Unions Horizon 2020 research and innovation program (Grant Agreement No. 787539). The computational resources of this work were supported by the National Sciences and Engineering Research Council of Canada (NSERC) and Compute Canada (Job Grant No. pve-323-ac). Part of the experimental system is covered by a patent (No. 1000183285, 2013, INPI-France). The FLASH software used was developed, in part, by the DOE NNSA ASC- and the DOE Office of Science ASCR-supported Flash Center for Computational Science at the University of Chicago. We thank J. L. Dubois for providing us EOS and opacities. The research leading to these results is supported by Extreme Light Infrastructure Nuclear Physics (ELI-NP) Phase II, a project co-financed by the Romanian Government and the European Union through the European Regional Development Fund, and by the Project No. ELI-RO-2020-23 funded by IFA (Romania). JIHT RAS team members are supported by the Ministry of Science and Higher Education of the Russian Federation (State Assignment No. 075-00460-21-00). The study reported here was funded by the Russian Foundation for Basic Research, Project No. 19-32-60008.

The authors declare no conflicts of interest.

All data needed to evaluate the conclusions in the paper are present in the paper. Experimental data and simulations are archived on the servers at LULI and LERMA laboratories, respectively, and are available from the corresponding author upon reasonable request.