Kinetic study of shock formation and particle acceleration in laser-driven quasi-parallel magnetized collisionless shocks

Magnetized collisionless shocks are ubiquitous in the universe and are presumed to be the source of some of the highest energy cosmic rays. Fermi acceleration was first proposed to account for the origin of the inverse power law spectral distribution of these cosmic rays.¹ In Fermi's theory, interstellar particles interact with wandering magnetic fields and gain energy from experiencing more head-on collisions than overtaking collisions. The stochastic mechanism is often known as “magnetic mirrors,” or second-order Fermi acceleration. Fermi's theory was further developed and referred to as diffusive shock acceleration (DSA) in the late 1970s,^2–5 which is associated with particle acceleration in strong shock waves, such as supernova explosions.⁶ DSA requires particles to have large enough gyro-radii to pass through the shock without strong deflection in order to be accelerated. However, how particles acquire the initial energy required to do this remains an open question (the “injection” problem). It is generally believed that ions can fulfill the condition for DSA more easily than electrons because of their larger Larmor radii compared to the shock thickness.⁷ Although DSA theory has been used to account for observations, direct evidence remains scarce.

In magnetized shocks, the angle $θn$ between the shock normal and the external magnetic field is a fundamental feature that dictates key shock characteristics. The angle $θn$ is used to categorize magnetized shocks into quasi-perpendicular (⁠ $θn>45∘$⁠) and quasi-parallel (⁠ $θn<45∘$⁠) geometries. Unlike quasi-perpendicular shocks, which can be formed on timescales less than an ion gyro-period,^8–11 quasi-parallel shocks rely on ion–ion instabilities that grow slowly and therefore are much more difficult to form. Consequently, while quasi-perpendicular shocks have been successfully created on several laser facilities,^12–15 quasi-parallel shocks have yet to be created experimentally, due to the high laser energy required to enable a sufficiently large plasma, among other constraints.¹⁰ Quasi-parallel shocks are believed to be more efficient accelerators of particles than quasi-perpendicular shocks.¹⁶ 

As a complement to astronomical observations, laboratory shock experiments driven by lasers can provide a repeatable, and systematic way of studying collisionless shock physics,¹⁷ and advanced diagnostic systems provide additional measurements and characterizations that are not possible in space. Laboratory experiments can access early stages of shock formation, while astrophysical shocks are typically already fully formed. The laboratory formation of quasi-parallel collisionless shocks requires the combination of large magnetized system size and strongly driven flows so that the shocks have sufficient space and time to develop and accelerate particles while still remaining collisionless. For example, the density of the laser-produced plasma must be low enough that Coulomb collisions remain negligible, but high enough to fit hundreds of ion inertial lengths that allow the shock formation in the experimental apparatus; the applied magnetic field must be strong enough to magnetize the plasma, but weak enough to satisfy the Mach number requirement.¹⁰ The experimental and physical constraints rule out most of the existing laser facilities. Our experiments are currently underway to achieve these conditions by utilizing the National Ignition Facility (NIF), the world's largest laser facility.^18–22 

Earlier NIF experiments have observed unmagnetized collisionless shocks and associated energized electrons,²³ and previous quasi-perpendicular shocks have observed proton energization.¹⁵ Some of the beam instabilities that contribute to parallel shock formation have been observed at the LAPD.^24,25 However, there is as yet no direct evidence of shock-accelerated particles from quasi-parallel laser-driven experiments available, which are more relevant to astrophysical observations. Numerical simulations can play a key role in facilitating the design of laboratory experiments on quasi-parallel shocks. Previous hybrid simulations^16,26,27 revealed much information on shock formation and particle acceleration, including acceleration efficiency, magnetic field amplification, and particle diffusion. However, hybrid studies utilize the kinetic–proton and fluid–electron model intrinsically and assume a Maxwellian thermodynamically equilibrium for electrons. Recent works^28–30 consider energy exchange between the ions and electrons in collisional plasma. In contrast, full particle-in-cell (PIC) simulations, which treat both ions and electrons kinetically, allow a self-consistent collisionless kinetic ion–electron energy exchange, which will be shown to be important under experiment-relevant conditions (see Sec. IV B).

In this paper, PIC simulations are used to study the formation of quasi-parallel collisionless shock and particle acceleration, relevant to NIF experimental parameters. With PIC simulations, not only are shock formation and ion acceleration studied, but also electron spectra are obtainable. Upstream waves generated by reflected high-energy particles are identified as Alfvén waves and magnetosonic waves. The downstream ion–electron energy partition is found to be on the order of unity, which is far from independent heating, implying a significant ion–electron energy exchange in the downstream.

