PARTICLE ACCELERATION IN RELATIVISTIC MAGNETIZED COLLISIONLESS PAIR SHOCKS: DEPENDENCE OF SHOCK ACCELERATION ON MAGNETIC OBLIQUITY

Unmagnetized relativistic shocks are mediated by Weibel-like filamentation instabilities (Weibel 1959; Medvedev & Loeb 1999; Gruzinov & Waxman 1999). Incoming particles are scattered by the Weibel-generated fields and the incoming flow is isotropized, decelerated and compressed (e.g., Spitkovsky 2005). The same mechanism is likely to mediate moderately magnetized (σ ≲ 1) parallel shocks. In fact, we find that the width of the transition layer in σ = 0.1 parallel shocks, as evaluated from the density profile (Figure 3(b) at x_sh ≈ 730 c/ω_p), is of the order of a few tens of skin depths, as observed in unmagnetized shocks.³ Also, the magnetic energy density peak at the shock (reaching ∼30% of the upstream kinetic energy density, see Figure 3(d)) is a characteristic signature of magnetic field amplification by a Weibel-like process. The field then decays to the plateau predicted by MHD calculations on a length scale ∼50 c/ω_p.

Oblique shocks (θ ≠ 0°) may be triggered by another process, namely, magnetic reflection of the incoming flow from the shock-compressed magnetic field (Alsop & Arons 1988), as reported for perpendicular shocks by previous 1D (e.g., Langdon et al. 1988; Gallant et al. 1992) and 3D (Spitkovsky 2005) PIC simulations. The synchrotron maser instability, excited by coherent gyration of the incoming particles in the shock magnetic field, mediates the thermalization of the upstream plasma (Hoshino & Arons 1991). Indeed, the peak in the transversely averaged density profile observed for θ = 45° at the shock location (see Figure 9(b) at x_sh ≈ 920 c/ω_p) is due to Larmor gyration of the incoming particles in the compressed shock fields, which decelerates the upstream flow and therefore increases its density. The density rise is as sharp as a few Larmor radii in the compressed field (i.e., a few plasma skin depths), much thinner than for parallel shocks (compare Figure 9(b) with Figure 3(b)). In just a few Larmor scales downstream from the shock, the density profile flattens to a plateau in reasonable agreement with MHD jump conditions. Shocks with larger obliquities (θ ≳ 45°) are also triggered by magnetic reflection.

Both Weibel instabilities and magnetic reflection play a role in mediating shocks with intermediate obliquities (0° < θ ≲ 45°). Filamentation instabilities dominate for lower obliquities, such as θ = 15°, whose fluid structure resembles θ = 0° in the magnetic overshoot at the shock (Figure 5(d) at x_sh ≈ 720 c/ω_p). Instead, such Weibel-induced magnetic energy excess is not seen for θ = 30° (Figure 7(d) at x_sh ≈ 770 c/ω_p) or larger obliquities; rather, investigation of early times (ω_pt = 225) for θ = 30° shows that this shock started out with the prominent density peak characteristic of highly oblique shocks, which suggests that this shock is mainly triggered by magnetic reflection. For θ = 30°, the sharp density peak observed at early times was later smoothed by the increased pressure in shock-accelerated particles.

For all oblique magnetic configurations (0° < θ < 90°), a fluid structure with two shocks (a strong "fast" shock and a weak "slow" shock) is seen in our PIC simulations. The presence of two shocks is predicted for the parameters (γ₀ = 15 and σ = 0.1) that we adopt (Landau et al. 1984) and has been confirmed by MHD simulations of the collision of a magnetized shell against a wall (S. Komissarov 2007, private communication). The particle number density rises across both shocks; however, while the fast shock involves an amplification of the magnetic field, across the slow shock the magnetic energy decreases. The slow shock propagates into the downstream medium of the fast shock. Its transition layer in our PIC simulations is not as narrow as MHD would predict, due to insufficient dissipation along the magnetic field. Rather, the slow-shock transition region stretches self-similarly with time; the jump in fluid quantities between the wall and the fast-shock downstream stays roughly constant in time. We observe that the slow shock is slower and weaker as θ approaches either 0° or 90°, and relatively stronger and faster for intermediate obliquities (θ ≈ 30°). For oblique magnetic configurations, we shall refer to the fast shock simply as "the shock" and to its downstream region as "the downstream medium."

For the main (i.e., fast) shock, the velocity and jump conditions extracted from our PIC simulations agree within a few percent with the solution of relativistic MHD equations (Appl & Camenzind 1988) supplemented by the 3D Synge (1957) equation of state, provided that MHD results are Lorentz-transformed from the standard shock frame into the simulation frame.⁴ In Table 1, we compare the results from PIC simulations and MHD calculations for the shock velocity β_sh and the jump across the shock in particle number density (n_d/n_u), magnetic energy density (B²_d/B²_u), transverse magnetic field (B_z,d/B_z,u), and transverse electric field (E_y,d/E_y,u). We remark that, for oblique shocks (0° < θ < 90°), the background motional electric field (last column in Table 1) does not drop to zero after the shock since, as discussed above, the downstream flow preserves a small transverse velocity β_z,d < 0 which gives a motional field E_y,d = −β_z,dB_x,d > 0.

As shown in Table 1, for θ ≲ 30° our kinetic calculations are in remarkable agreement with 3D MHD results, suggesting that for low magnetic obliquities the downstream plasma in our 2D simulation box has a 3D effective adiabatic index, as anticipated above. For θ ≳ 45°, the adiabatic index of the downstream fluid in our simulations tends to that of a 2D gas, which explains the discrepancy between PIC simulations and 3D MHD calculations for θ = 45° in Table 1.⁵

3.1.1. Comparison with Previous 1D Simulations

The internal structure at time ω_pt = 2250 of shocks with magnetic obliquities θ = 0° (Figures 3 and 4), 15° (Figures 5 and 6), and 30° (Figures 7 and 8) shows some features common to all subluminal magnetic configurations. Subluminal shocks are characterized by the presence of a diffuse stream of high-energy shock-accelerated particles that propagate ahead of the shock and trigger the growth of oblique waves in the upstream medium. In turn, as we explain in Section 5, such waves are of fundamental importance for particle acceleration.

The longitudinal phase space plots of positrons in panels (l)–(n) of Figures 3, 5, and 7 show the injected plasma as a cold dense beam propagating with γβ_x ≈ −15 and broadening while it approaches the shock, which is located at x_sh ≈ 730 c/ω_p for θ = 0°, at x_sh ≈ 720 c/ω_p for θ = 15°, and at x_sh ≈ 770 c/ω_p for θ = 30°. In the shock transition layer, the incoming plasma is isotropized and thermalized. A population of high-energy particles moving ahead of the shock appears in the upstream region as a hot diffuse beam with γβ_x > 0 (panel (l) in Figures 3, 5, and 7); as we show in Section 5, these particles have been accelerated at the shock. Additional insight is provided by the positron γβ_x − γβ_z momentum space and electron (blue) and positron (red) spectra in Figures 4, 6, and 8, extracted from three different locations across the flow as marked by the arrows at the bottom of Figures 3, 5, and 7, respectively. Upstream slabs (panels (c) and (f) in Figures 4, 6, and 8) show two components: the cold incoming plasma appears in particle spectra (panel (f) in Figures 4, 6, and 8) as a relatively narrow bump peaking at low energies, around γ ≈ 15, which broadens as the flow approaches the shock; the returning diffuse beam of shock-accelerated particles populates the high-energy bump, which is more and more depleted (especially at intermediate energies) as we move farther upstream. Far downstream slabs (panels (a) and (d) in Figures 4, 6, and 8) show roughly isotropic particle distributions with thermal spectra supplemented by high-energy tails. For θ = 30°, the downstream population of energetic particles significantly varies along x, with regions closer to the shock showing more particles of higher energies than slices farther downstream (Figure 7(l)–(n), and compare panels (a) and (d) with panels (b) and (e) in Figure 8). Since fluid elements which are farther downstream passed through the shock at earlier times, this suggests that the shock acceleration efficiency increased over time, as we discuss in Section 4.2.

For oblique magnetic configurations (Figure 5 for θ = 15° and Figure 7 for θ = 30°), the phase space plots x − γβ_y (panel (m)) and x − γβ_z (panel (n)) are not symmetric around γβ_i = 0 (i = y,  z), when focusing on the high-energy particles located in the vicinity of the shock (see also panels (b) and (c) in Figures 6 and 8). The high-energy positrons with γβ_z > 0 in the shock and upstream regions are shock-accelerated particles whose gyro-center is heading upstream along the oblique magnetic field (we remind that B_x, B_z > 0).⁶ With increasing magnetic obliquity from θ = 15° to 30°, the asymmetry around γβ_z = 0 becomes more evident and the returning high-energy positrons show increasingly larger γβ_z than γβ_x (compare panels (l) and (n) between Figures 5 and 7). The electron phase space x − γβ_z (not plotted) shows the same asymmetry.

The excess of high-energy positrons with γβ_y < 0 observed in the vicinity of the shock (panel (m) of Figures 5 and 7) is a natural consequence of the ∇B drift experienced by particles accelerated at the shock. For positively charged particles, the drift is directed along $-\hat{\mathbf {y}}$, which explains the sign of the asymmetry in Figures 5(m) and 7(m). Since the magnetic field jump at the shock is larger for increasing obliquities (Table 1), and the drift velocity is then higher, we expect the asymmetry to be more significant for θ = 30° than for 15°, as observed. Moreover, since electrons and positrons move in opposite directions under ∇B drift, the electron x − γβ_y phase space (not shown) exhibits a reversed asymmetry, with an excess of particles having γβ_y > 0.

The returning beam of high-energy particles triggers in the upstream medium the generation of waves whose wavevector is oblique to the background magnetic field. Such waves have an electromagnetic component (B_z in panel (e) and E_y in panel (g) of Figures 3, 5, and 7) and a smaller electrostatic component (E_x in panel (i) of Figures 3, 5, and 7). The electrostatic component also appears as upstream oscillations in the transversely averaged 〈E_x〉 (panel (j) of Figures 3, 5, and 7) and in the electric pseudo-potential −∫〈E_x〉 dx (panel (k) of Figures 3, 5, and 7). The upstream waves grow in amplitude while propagating toward the shock, where the regular oblique pattern observed far upstream (with characteristic wavelength ∼20 c/ω_p at time ω_pt = 2250) changes to a more complicated structure (e.g., compare x ≳ 1000 c/ω_p with 750 c/ω_p ≲ x ≲ 900 c/ω_p for θ = 15° in Figure 5(e), (g), and (i)). The ratio between the magnetic energy of the oblique waves (mostly contributed by their B_z) and of the background upstream field is of the order of a few tens of percent, roughly constant with magnetic obliquity.