This paper is organized as follows: Sec. II presents the simulation setup. The formation of quasi-parallel shocks within the parameter range achievable on the NIF at different oblique angles $θn$⁠, shock structures, and ion acceleration are discussed in Sec. III. This section also gives a comparison of the acceleration capabilities of quasi-parallel and quasi-perpendicular shocks. Section IV briefly discusses upstream wave excitation, and ion–electron energy exchange observed in the downstream. Conclusions are given in Sec. V.

The simulations are designed to model experiments currently underway on the NIF. In those experiments, a laser-driven piston plasma expands supersonically into an ambient magnetized plasma. The piston acts to sweep out magnetic flux and ambient particles, accelerating them to super-magnetosonic speeds and driving a collisionless shock.^13,14,31 To achieve quasi-parallel shocks, the piston is designed to expand at an angle of $θn=30∘$ relative to an externally applied background magnetic field. The size of the ambient plasma extends for several hundred ion inertial lengths, allowing the shock to propagate for several tens of ion gyro-times. The schematic of the experimental setup is shown in Fig. 1(a). Additional details on the experimental setup will be reported in a future paper.

The simulations are performed using OSIRIS 4,³² a fully kinetic, relativistic, parallel PIC code. Our simulations, shown schematically in Fig. 1(b), utilize a reflective wall (ideal piston) boundary condition. Instead of simulating a piston plasma, a fixed reflective wall is implemented on the right boundary, from which a clear downstream, free of the influence of the piston, is obtained. This setup simplifies piston–ambient coupling^11,13,14,31 via ambient–wall collision and focuses on the shock itself. Magnetized plasma drifting against the wall (⁠ $+x$⁠) with a bulk velocity of $Vd$ is reflected off the wall and forms a shock that propagates to the left (⁠ $−x$⁠) in the simulation. This ideal wall is a good analog for the solid-density piston in the experiments.⁸ The simulations run in a piston/wall frame of reference, in which most of the figures in this paper are shown, except where the laboratory frame is explicitly stated in the captions. The velocity transformation between the laboratory and simulation frame simply follows $vlab.=vsim.−Vdex$⁠, where the superscripts denote the frames of reference, and $ex$ represents the unit vector in the x-direction.

Unlike perpendicular shock formation, where the ion–electron mass ratio is critical,⁸ the formation and evolution of quasi-parallel shocks are dominated by ion-scale physics.^10,16,33,34 A reduced ion mass of $mi/me=50$ is thus used in most of our two-dimensional (2D) PIC simulations. Some key initial upstream parameters are summarized in Table I. To better describe the shock, the subscripts 1 and 2 are hereinafter used to denote the upstream and downstream, respectively. Our simulations use a rectangular simulation domain in the $x−y$ plane, with a minimum size of 883  $c/ωpi1$ in the x-direction and 21  $c/ωpi1$ in the y-direction. The shock propagates in the negative x-direction, and the external B field is uniform in the $x−y$ plane, at shock angle $θn$ to the shock normal. In the y-direction, a periodic boundary condition is applied for particles and fields. The grid size is 0.5  $c/ωpe1$⁠, and the timestep is 0.35  $ωpe1−1$ to ensure convergence. We used 12 × 12 macro-particles per cell in the simulations to maintain a low numerical collision effect. The total energy is conserved within 0.1% in the simulations. The expected mean free path³⁵ between the piston ions and the ambient ions $λii≈2.5$ cm is comparable to the system size in the experiments; therefore, the physics involved is fundamentally collisionless. In the simulations, the collision module in OSIRIS is turned off.

Quasi-perpendicular shocks in strongly magnetized cases relevant to kilojoule, terawatt laser experiments like OMEGA EP,^36–38 have been found to be mediated by a modified two-stream instability (MTSI) between incoming and reflected ions, so can be formed within an ion gyro-period.⁸ The MTSI relies on magnetically tied electrons and unmagnetized (or weakly magnetized) ions, similar to the lower hybrid mode, and has a growth rate much larger than the ion gyro-frequency,^8,39,40 thus resulting in a shock formation time much less than an ion gyro-period $Tci(=2π/ωci)$⁠. In the quasi-parallel regime, however, the absence of an effective MTSI causes a much longer shock formation time, on the order of ion gyro-periods.