By comparing panels (e), (g), and (i) with panels (l)–(n) of Figure 7 for θ = 30°, the growth of the upstream oblique modes appears to be spatially coincident (at x ≈ 1000 c/ω_p) with the leading edge of the population of returning particles. The same holds true for smaller obliquities, but for θ = 0° and 15° the returning particles that are farthest upstream are outside the x-range plotted in Figures 3 and 5. This suggests that interplay between the unperturbed incoming plasma and the returning high-energy beam may be responsible for triggering the oblique waves.

The growth of the upstream waves is neither an unphysical consequence of our jumping injector nor an effect of the numerical heating associated with the propagation of a cold beam through the grid. In fact, if we replace the reflecting wall at x_wall = 0 with an open boundary, everything else being equal, we observe that no upstream waves are produced, nor is there a shock. In order to test the role played by the returning high-energy particles in generating the upstream waves, we have performed a suite of simulations in which we artificially suppress particle acceleration. In these runs, we introduced cooling, so that particles with γ>γ_cool lose a random fraction of their nonthermal energy (a similar technique was used in Keshet et al. 2009). By varying γ_cool (we tried γ_cool = 50,  80,  140 for θ = 0°), we conclude that particle cooling suppresses the growth of the upstream waves, and by decreasing γ_cool (i.e., by forcing more particles to cool) both the wave amplitude and the distance from the shock where they are triggered decrease.

In summary, the oblique waves seen in the upstream medium are generated by the returning beam of high-energy particles; they may result from a mixed-mode streaming instability which couples the electrostatic two-stream instability with the electromagnetic filamentation (Weibel) instability, as envisioned by Bret & Dieckmann (2008). A complete analysis of the wave dispersion relation and generation mechanism is deferred to a future paper.

Finally, we point out that the electrostatic component of the upstream waves induces oscillations in the γβ_x-velocity of the incoming plasma, as seen in the positron phase space of Figure 3(l). The electron x − γβ_x phase space (not shown) presents oscillations of opposite sign. The difference at low energies between electron and positron spectra in the upstream medium (panel (f) in Figures 4, 6, and 8) can be correlated with the sign of the electric pseudo-potential at the same longitudinal location (panel (k) in Figures 3, 5, and 7, respectively), since a net electric potential accelerates incoming charges of one sign and decelerate charges of the opposite sign. We further comment on the effects on particle spectra of the electrostatic component of upstream oblique waves in Section 4.2.

Figures 9 and 10 show the shock structure, phase space plots and particle spectra for a representative superluminal shock (θ = 45°) at time ω_pt = 2250. No significant particle acceleration occurs in superluminal shocks, and the upstream plasma is not appreciably perturbed by the presence of the shock.

In the phase space plot of Figure 9(l), the incoming plasma is seen as a cold beam with γβ_x ≈ −15, which is slightly heated while approaching the shock and then sharply thermalized in the shock transition layer at x_sh ≈ 920 c/ω_p. As discussed by Hoshino & Arons (1991) and anticipated in Section 3.1, the incoming plasma thermalizes due to the absorption of cyclotron waves produced by coherent particle gyration in the compressed shock fields. Thermalization proceeds with time in the downstream medium: panels (l) and (m) of Figure 9 show that the plasma farther downstream, which had more time to thermalize after passing through the shock, has a broader x − γβ_x and x − γβ_y distribution than the fluid closer to the shock.

Downstream particles preferentially move in a plane perpendicular to the downstream oblique field; indeed, the oblique boxy shape of their downstream γβ_x − γβ_z momentum space (Figures 10(a), 10(b) for two different downstream slabs) is suggestive of particle gyro-motion around the oblique field. As discussed at the beginning of Section 3, for magnetic obliquities θ ≳ 45° it is harder for downstream particles to isotropize along the magnetic field rather than in the plane orthogonal to it. Since for high obliquities such plane approaches the xy simulation plane, this should result in the downstream thermal width of the x − γβ_z phase space being smaller than the x − γβ_x or x − γβ_y distributions, as observed (compare panel (n) with panels (l) and (m) in Figure 9).

The net fluid velocity β_z,d < 0 behind the shock in Figure 9(n) is responsible for the residual motional electric field E_y,d > 0 observed in the simulation frame downstream from the shock (hard to see in Figure 9(h), but see Table 1), a feature common to all oblique magnetic configurations, as discussed in Section 3. This residual fluid velocity vanishes farther downstream, in regions already overtaken by the slow shock. The slow shock appears as a smooth increase in density (Figure 9(b)) and decrease in magnetic energy (Figure 9(d)) and B_z (Figure 9(f)) behind the main shock toward the wall. Similar properties for the slow shock are observed for θ = 30° in Figure 7, whereas the slow shock for θ = 15° is outside the x-range shown in Figure 5.

Most notably, no evidence for particle acceleration is present for θ = 45°. No high-energy particles are observed to propagate upstream from the shock (Figure 9(l)–(n)), and the particle energy spectrum in upstream slices is just the thermal distribution of the injected plasma Lorentz-boosted along $-\bf {\hat{x}}$ (Figure 10(f)). The range of downstream four-velocities is much smaller than for subluminal shocks (compare panels (l)–(n) in Figure 9 with the same panels in Figures 3, 5, and 7), and the downstream energy spectrum (Figures 10(d) and 10(e)) resembles a purely thermal distribution without any suprathermal tail. The lack of returning particles explains why no oblique waves are seen in the upstream medium of superluminal shocks (compare panels (e), (g), and (i) in Figure 9 with the same panels in Figures 3, 5, and 7).

We remark that the main features discussed here for θ = 45°, and most importantly the absence of shock-accelerated particles, are common to all superluminal configurations up to the perpendicular case θ = 90° investigated via PIC simulations by, e.g., Gallant et al. (1992) and Spitkovsky (2005).

Figure 11 shows the main result of our work, namely, the dependence of the acceleration efficiency on magnetic obliquity: downstream spectra for θ = 0° (blue), 15° (green), 30° (red), and 45° (black) are shown at the final time ω_pt = 9000 of our simulation runs. Subpanel (a) shows the slope p of the fitting power law as a function of magnetic obliquity, while subpanels (b) and (c) present how the fraction of particles and energy stored in the suprathermal tail depends on the magnetic inclination.

Zoom In Zoom Out Reset image size

For subluminal shocks, the downstream energy spectrum consists of a low-energy thermal bump and a high-energy suprathermal tail. The low-energy part, which accounts for most of the energy, is populated by the particles that are directly transmitted downstream, where they are isotropized and thermalized. A fraction of the incoming particles may be scattered or reflected back upstream when they first encounter the shock; they remain in the vicinity of the shock for longer times and may be significantly accelerated (see the high-energy population at the shock in panels (l)–(n) of Figures 3, 5, and 7). If they finally escape upstream, they may contribute to the beam of returning particles which trigger the generation of the upstream oblique waves; if transmitted downstream, they populate the high-energy suprathermal tail of downstream particle spectra.

From Figure 11 it is apparent that the suprathermal tail in subluminal shocks (blue, green, and red lines for θ = 0°, 15°, and 30°, respectively) is much more pronounced than in superluminal shocks (black line for θ = 45°), whose particle spectrum does not significantly deviate from a thermal distribution. The suprathermal tail becomes flatter and more populated as the magnetic obliquity increases within the subluminal range, meaning that shocks close to the superluminality threshold θ_crit ≈ 34° accelerate more efficiently than quasi-parallel shocks. As θ increases from 0° to 30°, the slope of the power-law tail increases from −2.8 ± 0.1 to −2.3 ± 0.1 (subpanel (a)), the fraction of particles in the tail grows from ∼1% to ∼2% (subpanel (b)) and the energy fraction they account for increases from ∼4% to ∼12% (subpanel (c)); at the end of our simulations the spectrum extends up to γ_max ≈ 500 for θ = 0° and reaches γ_max ≈ 2000 for θ = 30°. The negligible suprathermal tail observed for θ = 45° is much steeper than for subluminal shocks; as a conservative upper limit, the fraction of particles and energy it accounts for is smaller by 2 orders of magnitude than in subluminal configurations.

The dependence of the acceleration efficiency on magnetic obliquity should affect the temperature Δγ_mb of the low-energy Maxwellian; we expect that, the larger the fraction of energy stored in the suprathermal tail, the cooler the low-energy bump, if the mean downstream particle Lorentz factor $\bar{\gamma }$ were constant with magnetic obliquity. However, with increasing θ, a larger fraction of the upstream kinetic energy is converted to downstream magnetic energy (see Table 1), so that the average downstream particle energy monotonically decreases with θ. In summary, we expect that, within the subluminal range, Δγ_mb should systematically decrease with θ, both because $\bar{\gamma }$ decreases and because a larger fraction of downstream particle energy is stored in the suprathermal tail; across the boundary between subluminal and superluminal configurations, Δγ_mb should sharply increase, if the drop in acceleration efficiency prevails over the decrease of $\bar{\gamma }$ with θ. Indeed, the thermal spread of the fitting Maxwellian decreases from ∼4.6 for θ = 0° to ∼4.2 for θ = 30°, and then jumps back to ∼4.5 for θ = 45°; in Figure 11, such trend can be seen in the position of the peak of the low-energy bump.

In summary, the high-energy tail becomes harder and more populated as the magnetic obliquity increases within the subluminal range, in agreement with semianalytical (e.g., Kirk & Heavens 1989; Ballard & Heavens 1991) and Monte Carlo (e.g., Bednarz & Ostrowski 1998; Niemiec & Ostrowski 2004; Ellison & Double 2004) calculations. For superluminal shocks, we find that self-generated turbulence is not sufficiently strong to trap a significant fraction of particles in the vicinity of the shock for efficient acceleration. As we show in Section 5, the upstream oblique waves discussed in Section 3.2 regulate injection into the acceleration process for θ ≲ θ_crit. Since no waves are present in the upstream medium of superluminal shocks, due to the absence of the stream of returning particles that trigger their growth, a negligible fraction of the incoming particles enters the acceleration process, and the acceleration efficiency is correspondingly small. The suprathermal tail barely seen for θ = 45° in Figure 11 is probably populated by particles accelerated during a single shock crossing from upstream to downstream, and then advected downstream. We further explore the transition between subluminal and superluminal shocks with regards to downstream particle spectra in Appendix B, where we show that the drop in acceleration efficiency across θ_crit ≈ 34° is rather abrupt.