The dependence of shock formation on the shock angle $θn$ in a quasi-parallel geometry has been studied by 2D PIC simulations for $θn=0∘$⁠, 15°, and 30°. Although shocks with smaller $θn$ take longer to form, we found that supercritical shocks where reflected ions are present and important^34,41 can all be formed within a few ion gyro-periods, consistent with previous studies.⁴² The 2D density maps of three simulations are presented in Fig. 2, with y-averaged density profiles (solid lines) plotted as well. The three rows of Fig. 2 correspond to $t=$ 3.15, 6.30, 12.59  $ωci1−1$⁠, respectively. An effective density jump^43, $r(≡n2/n1)≳4$⁠, which indicates the formation of the shocks, is observed in Figs. 2(c), 2(e), and 2(g) at $θn=0∘$⁠, 15°, and 30°, respectively, despite different shock profiles. The 3D adiabatic limit for the density jump is 4.^43,44 The observed higher density jump is an overshoot in this transition region. It relaxes toward 4 in the downstream region. The density compression builds up faster at $θn=30∘$⁠, within one ion gyro-period, while at smaller $θn$ a similar compression takes longer to form. The shock foot, which consists of ions reflected from the shock front, is larger for smaller oblique angles because the ion trajectories are impeded by the perpendicular component of the magnetic field, namely, $By$⁠. Those reflected ions then gyrate upstream along the B field lines, excite waves, and deposit their energy to the ambient plasma, which is called the “diffusion” process. Ultimately, the upstream region reaches a finite length^16,26 and sets an energy threshold, above which ions can escape from the shock. Therefore, as can be clearly distinguished in Figs. 2(c), 2(f), and 2(i), the wave pattern in the foot shows a close correlation with the B field direction.

Therefore, based on the above simulation results, we will use $θn=30∘$ as our main setup following the experiments. In this case, the quasi-parallel shock formation time can be as short as $∼3ωci1−1$ (∼1 ns in the experiments). This setup gives faster shock evolution than $θn=0∘$⁠, while preserving features of particle acceleration similar to $θn=0∘$⁠. Most importantly, as will be shown in Sec. III B, the energy spectra of energized ions remain similar between the $θn=30∘$ case and smaller angles within the time accessible in the experiments.

In the theoretical picture of DSA, particles gain energy by interacting with the shock front. However, a one-time interaction can only increase the energy of a particle by a moderate amount. Particles interacting many times with the shock can achieve very high energies.⁴⁵ Repetitive acceleration by crossing the shock front, one of the key features of DSA theories, is seen in our simulations. However, the experimentally accessible time is far from the time required for DSA to fully develop a population of highly energized ions, which is on the order of hundreds or thousands of ion gyro-times. Therefore, even on the NIF, only the early stages of quasi-parallel shock formation and the DSA process are achievable. We therefore limit the total time in our PIC simulations to 40  $ωci1−1$ to determine experimentally relevant characteristics of accelerated particles. A longer simulation with a sustained shock potentially could yield similar power-law energetic particle distributions seen by the hybrid simulations.^16,26,27

The trajectories of some representative ions, in the $θn=30∘$ simulation, are plotted in $t−x$ [Fig. 3(a)] and $x−y$ [Fig. 3(b)] plots, in which their kinetic energies are indicated by color. Note that although the trajectories are given in the simulation frame for consistency, the particle energies are calculated in the laboratory frame throughout this paper. In Fig. 3(a), the y-averaged density profile is shown in the background in grayscale, and the shock front is plotted by the white dashed line. Tracked ions are repetitively reflected by the shock and gain energy by the Fermi process.^1,3 The maximum kinetic energy the ions can achieve within 40  $ωci1−1$⁠, comparable to the experiments, is ∼350 keV, which is ∼35× the shock energy (defined as $Esh=12miVsh2=9.6$ keV, in which $Vsh=43Vd$ is used to approximate the shock velocity for a strong shock⁴³). Figure 3(b) shows the 2D trajectories in the $x−y$ plane where the periodic boundary conditions in the y-direction are post-processed to show a continuous 2D space. Accelerated ions gyrate along the B field lines in the upstream, which would make them convenient to collect experimentally [refer to Fig. 1(a)].