Finally, we remind that in Figure 11 we are comparing downstream spectra computed in a slab whose distance behind the shock is the same for all magnetic obliquities. This may be misleading, since the acceleration process in subluminal shocks is not in steady state yet (see Section 4.2), and downstream particle spectra at a given time may vary among different longitudinal locations (e.g., compare panels (d) and (e) in Figure 8). The particle energy spectrum in a slab at fixed distance downstream from the shock is affected by the efficiency both in accelerating particles and in transmitting downstream the freshly accelerated particles. Since the latter decreases with increasing obliquity if cross-field diffusion is of minor importance, shocks of larger obliquities may appear, as evaluated from a slab at fixed distance behind the shock, to accelerate less efficiently just because particle diffusion downstream from such shocks is slower. By comparing particle spectra in a slab at constant post-shock distance (Figure 11), we find that subluminal shocks with higher obliquities accelerate more efficiently; this will hold true a fortiori when taking into account the obliquity dependence of the downstream diffusion probability.

Figures 12–14 follow the time evolution of downstream particle spectra for our representative sample of subluminal shocks (θ = 0°, 15°, and 30°, respectively), whereas Figure 15 refers to the superluminal shock θ = 45°. In the upper subpanel of each figure, we plot the time evolution of the electron (blue) and positron (red) maximum Lorentz factor in the downstream slab where the spectrum is computed. The lower subpanel shows the electron (blue) and positron (red) spectra at the final time of our simulations (ω_pt = 9000), together with the best-fitting low-energy Maxwellians (black dashed lines) and high-energy power laws (black solid lines).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

For both subluminal and superluminal magnetic configurations, the low-energy thermal bump does not significantly evolve with time. Instead, the high-energy tail grows with time as more and more particles are shock-accelerated, but its evolution is much more significant for subluminal shocks than for superluminal configurations. For θ = 45° (Figure 15), the particle spectrum for early times (black line) is remarkably thermal, and the minor suprathermal tail visible at later times stops evolving after ω_pt = 6650 (compare the green and red lines in Figure 15). The upper subpanel in Figure 15 confirms that the maximum particle Lorentz factor saturates at γ_max ≈ 120 after ω_pt ≈ 5000. We are then reasonably confident that for θ = 45° our simulation captures the full time evolution of the downstream particle spectrum, so that the claimed inefficiency of particle acceleration in superluminal shocks should still hold true with a longer simulation timespan.

For subluminal shocks, the high-energy tail stretches and flattens with time. Its time evolution is moderate for low obliquities (θ = 0° in Figure 12 and θ = 15° in Figure 13) and more dramatic for magnetic inclinations close to the superluminality threshold (θ = 30° in Figure 14). For θ = 0° the maximum particle Lorentz factor increases from γ_max ≈ 300 at ω_pt = 2250 to γ_max ≈ 500 at ω_pt = 9000, whereas for θ = 30° it grows much faster, from γ_max ≈ 300 to γ_max ≈ 2000. The fraction of particles and energy stored in the tail also increases during our simulations. Therefore, for subluminal obliquities (and especially for angles θ ≲ θ_crit, which show the most significant time evolution), the shock acceleration process is not in steady state yet, and the values we quote for the fractional contribution by the suprathermal tail to the downstream particle and energy census should be taken as lower limits. On the other hand, for late times (compare the green and red lines in Figures 12–14), the shape of the spectrum does not change significantly, apart from the linear increase with time of the upper energy cutoff. Since the power-law spectral index is always p > 2, this means that the values we quote for the acceleration efficiency at ω_pt = 9000 will not significantly increase at later times.

Concerning the time evolution of the maximum particle Lorentz factor, the upper subpanels of Figures 13 and 14 (for θ = 15° and 30°, respectively) show a roughly linear increase with time, whereas for θ = 0° it saturates at γ_max ≈ 500 after ω_pt ≈ 6000 (upper subpanel in Figure 12). By running simulations of parallel shocks in which the first backward jump of the injector happens at different times (we tried 30,000, 40,000, and 50,000 time steps), we observe that, the sooner the injector starts jumping, the quicker the acceleration process saturates. Each jump unphysically removes the returning particles which are farthest upstream. This has two consequences: first, we artificially prevent these particles from returning to the shock and possibly contributing to the downstream spectrum; also, since the stream of returning particles triggers the growth of the upstream oblique waves, which play a major role in the acceleration process (see Section 5), we are unphysically affecting the whole process of shock acceleration. In summary, the spectral time evolution for θ = 0° saturates due to our jumping injector, and the acceleration efficiency and maximum Lorentz factor we report for θ = 0° should be taken as lower limits. Instead, for higher obliquities (θ = 15° and 30°), no unphysical saturation is observed in the spectrum, since the returning particles remain closer to the shock, so that fewer particles are removed by our jumping injector. The jumping injector does not seem to cause other unphysical effects on particle acceleration.

Finally, we point out that the difference between the downstream spectra of positrons and electrons for θ = 30° (lower subpanel of Figure 14, at ω_pt = 9000), which is reminiscent of an analogous difference observed at earlier times (Figure 8(e), at ω_pt = 2250), is an artifact of the small transverse size of our simulations. For magnetic obliquities near θ = 30°, the stream of returning particles, which triggers the growth of the upstream oblique waves, does not propagate very far from the shock. As a consequence, the upstream pattern of electrostatic waves does not move far enough from the shock to let the shock electric potential fluctuate in time (Figure 7(k) for ω_pt = 2250). Particles of one charge will then enter the shock with energies systematically higher than particles of the opposite charge, and the spectral difference between electrons and positrons will persist over time. For smaller obliquities, the stream of returning particles recedes from the shock at higher speeds and the time average of the shock electric potential vanishes, so that no systematic difference is observed between the downstream spectra of electrons and positrons for θ = 0° and θ = 15° (lower subpanels in Figures 12 and 13). For θ = 30°, we have tested that the downstream spectral difference disappears when the simulation box becomes larger (more than 4096 transverse cells) than the characteristic transverse coherence length of the upstream electrostatic oscillations. Therefore, no spectral asymmetry between electrons and positrons is expected in any realistic astrophysical scenario. Moreover, the sum of electron and positron spectra does not vary among runs with different box width, proving that the shock electric potential does not play any appreciable role in the overall acceleration process.⁷

The orbit of a representative high-energy positron from the simulation of a parallel shock is plotted in Figure 16. For ω_pt ≲ 5800, the particle is approaching the shock from upstream with velocity $\mbox{\boldmath {$\beta $}}(t)\approx -\beta _0\,\hat{\mathbf {x}}$ (see x(t), the blue line in panel (c)). Since y(t) stays approximately constant before the particle encounters the shock, the E_y strips stacked in panel (b) are all centered at the same y-location for ω_pt ≲ 5800. Therefore, panel (b) shows for ω_pt ≲ 5800 how the E_y component of the upstream oblique modes (Figure 3(g) for ω_pt = 2250) evolves in time within a fixed longitudinal slice. It is thus evident that such waves are moving toward the shock and their phase velocity in the wall frame is smaller than the speed β₀ ≈ 1 of the incoming particles. The positron energy varies with time even before it encounters the shock (Figure 16(a)), in response to oscillations in γβ_x(t) (the blue line in Figure 16(d)) induced by the electrostatic component of the upstream waves.

From ω_pt ≈ 5800 to ω_pt ≈ 6600, the selected positron scatters between the shock and the upstream medium, and its energy is rapidly boosted up to a Lorentz factor γ ≈ 250 (Figure 16(a)). As the particle energy grows, its gyro-radius also proportionally increases, as seen in panel (c) for y(t) and z(t). A moderate overall energy gain occurs upstream from the shock between ω_pt ≈ 6600 and ω_pt ≈ 7500, and finally at ω_pt ≈ 7500 the positron is transmitted downstream, where its energy remains constant. Most of the accelerated particles end up downstream, but the small fraction that escapes upstream populates the returning stream which triggers the formation of the upstream oblique waves.

As seen in Figure 16(b), the acceleration process involves stochastic bouncings between the shock and the upstream medium; magnetic turbulence in the shock transition layer prevents the particle from being transmitted downstream, whereas the upstream scattering is provided by the oblique waves discussed in Section 3.2, whose amplitude is largest close to the shock. In particular, Figure 16(b) suggests that the positron gains energy whenever its y-velocity β_y(t) (red if negative, yellow if positive) has the same sign as the electromagnetic field E_y (black if negative, white if positive), and decelerates if they have opposite signs. The converse will be true for electrons (not shown). Overall boosting of the selected positron results from favorable scatterings β_y(t)E_y > 0 being more frequent than encounters in which β_y(t)E_y < 0.

The stochastic nature of such wave–particle interactions, which are ultimately responsible for particle acceleration, is a defining feature of DSA. So, particle energization in parallel shocks proceeds via DSA, at the very least because the absence of a uniform background motional electric field rules out SDA. However, in the traditional picture of first-order Fermi acceleration, accelerated particles scatter off magnetic turbulence embedded in the upstream and downstream media. Here, instead, high-energy particles bounce between the upstream medium and the shock front; moreover, the oblique modes responsible for the upstream scattering are not convected with the upstream plasma, as discussed above. Nevertheless, the stochastic character of the acceleration process points toward a diffusive Fermi-like acceleration mechanism. Indeed, as theoretically predicted for diffusive acceleration in relativistic magnetized shocks (e.g., Achterberg et al. 2001), we find that each acceleration cycle between the shock and the upstream medium involves an energy gain Δγ ∼ γ, at least for the limited number of shock crossings (≲10) observed in our simulation.

These results for a parallel shock should also apply to quasi-parallel magnetic configurations, in which the obliquity is small enough that the upstream motional electric field does not significantly affect the acceleration process. However, in Section 5.2 we show that SDA provides most of the energization for obliquity angles as low as θ = 15°, so that DSA is the dominant acceleration mechanism only for θ ≲ 10°.

Figure 17 presents the acceleration process for a representative positron in a shock with θ = 15°. The particle is reflected back upstream on its first interaction with the shock (ω_pt ≈ 5100) and then gains energy at each shock encounter (ω_pt ≈ 5200, ω_pt ≈ 5900, and ω_pt ≈ 8100). The selected positron eventually escapes upstream at time ω_pt ≈ 8200 with Lorentz factor γ ≈ 400, contributing to the beam of returning high-energy particles. However, due to the limited timespan of our simulation, we cannot exclude the possibility that the positron will finally be deflected back and advected downstream.

The time evolution of the positron Lorentz factor plotted in panel (a) suggests that energy gain is always associated with β_y(t) < 0 (red), i.e., when the positron y-velocity has the same sign as the uniform background electric field E_y,u < 0. Whenever the positron gyro-orbit is fully contained in the upstream region (e.g., from ω_pt ≈ 5300 to ω_pt ≈ 5800, or from ω_pt ≈ 6000 to ω_pt ≈ 8000, as well as after the selected positron escapes upstream), the energy gain within the "red" half of the gyro-period is balanced by an approximately equal energy loss during the "yellow" half, since the background motional electric field E_y,u stays approximately constant throughout the gyro-orbit. However, when the positron trajectory intersects the shock front, in the "upstream" portion of its clockwise Larmor gyration around the oblique field the positron y-velocity has always the same sign as the background motional electric field, thus resulting in systematic acceleration. Instead, the "downstream" portion of the orbit does not result in deceleration because in the simulation frame the motional electric field downstream is much smaller than in the upstream.⁸ For electrons (not shown), whose Larmor gyration around the magnetic field is counterclockwise, the y-velocity in the upstream medium will have opposite sign than E_y,u, resulting again in net energization.