Figure 4 shows the history of particle energization, in which the energy spectra at different timesteps are plotted in different colors. Figures 4(a)–4(c) are the ion spectra for $θn=0∘$⁠, 30°, and 60°, respectively; the electron spectra are shown in Figs. 4(d)–4(f). Nonthermal particle populations are created by both the quasi-parallel and quasi-perpendicular shocks. However, spectra that extend further indicate that even at this early stage in their formation, quasi-parallel shocks [Figs. 4(a) and 4(b)] are capable of energizing more particles to higher energies than quasi-perpendicular shocks [Fig. 4(c)]. $E=10 Esh$ is marked by the vertical lines in Fig. 4. In Fig. 4(c), the ions are accelerated to higher energies faster by the quasi-perpendicular shock but achieve a nearly stable spectrum with a steep energy cutoff at $∼10 Esh$⁠. In contrast, as shown in Figs. 4(a) and 4(b), the energy spectra produced by quasi-parallel shocks continue to grow further to higher energies. Only quasi-parallel shocks can accelerate a non-negligible population of ions to $≳10 Esh$⁠,^15,16 meaning that any significant population of particles measured above this range in the experiment can only be produced by the quasi-parallel shock. We then define an ion acceleration efficiency as $ηi=NE>10 Esh/NDS$⁠, where $NE>10Esh$ is the number of ions whose energy is above $10 Esh$⁠, and $NDS$ is the number of ions in the downstream. For $θn=0∘$⁠, 30°, and 60° shocks, the efficiencies of $ηi≈0.14%$⁠, $0.06%$⁠, and $0.0004%$⁠, respectively, are obtained from the distributions at $t=40ωci1−1$⁠. Therefore, the parallel shocks are expected to be more efficient in particle acceleration. While ions can be accelerated to higher energies by quasi-parallel shocks, the experimentally relevant timescales explored here are still too short to obtain a power-law spectrum³ that is usually seen in astrophysical observations of energized particles. Figures 4(d)–4(f) indicate significant electron energization by the shocks, which will be discussed in Sec. IV B.

Nevertheless, the simulations show that this experimental achievable early-stage quasi-parallel shock has already (1) formed a density compression $r∼4$ that Rankine–Hugoniot conditions⁴³ predict for a strong shock (Fig. 2); (2) created an environment that allows ions to repetitively cross the shock front and gain energy (Fig. 3); and (3) energized a number of ions that exceed $10 Esh$ (Fig. 4).

In quasi-parallel shocks, high-energy ions that are accelerated by the first-order Fermi process gyrate along the B field in the upstream. They also interact with ambient plasma and excite waves^5,46 as they escape upstream and travel in the ambient plasma. Figure 5(a) shows the $Bz$ field profile at $t=12.59ωci1−1$ in the $θn=0∘$ simulations. Unlike downstream turbulence, apparent wave patterns are seen in the far upstream region (weakly perturbed, $x<790cωpi1−1$⁠). In Fig. 5(b), the Fourier spectrum of the upstream $Bz$ field shows an X-shaped mode, whose arms are at an oblique angle of 65° to the $kx$ axis. In the weakly perturbed upstream region, we can use the upstream initial parameters in Table I to estimate the magnetohydrodynamic (MHD) wave velocities, $vms1≈2.17vA1$⁠, where $vms1=vA12+cs12$ is the upstream magnetosonic wave velocity, $vA1=B1/4πmini1$ is the upstream Alfvén velocity, $cs1=γkB(Te1+Ti1)/mi$ is the upstream sound velocity, $γ=5/3$ is the polytropic index, and $kB$ is Boltzmann constant. As the upstream magnetosonic waves and Alfvén waves propagate at their own characteristic velocities, the wavefront of the combined wave would show a propagation at an oblique angle of $arctan(vms1/vA1)≈65∘$⁠, which agrees with the X-shaped mode in Fig. 5(b). Similarly, in Fig. 5(c), the Fourier spectrum of the upstream waves in the $θn=30∘$ simulation also shows an X-shaped mode but with a tilt angle of 30°. Therefore, the far upstream waves are composed of two MHD waves, magnetosonic waves perpendicular to the B field, and Alfvén waves parallel to the B field. These MHD waves are excited by the accelerated energetic ions that travel through the upstream ambient plasma and deposit their energies. The ions are diffused upstream, slowed down, and caught up by the shock, then enter downstream. Since the magnetic fields dominate the excited waves upstream, the ions are likely reflected by the magnetic fields. However, this needs further study. The downstream turbulence acts as another diffusive mechanism that reflects the ions back to the upstream. This process repetitively reflects and accelerates the ions until they gain enough energy and finally escape.^3,5 Observing these waves in the simulation shows that the formation of the upstream, a critical stage in the formation of a quasi-parallel shock, is expected to be achievable in the NIF experiment.

The energy partition result implies a significant electron heating due to ion–electron energy exchange, which also agrees with planetary and astrophysical observations⁴⁷ that indicate $Ti/Te≈0.5−30$ under different interstellar conditions. As a comparison, adiabatic electron heating, which is assumed in hybrid simulations, would give a downstream electron temperature $Te2=Te1·rγ−1≈250$ eV. Independent thermalization of each species would predict a temperature increase $ΔTi/ΔTe≈mi/me=50$⁠, i.e., $ΔTe≈70$ eV. Both estimates are much smaller than the simulation result.