In summary, particles are systematically accelerated by the upstream motional electric field whenever their gyro-orbit intersects the shock surface, as expected for the SDA process. SDA is inevitably associated with a ∇B drift along (or opposite, according to the particle charge) the background electric field, due to the magnetic field jump at the shock. This is seen, for the positron in Figure 17, in the slow drift along $-\bf {\hat{y}}$ of y(t) in panel (c). Such a drift will be more significant for higher obliquities, since the magnetic field jump at the shock is larger.

A comparison between x(t) and z(t) in Figure 17(c) suggests that our representative positron is heading upstream along the oblique magnetic field. As compared to a parallel shock, it will take longer to escape ahead of the shock, the more so for higher obliquities. Many gyro-orbits could then intersect the shock during a single shock-crossing of the positron gyro-center, thus resulting in multiple energizations (as happens at ω_pt ≈ 5200 and ω_pt ≈ 5900). The particles that escape upstream are not systematically shock-drift accelerated any longer, until they return to the shock (ω_pt ≈ 8000 in Figure 17), due to deflection by the upstream background fields or scattering off the upstream oblique waves. In that respect, even though the upstream motional field E_y,u is the primary agent for particle energization, DSA cooperates with SDA to trap the accelerated particles at the shock, thus increasing the overall acceleration efficiency. Making use of Equation (2) with a sample of ∼1000 high-energy particles, we are able to confirm that the ratio Δγ_sda/Δγ is independent of the total energy gain Δγ (as predicted by Jokipii 1982), and we estimate the fractional contribution by SDA to be ∼60%. Therefore, for obliquity angles as small as θ = 15°, SDA is already the main acceleration mechanism.

The acceleration process for a representative positron in a shock with θ = 30° is shown in Figure 18. As the obliquity angle approaches the critical boundary between subluminal and superluminal configurations (θ_crit ≈ 34°), particles reflected upstream can hardly escape ahead of the shock. They are effectively trapped at the shock while heading upstream along the oblique magnetic field, since for θ ≲ θ_crit their gyro-center x-velocity is comparable to the shock speed. Before escaping upstream, their orbit may intersect the shock many times (many more than for θ = 15°), thus resulting in repeated energizations via SDA.

As discussed in Section 5.2, for a particle gyrating across the shock any excursion into the upstream region results in efficient SDA by the upstream motional electric field E_y,u (compare panels (a) and (b) in Figure 18), whereas in the downstream medium its energy stays roughly constant since E_y,d ≪ E_y,u in the simulation frame. The corresponding time evolution of the Lorentz factor γ(t) resembles a stairway (see Figure 18(a)), with steps that grow linearly with γ since the time spent upstream is proportional to the particle energy. When a particle completes a full gyro-orbit in the upstream medium (between ω_pt ≈ 7200 and ω_pt ≈ 7600 for the positron in Figure 18), its energy gain is balanced by an equal energy loss. Due to scattering off the incoming oblique waves or regular deflection by the upstream background fields, the particle may be recaptured by the shock, where its acceleration via SDA resumes. The representative positron in Figure 18, as most of the high-energy particles extracted from our θ = 30° simulation, is still being shock-accelerated at the final simulation time ω_pt = 9000; its final Lorentz factor γ ≈ 1000, which is already much larger than in shocks with smaller obliquities, is therefore just a lower limit.

Using Equation (2) with a sample of ∼1000 high-energy particles, we find that on average SDA contributes a fraction ∼90% to the overall energy gain, so that diffusive scattering plays a minor role in the acceleration process for θ = 30°. Thus, with increasing magnetic obliquity in the subluminal range, the acceleration process switches from DSA, which is the only mechanism acting for parallel shocks, to SDA, which for θ ≲ θ_crit operates in an essentially scatter-free regime.

For θ = 30° the upstream oblique waves are important to regulate injection into the acceleration process: all the particles that end up gaining the highest energies have been deflected from their upstream straight path before encountering the shock (see Figure 18(b) at ω_pt ≈ 6500). We have tested that incoming particles enter the acceleration process only if their first shock encounter results in a reflection back upstream, whereas transmitted particles populate the low-energy bump of downstream spectra. In agreement with Ballard & Heavens (1991), we have verified that such reflected particles had encountered the shock with pitch angle cosine μ_ht in the de Hoffmann–Teller frame (de Hoffmann & Teller 1950) which satisfies

$$\begin{equation} |\mu _{{\mbox{\scriptsize {\textsc {ht}}}}}|<\sqrt{1-\frac{B_{\rm u,{\mbox{\scriptsize {\textsc {ht}}}}}}{B_{\rm d,{\mbox{\scriptsize {\textsc {ht}}}}}}}, \end{equation} \tag{ 3 }$$

where B_d,ht/B_u,ht is the magnetic field jump in the de Hoffmann–Teller frame. In other words, only the incoming particles that are deflected by the upstream waves, such that their pitch angle satisfies Equation (3) when encountering the shock, will be reflected back upstream and participate in the SDA process.