We note that this conclusion may only be valid for the low Mach number, quasi-parallel regime, in contrast to astrophysical observations where the Alfvénic Mach number is of the order of 100.²³ The isotropization process is likely the result of a complicated mixture of micro-instabilities, until turbulence is sufficiently developed. The relevant instabilities will be discussed in future publications.

This paper focused on ion acceleration, which is expected to be a key signature of the planned NIF experiments. However, non-thermal electrons were also observed in Figs. 4 and 6, indicating electron acceleration. In addition, collisionless ion–electron energy transfer could also play a role. Possible mechanisms are not currently understood and will be a future research topic.

Quasi-parallel shock formation relevant to laser-driven experiments on the NIF has been studied with 2D PIC simulations. It is found that an early stage of quasi-parallel shock formation can be achieved on the NIF, and that energetic ions accelerated by the Fermi process are observable. Even so, the limited experimental extent in time and distance means that a power law energy spectrum does not have time to develop. However, one of the key features of DSA, ion repetitive crossing of the shock front by scattering from waves and turbulence is seen in the simulations. We expect that the shock formation and structure are essentially captured by these 2D simulations, based on the observed similarity between our 1D and 2D simulations. This is further supported by the observed similarity between 2D and 3D hybrid simulations.¹⁶ However, whether particle acceleration or ion–electron coupling may change in 3D remains to be studied in the future.

More importantly, the simulations confirm that quasi-parallel shocks are capable of energizing more particles and creating energy spectra that extend further than quasi-perpendicular shocks. Under experimental conditions, our simulations indicate that a population of ions with energies above $10Esh$ can only be effectively produced by quasi-parallel shocks. Upstream waves excited by these accelerated ions are identified as parallel Alfvén waves and perpendicular magnetosonic waves. In the downstream, ions and electrons reached an energy partition $Ti2/Te2≈1.3$⁠. These findings, including density compression, high energy ions, turbulent structures produced by the shock, and the temperatures of ions and electrons, provide experimental observables that will be tested in the experiments on the NIF. Most importantly, energetic particles accelerated by quasi-parallel shocks have never been experimentally observed before. This NIF platform aims to validate the experimental capability of quasi-parallel shock studies.

The authors would like to thank Dr. Damiano Caprioli for the useful discussions. This material is based upon work supported by the Department of Energy (National Nuclear Security Administration) University of Rochester “National Inertial Confinement Fusion Program” under Award No. DE-NA0004144, the Department of Energy under Award Nos. DE-SC0020431 and DE-SC0024566, and the resources of the National Energy Research Scientific Computing Center (NERSC), a U.S. Department of Energy Office of Science User Facility located at Lawrence Berkeley National Laboratory. The authors thank the UCLA-IST OSIRIS consortium for the use of OSIRIS. This report was prepared as an account of work sponsored by an agency of the U.S. Government. Neither the U.S. Government nor any agency thereof, nor any of their employees, makes any warranty, express or implied, or assumes any legal liability or responsibility for the accuracy, completeness, or usefulness of any information, apparatus, product, or process disclosed, or represents that its use would not infringe privately owned rights. Reference herein to any specific commercial product, process, or service by trade name, trademark, manufacturer, or otherwise does not necessarily constitute or imply its endorsement, recommendation, or favoring by the U.S. Government or any agency thereof. The views and opinions of authors expressed herein do not necessarily state or reflect those of the U.S. Government or any agency thereof.

The authors have no conflicts to disclose.

Yu Zhang: Conceptualization (equal); Data curation (lead); Formal analysis (lead); Investigation (lead); Visualization (lead); Writing – original draft (lead); Writing – review & editing (equal). Peter V. Heuer: Formal analysis (equal); Supervision (equal); Writing – review & editing (equal). Jonathan R. Davies: Conceptualization (equal); Formal analysis (equal); Funding acquisition (lead); Supervision (equal); Writing – review & editing (equal). Derek B. Schaeffer: Conceptualization (equal); Writing – review & editing (equal). Han Wen: Formal analysis (supporting); Writing – review & editing (supporting). Fernando García-Rubio: Formal analysis (supporting); Writing – review & editing (supporting). Chuang Ren: Conceptualization (equal); Formal analysis (equal); Resources (equal); Software (equal); Supervision (equal); Writing – review & editing (equal).

The data that support the findings of this study are available from the corresponding author upon reasonable request.