In summary, while for quasi-parallel shocks (θ ≲ 10°) the upstream oblique waves generated by the beam of returning particles provide the upstream scattering for DSA, for magnetic obliquities close to θ_crit the same waves mediate the injection process into SDA.

Figures 19 and 20 show the shock structure and positron phase space at time ω_pt = 2250 for in-plane magnetic configurations with obliquity θ = 15° and θ = 30°, respectively. Their out-of-plane counterparts are plotted in Figures 5 and 7, respectively. We find that the shock internal structure does not significantly change by varying the orientation of the field with respect to the simulation plane.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Out-of-plane and in-plane configurations with the same obliquity show a comparable shock velocity and similar transversely averaged density (panel (b)) and magnetic energy (panel (d)) profiles, as relates to the width of the shock transition region (∼10 c/ω_p), the jump across the shock (in agreement with 3D MHD calculations) and the presence of a fast and a slow shock (compare Figures 7 and 20 for θ = 30°). For θ = 15°, the y-averaged magnetic energy profile peaks in the shock layer at ∼30% of the upstream kinetic energy density, both for out-of-plane (Figure 5(d)) and in-plane (Figure 19(d)) configurations, and in both cases it decays downstream on a length scale ∼50 c/ω_p. For θ = 30°, the overall fluid structure does not significantly change between in-plane and out-of-plane configurations, but a few minor differences are present: for the in-plane case, the magnetic energy in the shock layer slightly exceeds the downstream plateau (Figure 20(d)), which was not observed in its out-of-plane counterpart (Figure 7(d)); also, the speed of the slow shock, tracked from the point where the magnetic energy starts to decrease downstream from the main shock, seems to be higher for in-plane than out-of-plane configurations. For both θ = 15° and θ = 30°, the motional electric field (panel (h) in Figures 5 and 7, not shown in Figures 19 and 20) does not drop to zero behind the main shock due to a residual downstream bulk velocity transverse to the upstream flow (β_y,d < 0 and β_z,d < 0 for in-plane and out-of-plane configurations, respectively). The transient electromagnetic precursor wave generated by coherent Larmor bunching at the shock, which showed up in B_z and E_y for out-of-plane oblique configurations, is now seen in the corresponding in-plane fields B_y and E_z with comparable wavelengths and amplitude (further upstream than Figures 19 and 20 show).

The phase space plots in panels (l)–(n) of Figures 19 and 20 are similar to their counterparts in Figures 5 and 7, respectively, provided that out-of-plane γβ_z and γβ_y are compared with in-plane γβ_y and −γβ_z, respectively. As observed for out-of-plane configurations, the returning high-energy particles slide along the upstream magnetic field, which means for the in-plane case that they show large and positive γβ_x and γβ_y (panels (l) and (m) in Figures 19 and 20). The in-plane excess of positrons with γβ_z > 0 in the upstream region close to the shock (panel (n) in Figures 19 and 20) is due to ∇B drift, and it corresponds to the asymmetry toward negative values of γβ_y (panel (m) in Figures 5 and 7, respectively) observed for out-of-plane configurations.

Aside from the similarities between in-plane and out-of-plane results, Figures 19 and 20 improve our knowledge of the internal structure of oblique shocks. The 2D plots of magnetic energy density (panel (c)) and B_z (panel (e)) of Figures 5 and 7 show islands of low magnetic field in the downstream region of the strong shock. Fewer islands are produced at the shock at later times; when formed, each island slowly decays with time at constant x-location, and it disappears when caught up by the slow shock. The corresponding panels in Figures 19 and 20 clarify that such islands are actually tubes of low magnetic energy stretching along the downstream magnetic field, as confirmed also by 3D simulations with 256 cells along each transverse dimension.

Most importantly, the electromagnetic component of the upstream oblique waves, which was mostly seen in B_z and E_y for the out-of-plane case (panels (e) and (g) in Figures 5 and 7), shows up mainly in the same field components also for in-plane configurations. In other words, the wave magnetic field is still mostly perpendicular to the simulation plane, regardless of the orientation of the upstream background magnetic field with respect to the simulation plane. A further study of such modes will therefore have to be fully 3D to capture their true physical nature. What Figures 19 and 20 (panels (e), (g), and (i)) add to our understanding of such waves is that, in the region far upstream where they are triggered, the surface of constant phase is roughly a cone of constant opening angle around the upstream background field.

Figure 21 shows at ω_pt = 9000 a comparison of downstream particle spectra between in-plane (black line) and out-of-plane (red line) configurations: for θ = 15° (panel (a)), the spectrum is remarkably similar between the two cases, whereas a significant difference is seen for θ = 30° (panel (b)), with the high-energy tail for the out-of-plane configuration being much more pronounced than in its in-plane counterpart (and the peak of the low-energy Maxwellian being correspondingly cooler). However, the out-of-plane spectrum for θ = 30° is in good agreement with a 3D simulation having 64 cells along each transverse dimension (green line in panel (b) of Figure 21), which suggests that results from our 2.5D out-of-plane simulations may be confidently applied to realistic 3D scenarios.

Zoom In Zoom Out Reset image size

Concerning the fact that the θ = 30° spectrum at ω_pt = 9000 shows a significant difference at high energies between in-plane and out-of-plane configurations, we report that the spectrum computed in the shock transition layer (not shown) is remarkably similar between the two cases, suggesting that it is not the intrinsic acceleration efficiency to differ but rather the transmission probability of shock-accelerated particles into the downstream medium.