Hostname: page-component-78c5997874-xbtfd Total loading time: 0 Render date: 2024-11-05T04:14:59.927Z Has data issue: false hasContentIssue false

Density jump for oblique collisionless shocks in pair plasmas: physical solutions

Published online by Cambridge University Press:  16 April 2024

Antoine Bret*
Affiliation:
ETSI Industriales, Universidad de Castilla-La Mancha, 13071 Ciudad Real, Spain Instituto de Investigaciones Energéticas y Aplicaciones Industriales, Campus Universitario de Ciudad Real, 13071 Ciudad Real, Spain Center for Astrophysics — Harvard & Smithsonian, Harvard University, 60 Garden St., Cambridge, MA 02138, USA Black Hole Initiative at Harvard University, 20 Garden St., Cambridge, MA 02138, USA
Colby C. Haggerty
Affiliation:
Institute for Astronomy, University of Hawaii, Manoa, 2680 Woodlawn Dr., Honolulu, HI 96822, USA
Ramesh Narayan
Affiliation:
Center for Astrophysics — Harvard & Smithsonian, Harvard University, 60 Garden St., Cambridge, MA 02138, USA Black Hole Initiative at Harvard University, 20 Garden St., Cambridge, MA 02138, USA
*
Email address for correspondence: [email protected]

Abstract

Collisionless shocks are frequently analysed using the magnetohydrodynamics (MHD) formalism, even though MHD assumes a small mean free path. Yet, isotropy of pressure, the fruit of binary collisions and assumed in MHD, may not apply in collisionless shocks. This is especially true within a magnetized plasma, where the field can stabilize an anisotropy. In a previous article (Bret & Narayan, J. Plasma Phys., vol. 88, no. 6, 2022b, p. 905880615), a model was presented capable of dealing with the anisotropies that may arise at the front crossing. It was solved for any orientation of the field with respect to the shock front. Yet, for some values of the upstream parameters, several downstream solutions were found. Here, we complete the work started in Bret & Narayan (J. Plasma Phys., vol. 88, no. 6, 2022b, p. 905880615) by showing how to pick the physical solution out of the ones offered by the algebra. This is achieved by 2 means: (i) selecting the solution that has the downstream field obliquity closest to the upstream one. This criterion is exemplified on the parallel case and backed up by particle-in-cell simulations. (ii) Filtering out solutions which do not satisfy a criteria already invoked to trim multiple solutions in MHD: the evolutionarity criterion, that we assume valid in the collisionless case. The end result is a model in which a given upstream configuration results in a unique, or no downstream configuration (as in MHD). The largest departure from MHD is found for the case of a parallel shock.

Type
Research Article
Copyright
Copyright © The Author(s), 2024. Published by Cambridge University Press

1. Introduction

Shock waves are fundamental processes in fluids. They have been the subject of numerous studies for nearly two centuries (Salas Reference Salas2007). When the frequency of collisions between particles is high, the thickness of the shock front is of the order of a few mean free paths, since binary collisions are ultimately the only microscopic mechanism capable of transferring some of the kinetic energy of the upstream medium into heat in the downstream (Zel'dovich & Raizer Reference Zel'dovich and Raizer2002).

However, in situ observations of the bow shock of the Earth's magnetosphere in the solar wind have shown that its front is approximately a hundred kilometres thick, while the mean free path at the same location is of the order of the Sun–Earth distance (Bale, Mozer & Horbury Reference Bale, Mozer and Horbury2003; Schwartz et al. Reference Schwartz, Henley, Mitchell and Krasnoselskikh2011). Such a shock cannot be mediated by collisions. It is mediated by collective plasma electromagnetic effects (Sagdeev Reference Sagdeev1966). This type of shock is known as a ‘collisionless shock’.

Strictly speaking, collisionless shocks should be studied at the kinetic level, using the Vlasov equation, since the absence of collisions can even make it difficult to define a local velocity, as is the case, for example, in counter streaming systems. Due to the complexity involved in solving the Vlasov equation, collisionless shocks, and in particular the density, temperature or velocity jumps they present, are often interpreted via magnetohydrodynamics (MHD).

Yet MHD relies on hydrodynamics and as such entails the same hypothesis of small mean free path (Gurnett & Bhattacharjee Reference Gurnett and Bhattacharjee2005; Goedbloed, Keppens & Poedts Reference Goedbloed, Keppens and Poedts2010; Thorne & Blandford Reference Thorne and Blandford2017). This hypothesis of small mean free path has 2 consequences that are important for the calculation of the density jump around a shock. The first consequence is that pressure is isotropic, both before and after the shock. In fact, even if a fluid is subjected to pressure anisotropy during its transport from the upstream to the downstream, binary collisions will restore isotropy on a time scale of the order of the collision frequency, well below the macroscopic times involved. The second consequence is that all of the upstream fluid passes into the downstream, along with the matter and momentum it carries.

It turns out that, in a collisionless shock, these 2 consequences can be invalidated (Bret Reference Bret2020). The first, because in the absence of collisions, a plasma can maintain a stable anisotropy in the presence of an ambient magnetic field (Hasegawa Reference Hasegawa1975; Gary Reference Gary1993). The second consequence, because due to the large mean free path, plasma particles can bounce off the shock front or even travel upstream from the downstream, as is the case with accelerated particles (Drury Reference Drury1983; Blandford & Eichler Reference Blandford and Eichler1987).

This article proposes a remedy to the first consequence: How to correct the MHD jump equations so that they can account for an anisotropy in the plasma?

Notably, in the absence of an ambient magnetic field, the Weibel instability ensures isotropy of a collisionless plasma (Weibel Reference Weibel1959; Silva, Afeyan & Silva Reference Silva, Afeyan and Silva2021). Therefore, the MHD jump equations only need to account for anisotropies if a magnetic field is present.

Several authors have derived the MHD jump equations for the non-isotropic case (Hudson Reference Hudson1970; Karimabadi, Krauss-Varban & Omidi Reference Karimabadi, Krauss-Varban and Omidi1995; Erkaev, Vogl & Biernat Reference Erkaev, Vogl and Biernat2000; Gerbig & Schlickeiser Reference Gerbig and Schlickeiser2011). But in all of these works, while the anisotropy of the upstream is considered a free parameter, so is that of the downstream. These equations are therefore unable, on their own, to derive the density jump of a shock whose downstream is not isotropic, because they lack precisely one parameter: the anisotropy of the downstream.

In a recent series of papers, we developed a model that precisely fills this gap. It was first applied to the case of a parallel shock (Bret & Narayan Reference Bret and Narayan2018), i.e. a shock moving parallel to the ambient magnetic field. The assumptions made and the results obtained were confirmed by numerical simulations (Haggerty, Bret & Caprioli Reference Haggerty, Bret and Caprioli2022).

The model was then applied to the case of a perpendicular shock (Bret & Narayan Reference Bret and Narayan2019), and finally to the cases of a switch-on and of an oblique shock (Bret & Narayan Reference Bret and Narayan2022a,Reference Bret and Narayanb). The latter case, that of an oblique shock, is analytically much more complicated than the parallel and perpendicular cases, due to the complexity of the MHD jump equations for an oblique field and anisotropic temperatures.

In Bret & Narayan (Reference Bret and Narayan2022b), hereafter referred to as Paper 1, the algebra of these equations was solved, but the solutions were left unfiltered. As a result, several coexisted for some combinations of the upstream parameters.

Here, Paper 1 is completed by filtering the algebraic solutions so that, for a given combination of upstream variables, there is no more than 1 solution for the downstream.

2. Method

The system under scrutiny is represented in figure 1. It is identical to that of Paper 1. Although the conservation equations we shall use are valid for any upstream temperatures, we here, as in Paper I, only treat the strong sonic case $T_{\perp 1}=T_{\parallel 1}=0$.

Figure 1. System under scrutiny. It is identical to that of Paper 1. Although the conservation equations used are valid for any upstream temperatures, we here, as in Paper 1, only treat the strong sonic case $T_{\perp 1}=T_{\parallel 1}=0$. We work in the frame of reference where $\boldsymbol {v}_1$ is normal to the front.

2.1. Summary of previous works

As previously said, the basic caveat of MHD is that, if a collisionlessly stable anisotropy develops at the front crossing, MHD itself cannot derive it. The jump of quantities like the density is therefore under-determined.

For completeness, we now briefly recall the results obtained in previous works.

In Bret & Narayan (Reference Bret and Narayan2018) we presented a model capable of solving this issue for a parallel shock. We reasoned that, as it crosses the front, the plasma is compressed in the direction parallel to the motion. As a consequence, its parallel temperature increases while its perpendicular temperature remains constant. The 3 MHD conservation equations (matter, momentum, energyFootnote 1 ) are therefore completed by

(2.1)\begin{equation} T_{{\perp} 2} = T_{{\perp} 1}, \end{equation}

allowing us to derive the 4 downstream unknowns ($n_2,v_2,T_{\perp 2},T_{\parallel 2}$), in terms of the upstream variables. Note that, here, the perpendicular direction is common to the flow and the field.

Now, the state of the downstream resulting from the conservation of $T_\perp$ may be stable, or not. If it is stable, then this is the end state of the downstream. If it is unstable, the plasma migrates to the instability threshold.Footnote 2 Imposing marginal stability then provides again a fourth equation allowing us to fully determine the properties of the downstream.

Bret & Narayan (Reference Bret and Narayan2018), as well as every subsequent works of ours on the same model, is limited to pair plasmas. The reason for this is that the equality of the mass of the species allows us to consider only one parallel and one perpendicular temperature. In an electron/ion plasma where electrons and ions are heated differently in the shock, a 4 temperature model would be required (Guo, Sironi & Narayan Reference Guo, Sironi and Narayan2017, Reference Guo, Sironi and Narayan2018). Yet, since the model eventually relies on macroscopic physics, it should also apply to electron/ion plasmas, as preliminary particle-in-cell (PIC) simulations seem to indicate (Shalaby Reference Shalaby2024).

The model predicted, for a strong sonic parallel shock, a density jump going from 4 to 2 in the high field regime, a prediction successfully confirmed by PIC simulations in Haggerty et al. (Reference Haggerty, Bret and Caprioli2022). Such a result stands in contrast with MHD, where the density jump should always be 4, regardless of the field strength.Footnote 3

The perpendicular case was treated in Bret & Narayan (Reference Bret and Narayan2019). There, the direction perpendicular to the flow is eventually parallel to the field, so that the counterpart of (2.1) is

(2.2)\begin{equation} T_{{\parallel} 2} = T_{{\parallel} 1}. \end{equation}

The switch-on shocks, where the field is oblique in the downstream only, was treated in Bret & Narayan (Reference Bret and Narayan2022a). The model has also been solved for a parallel or a perpendicular shock, with an anisotropic upstream (Bret Reference Bret2023a,Reference Bretb).

Finally, the general case where the field may be oblique in both the upstream and the downstream was treated in Bret & Narayan (Reference Bret and Narayan2022b), namely, Paper 1.

In Bret & Narayan (Reference Bret and Narayan2022a,Reference Bret and Narayanb), the closure of the MHD jump equations was achieved through an ansatz interpolating between (2.1) and (2.2). In the limit of a cold upstream with $T_1=0$, which is the regime treated in Paper 1 and hereafter, the ansatz reads

(2.3)\begin{equation} \left.\begin{array}{c@{}} T_{{\parallel} 2} = T_e\cos^2\theta_2, \\ T_{{\perp} 2} = \dfrac{1}{2}T_e\sin^2\theta_2, \end{array}\right\} \end{equation}

where $T_e$ is a parameter determined when solving the equations and $\theta _2$ is the angle of the downstream field with the shock normal (see figure 1). Equations (2.3) correctly reduce to (2.1) and (2.2) in their respective limits since $\theta _2=0$ for the parallel case, while $\theta _2={\rm \pi} /2$ for the perpendicular one. Such a scheme guarantees that both perpendicular temperatures are equal, as required by the Vlasov equation (Landau & Lifshitz Reference Landau and Lifshitz1981, § 27). Also, the total thermal energy in the 3 directions sums up to $k_BT_e$, where $k_B$ is the Boltzmann constant.

In summary, and since the 2 instabilities involved are the firehose and the mirror instabilities (see § 4), our model can be stated as follow:

  1. As the plasma goes through the shock front, its temperature normal to the flow is conserved for the parallel and perpendicular cases. This translates directly to (2.1) and (2.2), respectively. For the oblique case with cold upstream $T_1=0$, (2.3) interpolate between these 2 extremes.

  2. The resulting state of the plasma in the downstream is called ‘stage 1’.

  3. If the downstream field is strong enough to stabilize stage 1, then this is the end state of the downstream.

  4. If the downstream field is too weak to stabilize stage 1, then

    1. If stage 1 is firehose unstable, it migrates to the firehose instability threshold. This is ‘stage 2 – firehose’.

    2. If stage 1 is mirror unstable, it migrates to the mirror instability threshold. This is ‘stage 2 – mirror’.

2.2. Present work

Paper 1 has 3 kinds of limitations:

  1. (a) It is restricted to strong sonic shocks, namely $T_1=0$, and to non-relativistic pair plasmas. These restrictions are still considered here.

  2. (b) It considers the simplest expressions for the Alfvén velocity and the stability criterion of the instabilities involved in our model. Yet, more accurate expressions are required in an anisotropic plasma. The present work accounts for one of them.

  3. (c) It only presents the allowed solutions to the conservation equations, plus (2.3). It does not filter these solutions according to their physical relevance. Such a filtering is the main goal of the present work.

Our purpose here is to deal with points (b) and (c) above.

Besides the variables explained on figure 1, we shall use the following dimensionless parameters:

(2.4)\begin{equation} \left.\begin{array}{c@{}} \displaystyle r = \dfrac{n_2}{n_1}, \\ \displaystyle A_i = \dfrac{T_{{\perp} i}}{T_{{\parallel} i}}, \\ \displaystyle M_{A,i} = \dfrac{v_i}{v_{A,i}}, \\ \displaystyle \sigma = \dfrac{B_1^2/8{\rm \pi}}{\dfrac{1}{2}n_1 m v_1^2} \equiv \dfrac{1}{M_{A,1}^2}, \end{array}\right\} \end{equation}

where $v_{A,i}$ is the Alfvén velocity

(2.5)\begin{equation} v_{A,i} = \frac{B_i}{\sqrt{4 {\rm \pi}n_i m}}. \end{equation}

The parameter $\sigma$ is frequently used in simulations of collisionless shocks like Haggerty et al. (Reference Haggerty, Bret and Caprioli2022), while the Alfvén Mach number is common in the MHD shock literature.

In addition to the Alfvén Mach number defined above, we shall often use in the sequel its following variant:

(2.6)\begin{equation} M_{Aix} \equiv \frac{v_i \cos\xi_i}{v_{A,i}\cos\theta_i}. \end{equation}

It compares the projection of the flow velocity along the shock normal ($x$ axis) with the projection of the Alfvén velocity, still along the shock normal.

Since the road map for solving our model is eventually the one already used in MHD, we start by reminding the reader how it applies there.

3. Isotropic MHD results and evolutionarity

One of the criteria used here to filter out some solutions is the so-called ‘evolutionarity criterion’, already present in isotropic MHD. It is therefore convenient to present how it operates there.

The MHD conservation equations for isotropic temperatures are reported in Appendix A. They can be used to derive an expression of the downstream field angle $\theta _2$ in terms of the upstream quantities only. Namely, the quantity

(3.1)\begin{equation} T_2 \equiv \tan\theta_2, \end{equation}

is a root of the polynomial (A9) in Appendix A. The MHD density jump is then given by

(3.2)\begin{equation} r = \frac{M_{A1}^2 T_2}{M_{A1}^2 \tan \theta_1+T_2 \cos ^2\theta_1-\sin \theta_1 \cos \theta_1}, \end{equation}

in terms of the upstream Alfvén Mach number (2.4).

Although such a derivation of the MHD density jump is uncommon in the literature, it mimics the derivation of the jump in our model. It is therefore pedagogically valuable for the present work.

Figure 2 presents the solutions for 3 different angles $\theta _1$. The upper row shows all the possible solutions. For $\theta _1=0$ and $M_{A1x}\in [1,2]$, there are 2 MHD branches, the lower one pertaining to the switch-on solutions. For $\theta _1=0.175 {\rm \pi}/2$ there are even 3 solutions for $M_{A1x}\in [1, 1.34]$.

Figure 2. Upper row: density jump $r$ for different angles $\theta _1$ in terms of $M_{A1x}=M_{A1}/\cos \theta _1$, defined by (2.6). Sometimes there is more than 1 solution. Second row: angle $\theta _2$ of the field $\textbf {B}_2$ with the $x$ axis (see figure 1). Third row: evolutionarity criterion. Some branches, namely the ones crossing the shaded areas, are to be excluded. Lower row: same as upper row, but without the branches excluded by the evolutionarity criterion. There is now but 1 solution for a given $M_{A1x}$, or none. The black circles for $\theta _1=0.5{\rm \pi} /2$ are the results of our MHD simulations (see § 3.1).

Which one will the shock choose? This question is at the heart of this work. Let us now see how it is solved in MHD.

The MHD answer relies on the notion of shock ‘evolutionarity’, which has been discussed several times in the literature (e.g. Kennel, Blandford & Wu Reference Kennel, Blandford and Wu1990; Farris et al. Reference Farris, Russell, Fitzenreiter and Ogilvie1994; Falle & Komissarov Reference Falle and Komissarov2001; Kulsrud Reference Kulsrud2005). For given upstream and downstream boundary conditions, the MHD Rankine–Hugoniot jump conditions give a unique solution for fast shocks, where the flow speed on both sides of the shock are super-Alfvénic, and also for slow shocks, where the fluid is throughout sub-Alfvénic. These two types of shocks are stable and are described as evolutionary. Four other potential shock types, each of which has super-Alfvénic upstream fluid and sub-Alfvénic downstream fluid, have no unique solutions to the Rankine–Hugoniot relations. In these shocks, the transverse magnetic field switches sign across the shock, and the fluid equations do not provide the correct number of Alfvén waves to handle the field flip (Falle & Komissarov Reference Falle and Komissarov1997; Kulsrud Reference Kulsrud2005). Such shocks will not arise naturally from generic initial conditions. Even if they are artificially set up, they will tend to deviate quickly from their initial configuration, typically splitting into two shocks, one fast and the other slow. These ‘forbidden’ shocks do not satisfy the evolutionarity condition. They are called ‘intermediate’ (Falle & Komissarov Reference Falle and Komissarov1997) or ‘extraneous’ (Kulsrud Reference Kulsrud2005) shocks.

The MHD solutions presented on the upper row of figure 2 need therefore to be filtered according to the aforementioned evolutionarity criterion. The upstream Alfvén Mach number to consider for the evolutionarity analysis is not the one given by (2.4), but its variant (2.6) instead. With this definition, the downstream Alfvén Mach number reads

(3.3)\begin{equation} M_{A2x}^2 = 1+\frac{\tan \theta_1 (M_{A1}^2 \sec ^2\theta_1-1)}{T_2}. \end{equation}

It is with definition (2.6) that the Alfvén Mach number of the MHD switch-on shocks, see left plot of third row in figure 2, is found equal to 1 (Goedbloed et al. Reference Goedbloed, Keppens and Poedts2010, p. 853). These modes, by definition, do not have $\theta _2=0$ (second row in figure 2) nor $\xi _2=0$. Corrections (2.6) to the Alfvén Mach number (2.4) are therefore important here.

The value of $M_{A2x}$ is plotted in terms of $M_{A1x}$ in the third row of figure 2. The forbidden, non-evolutionary, zones have been coloured in the plots. They feature the non-evolutionary transition just described, namely super-Alfvénic $\rightarrow$ sub-Alfvénic, together with the reverse transition sub-Alfvénic $\rightarrow$ super-Alfvénic. As a result, the MHD branches going through these regions are non-evolutionary, hence not physical solutions of the MHD problem.

The lowest row of figure 2 is the result of this evolutionary filter applied to the upper row. There is now but 1 solution for a given $M_{A1x}$, or none, as there are $M_{A1x}$-gaps where no solutions appear. Regardless of the initial conditions, the shock formed is never found with a $M_{A1x}$ inside such gaps like, for example, $M_{A1x}\in [1 , 1.4]$ for $\theta _1=0.5{\rm \pi} /2$ (see figure 2, bottom right plot).

We shall assume in the sequel that the evolutionarity criteria also applies in the collisionless case. This will have to be checked in future works.

3.1. Magnetohydrodynamics simulations

In order to illustrate these theoretical results, we ran some MHD simulations with the code KORAL (Sa̧dowski et al. Reference Sa̧dowski, Narayan, Tchekhovskoy and Zhu2013, Reference Sa̧dowski, Narayan, McKinney and Tchekhovskoy2014). KORAL is designed for multi-dimensional radiative MHD simulations in general relativistic space–times. However, if we turn off the radiation module as well as special and general relativity, the code reduces to a multi-dimensional non-relativistic MHD code. This version of the code was used here. Some of the results for $\theta _2$, $M_{A2x}$ and $r$ in terms of $M_{A1x}$ for upstream $\theta _1=0.5{\rm \pi} /2$, are pictured in figure 2 by the black circles. They perfectly line up with the theory and avoid the predicted gap in the solutions for $M_{A1x}\in [1 , 1.4]$.

Figure 3 shows the result of a simulation of a shock in a non-evolutionary case. At $t=0$, two cold fluids of identical density and opposite velocities $\pm 1.5 \boldsymbol {x}$ collide (blue lines). The magnetic field has modulus unity and $\theta _1=0.5{\rm \pi} /2$. At $t=1.5$, two shocks formed instead of one. The non-evolutionary case gives rise to 2 sub-shocks, as predicted in Kulsrud (Reference Kulsrud2005) for example.

Figure 3. Magnetohydrodynamics simulation in which the evolutionarity criterion results in the formation of 2 sub-shocks instead of one. The dashed blue lines show the initial fluid configuration at $t=0$, when two fluid columns collide at $x=0$. The solid red lines show the configuration at $t=1.5$. Each fluid column develops two shocks, a fast shock at $|x|=1.58$ and a slow shock at $|x|=0.34$. The little dip in density $n$ at $x=0$ for $t=1.5$ is a numerical artefact.

4. Instabilities involved and modified Alfvén velocity

As previously stated, the starting point for our model are the MHD conservation equations with anisotropic temperatures. They have been derived in Hudson (Reference Hudson1970) and are reported in Appendix B with the notations of Bret & Narayan (Reference Bret and Narayan2022a,Reference Bret and Narayanb).

For stage 1, we solve them imposing relations (2.3). This defines stage 1. Then, mirror or firehose stability of stage 1 has to be assessed. These instabilities are discussed here.

Now, evolutionarity involves the downstream Alfvén velocity (projected onto the shock normal). As long as the plasma is isotropic, the Alfvén velocity in given by (2.5). As a result, the downstream Alfvén Mach number for the MHD switch-on shocks is exactly 1.

Yet, in an anisotropic plasma, it was found in Abraham-Shrauner (Reference Abraham-Shrauner1967) that the Alfvén velocity readsFootnote 4

(4.1)\begin{equation} c_A ={\pm} v_A \sqrt{S_\perp{-} S_\parallel{+}1} \cos\theta, \end{equation}

where $v_A$ is still given by (2.5) and

(4.2)\begin{equation} S_{{\perp},\parallel} = \frac{n k_B T_{{\perp},\parallel}}{B^2/4{\rm \pi}}. \end{equation}

In (4.1), the $\pm$ sign refers to waves propagating along the field, or in the opposite direction. As shall be checked in the sequel (§§ 5.1 and 5.2 and figure 5), with this correction to the Alfvén velocity, the switch-on solutions of our model have $M_{A2x}=1$, exactly like in MHD.

In addition, the quantity below the square root can become negative. Such a situation indicates that the Alfvén waves become unstable, which corresponds to the firehose instability. For this to happen, $S_\perp - S_\parallel +1 < 0$ is required, which can be cast under the form

(4.3)\begin{equation} A \equiv \frac{T_\perp}{T_\parallel} < 1 - \frac{2}{\beta_\parallel}, \end{equation}

whereFootnote 5

(4.4)\begin{equation} \beta_\parallel{=} \frac{n k_B T_\parallel}{B^2/8{\rm \pi}}. \end{equation}

The criterion (4.3) is therefore the one we shall use in the sequel, in order to preserve the inner coherence with the switch-on solutions of our model having $M_{A2x}=1$. This criterion slightly differs, by the factor 2, from the one commonly used in the literature and in Paper 1.

The Solar Wind data are a key test for the threshold of the firehose instability. They are indeed compatible with both criteria, namely with or without the factor 2 (Hellinger et al. Reference Hellinger, Trávníček, Kasper and Lazarus2006; Bale et al. Reference Bale, Kasper, Howes, Quataert, Salem and Sundkvist2009; Maruca, Kasper & Bale Reference Maruca, Kasper and Bale2011).

The other instability considered in relation to the stability of stage 1, is the mirror instability. While the firehose instability occurs for too low an anisotropy $T_\perp /T_\parallel$, the mirror instability occurs for too high an anisotropy. The standard threshold given in the literature reads (Hasegawa Reference Hasegawa1975; Gary Reference Gary1993)

(4.5)\begin{equation} A > 1 + \frac{1}{\beta_\parallel}. \end{equation}

Yet, a different criterion is given in Abraham-Shrauner (Reference Abraham-Shrauner1967), namely

(4.6)\begin{equation} A > 1 + \frac{1}{A ~\beta_\parallel}. \end{equation}

However, while stage 2-firehose remains analytically tractable using (4.3), stage 2-mirror is not when using (4.6) instead of (4.5). Therefore, we adopt criterion (4.3) for the firehose instability and keep (4.5) for mirror. Note in this respect that (i) the Solar Wind data are compatible with both criteria and (ii), the inner coherence of the model does not impose a specific mirror criteria, as is the case for firehose.

Figure 4 pictures the various criteria commented here for the mirror and the firehose instabilities. Even if the corrected criterion for the mirror instability is functionally different than the one without the correction, the end result is qualitatively the same, and remains compatible with the Solar Wind data.

Figure 4. Stability diagram. The solid lines picture the criterion used here. The dashed lines pertain to other criteria discussed in § 4.

The two stability domains are disconnected, so that the plasma cannot be unstable to both at once. There is no possible competition between them.

5. Branch selection

Having specified the instabilities involved in the transition from stage 1 to stage 2, together with their respective thresholds, we can proceed to the filtering of the solutions our model offers. In this respect, it is instructive to single out the case $\theta _1=0$.

5.1. Case $\theta _1=0$

Figure 5(a) shows the solutions offered by our model for $\theta _1=0$, without any filtering. It displays the various solutions with stable stage 1, plus a branch, in green, corresponding to stage 2-firehose because stage 1 is firehose unstable for $\sigma \in [0 , 0.5]$. For $\sigma \in [0.5,1]$, there are up to 3 possible solutions: one has $r=2$ and $\theta _2 = 0$, and the other two pertain to the switch-on solutions in red, with $\theta _2 \neq 0$.

Figure 5. (a) All solutions offered by our model for $\theta _1=0$. (b) All the solutions fulfil the evolutionarity criterion. Note that this is different from the isotropic MHD problem, where some solutions for $\theta _1=0$ enter the shaded forbidden zone (figure 2, row 3). Note also that, with definition (4.1) of the Alfvén speed, $M_{A2x}=+\infty$ on the firehose threshold. This is why the green branch on the left is sent to $M_{A2x}=+\infty$ on the centre plot. Notably, $M_{A2x}=1$ for the switch-on shocks of our model (red segment on centre plot). For this analysis, the horizontal axis must be $\sigma _x=\sigma \cos ^2\theta _1$, which for $\theta _1=0$ makes no difference. (c) End result once the ‘$\theta _2$ closest to $\theta _1$’ filter has been applied. This is the result found in Bret & Narayan (Reference Bret and Narayan2018) and Haggerty et al. (Reference Haggerty, Bret and Caprioli2022).

Could the evolutionarity criterion help filtering them (the answer is ‘no’)? Figure 5(b) answers the question. Here is plotted the downstream Mach number $M_{A2x}$ of each solution, as defined by (2.6), where the Alfvén speed has been corrected according to (4.1). Our switch-on solutions have exactly $M_{A2x} = 1$, while the others also satisfy the evolutionary criterion.Footnote 6 In line with the adequate definition (2.6) of the Alfvén Mach number for evolutionarity analysis, we use on the horizontal axis the parameter

(5.1)\begin{equation} \sigma_x \equiv \sigma \cos^2\theta_1, \end{equation}

which, for the present case $\theta _1=0$, makes no difference.

Note that, with definition (4.1) of the Alfvén speed, $M_{A2x}=+\infty$ on the firehose threshold. This is the reason why the green branch on the left is sent to $M_{A2x}=+\infty$ on the centre plot. As a consequence, each time the system eventually settles in stage 2-firehose with $\sigma _x < 1$ (i.e. $M_{A1x}>1$), it is evolutionary.

As a conclusion, in the interval $\sigma \in [0.5,1]$, there are three solutions for stage 1, and all three satisfy the evolutionarity criterion. The evolutionarity criterion by itself is therefore not sufficient to trim the number of solutions down to 1, or even 0.

5.2. Particle-in-cell simulations for $\theta _1=0$

In order to check which solution the system chooses, we ran a series a PIC simulations for various $\sigma$ values, performed with the fully kinetic three-dimensional PIC code, TRISTAN-MP (Buneman Reference Buneman1993; Spitkovsky Reference Spitkovsky2005).Footnote 7 The surest way to tell whether the shock is of the switch-on type or not, is to plot the perpendicular component $B_{\perp 2}$ of the field in the downstream. According to (B2), with $\theta _1=0$ it simply reads

(5.2)\begin{equation} B_{{\perp} 2} = B_1 \tan \theta_2. \end{equation}

Figure 6(a) shows the expected value of $B_{\perp 2}/B_1$, in red for the switch-on solutions with $\theta _2\neq 0$, and in black for the solution with $\theta _2=0$. The circles show the values of $\sigma$ which have been simulated. Figure 6(b) shows the results of the PIC simulations, with the shock density profile on top, and the ratio $B_{\perp 2}/B_1$ below. Besides some perturbations around the shock front at low $\sigma$, $B_{\perp 2}/B_1$ does not depart from 0 in the downstream for any $\sigma$. These perturbations are due to particle acceleration and back-reaction in the precursor, and from the instabilities triggered by the interaction of the fast upstream flow with the front (Sironi, Spitkovsky & Arons Reference Sironi, Spitkovsky and Arons2013; Brown et al. Reference Brown, Juno, Howes, Haggerty and Constantinou2023).

Figure 6. (a) Expected value of $B_{\perp 2}/B_1$, in red for the switch-on solutions with $\theta _2\neq 0$, and in black for the solution with $\theta _2=\theta _1=0$. The circles show the values of $\sigma$ which have been simulated. (b) Results of the PIC simulations, with the shock density profile on top, and the ratio $B_{\perp 2}/B_1$ below. The $\sigma$ value for a curve is given by the circle of the same colour on the left plot. The same results are obtained for $\theta _1=2^\circ$. The flow is along the $x$ axis and $d_e=c/\omega _p$ is the electron inertial length.

Simply put, PIC simulations consistently discard the switch-on solutions. The same results are obtained for $\theta _1=2^\circ$, so that we are not witnessing a singular behaviour fruit of a perfect, hence unrealistic, symmetry. A similar pattern was observed in the relativistic regime in Bret et al. (Reference Bret, Pe'er, Sironi, Sa̧dowski and Narayan2017).

The situation for $\theta _1 \sim 0$ is here markedly different from the MHD case. In MHD, where isotropy is imposed, the evolutionarity criterion imposes switch-on shocks within some $\sigma$-range, as explained in § 3. In our model, where an anisotropy is driven by the field, the evolutionarity criterion does not impose switch-on solutions, while PIC simulations consistently choose the non-switch-on ones.

Arguably, this explains why so few detections of switch-on shocks have been made in the Solar System. Indeed, among the thousands of shocks observed (Farris et al. Reference Farris, Russell, Fitzenreiter and Ogilvie1994; Russell & Farris Reference Russell and Farris1995; Balogh & Treumann Reference Balogh and Treumann2013, § 2.3.6.) only one ‘possible’ detection of an interplanetary switch-on shock was reported in Feng et al. (Reference Feng, Lin, Chao, Wu, Lyu and Lee2009). Also, Russell & Farris (Reference Russell and Farris1995) reported the detection of only one switch-on shock among the International Sun-Earth Explorer data.Footnote 8

Still, what about these exceptions, since our collisionless scenario should rule them out? Several explanations are possible:

  1. The detections were faulty. Hence the adjective ‘possible’ associated with one of them.

  2. We only solve our model for a cold upstream, namely $T_1=0$. Maybe a finite $T_1$ would have our model allowing for some switch-on shocks.

  3. The ansatz (2.3) is not accurate enough, and a better version would allow for some switch-on shocks.

At any rate, PIC simulations and observations tell switch-on shocks in collisionless plasma are rare. We therefore propose the following criterion allowing us to choose between various solutions: the system chooses the solution with $\theta _2$ closest to $\theta _1$.

Note that this ‘$\theta _2$ closest to $\theta _1$’ criterion stems from our PIC simulations, as well as others of parallel shocks in pair (Bret et al. Reference Bret, Pe'er, Sironi, Sa̧dowski and Narayan2017; Haggerty et al. Reference Haggerty, Bret and Caprioli2022) or electron/ion (Niemiec et al. Reference Niemiec, Pohl, Bret and Wieland2012; Zeković Reference Zeković2019) plasmas, where the same, non-switch-on branch, is consistently chosen. This is why we used the verb ‘propose’. Its robustness on longer time scales, or other shock geometries, is beyond the scope of this work and should be assessed in further works.

Figure 5(c) eventually shows the end result once the ‘$\theta _2$ closest to $\theta _1$’ filter has been applied to the $\theta _1=0$ case. There is now but one solution for any $\sigma$, which indeed is the one found in Bret & Narayan (Reference Bret and Narayan2018) and Haggerty et al. (Reference Haggerty, Bret and Caprioli2022).

5.3. General algorithm for branches selection

We may now lay out the general algorithm for branches selection. The criteria used to eliminate some are, applied in this order:

  1. (a) Exclude stage 1 solutions with $r<1$. Then select the one with the $\theta _2$ closest to $\theta _1$.

  2. (b) In case stage 1 is unstable, is the anisotropy $A_2$ of a given stage 2 solution, negative? This was already implemented in Paper 1 and allows us to eliminate some stage 2 -firehose and -mirror solutions. This never happens with stage 1 since it has, by design from (2.3), $A_2 = \frac {1}{2}\tan ^2\theta _2$.

  3. (c) Does the resulting solution fulfil the evolutionarity criterion ?

The evolutionarity criterion is applied last because it operates on time scales related to the propagation of the shock,Footnote 9 whereas the other criteria operate on much shorter time scales, related to plasma instabilities.

The algorithm is eventually represented by the flowchart in figure 7. In case stage 1 is mirror unstable for a given pair $(\sigma,\theta _1)$, there can be a degeneracy in the choice of the stage 2-mirror states for the same pair $(\sigma,\theta _1)$. In such a case, we choose the stage 2-mirror with the $\theta _2$ closest to the $\theta _2$ of the stage 1 it comes from, as was already done in Bret & Narayan (Reference Bret and Narayan2022b). Such a situation never happens with stage 2-firehose.

Figure 7. Flowchart of the resolution of our model. In case stage 1 is mirror unstable for a pair $(\sigma,\theta _1)$, there can be a degeneracy in the choice of stage 2-mirror for the same pair $(\sigma,\theta _1)$. In such a case, we choose the stage 2-mirror with the $\theta _2$ closest to the $\theta _2$ of the stage 1 it comes from, as was already done in Bret & Narayan (Reference Bret and Narayan2022b). Such a situation never happens with stage 2-firehose.

According to the flowchart on figure 7, there is necessarily only 1 solution left before applying the evolutionarity criterion. There is therefore no need to apply the ‘$\theta _2$ closest to $\theta _1$’ filter in (c), since only 1 branch can make it to this stage. Yet, this does not mean the evolutionarity filter is useless, as it can simply forbid the existence of a solution in some $\sigma$ range, as is the case in MHD (see bottom row of figure 2).

6. Results

While the fruit of our algorithm has already been explained for $\theta _1=0$, it is interesting to detail how it unfolds for an oblique field, like for example $\theta _1 = 0.3{\rm \pi} /2$.

Figure 8(a) shows all possible density jumps for stage 1. Notably, $r < 1$ for $\sigma \in [1, 1.3]$. The corresponding branches are then eliminated in figure 8(b).

Figure 8. Detailed progression of the application of the filtering process described by the flow chart in figure 7, for $\theta _1=0.3{\rm \pi} /2$. In (d), the green colour indicates that these stage 1 solutions are firehose unstable. Then in (g), the green colour shows the density jump of stage 2-firehose. The red branch in ( f) with $A_2 \sim -1$ discards stage 2 firehose solutions with $r >4$, not shown in (e).

For these stage 1 solutions with $r > 1$, figure 8(c) then shows the corresponding values of $\theta _2$. The horizontal red line is at $\theta _2=\theta _1$. The result of the ‘$\theta _2$ closest to $\theta _1$’ selection rule is displayed in figure 8(d).

In figure 8(d), the green colour indicates that some stage 1 solutions at low $\sigma$ are firehose unstable. We therefore need to examine stage 2-firehose solutions. All of them are displayed in figure 8(e). Yet, figure 8f) shows that the lower branch has $A_2 < 0$ and needs to be eliminated. As a consequence, the only physical stage 2-firehose solution available when stage 1 is firehose unstable, is the upper branch.

Replacing then the firehose unstable stage 1 solutions of figure 8(d), by the corresponding stage 2-firehose solutions, gives figure 8(g).

We finally need to apply the evolutionarity test to the solutions of figure 8(g). This is done in figure 8(h) where the downstream Alfvén Mach number $M_{A2x}$ has been computed for the solution plotted in figure 8(g), with the forbidden zones shaded. Note that, for such an analysis, the horizontal scale has to be $\sigma _x=\sigma \cos ^2\theta _1$. Then and only then has the evolutionary analysis some branches passing exactly through the point $(1,1)$, as in figure 8(h).

The stage 2-firehose branch visible in green in figure 8(g) at small $\sigma$ passes the evolutionarity test since it has $\sigma _x < 1$ and $M_{A2x} = +\infty$. The analysis shows that some stable stage 1 solutions do not pass the evolutionarity test. As a result, figure 9 for $\theta _1=0.3{\rm \pi} /2$ displays a gap without solution for $\sigma \in [1,1.17]$, that does not show in figure 8(g).

Figure 9. Result of the filtering process elaborated in this work for various angles $\theta _1 \in [0, {\rm \pi}/2]$. All of them but $\theta _1=0.3{\rm \pi} /2$ are the ones considered in figure 8 of Bret & Narayan (Reference Bret and Narayan2022b). The grey dashed lines picture the MHD result, evolutionary filtered. The green solutions pertain to (stage 2-firehose, green), and the red ones to (stage 2-mirror, red).

The end result of such a filtering process is eventually shown in figure 9 for several angles $\theta _1$ between 0 and ${\rm \pi} /2$. For comparison, all of them but $\theta _1=0.3{\rm \pi} /2$, are the ones considered in figure 8 of Bret & Narayan (Reference Bret and Narayan2022b). The grey dashed line pictures the MHD result, evolutionary filtered.

Figure 9 eventually interpolates between the parallel case treated in Bret & Narayan (Reference Bret and Narayan2018), and the perpendicular one treated in Bret & Narayan (Reference Bret and Narayan2019). While § 5 showed that some filtering is needed to get to the end result for the parallel case, no filtering at all is required for the perpendicular case. Figure 9 for $\theta _1=0.9{\rm \pi} /2$ gives the result already found in Bret & Narayan (Reference Bret and Narayan2019), without any filtering but the $A_2 < 0$ for the mirror solutions. The main reason for this is that for the perpendicular case, there is but one stage 1 solution, which, when stable, fulfils the evolutionarity criterion.

7. Conclusion

Applying the MHD formalism to analyse collisionless shocks may be problematic in the presence of an ambient magnetic field, capable of stabilizing pressure anisotropies. In a series of recent papers, we developed a method allowing us to determine the downstream anisotropy of a collisionless shock (Bret & Narayan Reference Bret and Narayan2018, Reference Bret and Narayan2019, Reference Bret and Narayan2022a,Reference Bret and Narayanb; Bret Reference Bret2023a,Reference Bretb). The anisotropy can then be included in the MHD conservation equations for anisotropic pressures, and the modified density jump derived. The case of a parallel shock has been successfully tested by PIC simulations in Haggerty et al. (Reference Haggerty, Bret and Caprioli2022), confirming that for a strong enough field, the density jump can go from 4, the expected MHD value, to only 2, the anisotropy corrected one.

Once the density jump is found, all the other jumps such as pressure, temperature, entropy, etc. can be straightforwardly derived (Bret Reference Bret2021).

As can be seen in figure 9, our results differ from the isotropic MHD ones in 2 ways:

  1. The ranges of $\sigma$ without solutions differ.

  2. For $\sigma$ values with a solution, our density jump is usually, although not always, lower than the MHD one.

Overall, our results and the MHD ones bear a ‘family resemblance’, the largest discrepancy being found for the parallel case $\theta _1=0$. The predicted large reduction of the density jump could have important consequences for particle acceleration, since their index scales like $(r-1)^{-1}$ (Axford, Leer & Skadron Reference Axford, Leer and Skadron1977; Bell Reference Bell1978a,Reference Bellb; Blandford & Ostriker Reference Blandford and Ostriker1978; Caprioli, Haggerty & Blasi Reference Caprioli, Haggerty and Blasi2020; Haggerty & Caprioli Reference Haggerty and Caprioli2020).

Besides the strong sonic shock assumption of the present work, namely $P_1=0$, its main limitation may lie in the composition of the plasma, here a pair plasma. As stated in the introduction, such an assumption allows us to consider only one parallel and one perpendicular temperature. Yet, our formalism being eventually macroscopic, we expect our conclusions to hold for electron/ion plasmas as well. Such a conjecture is currently being tested through PIC simulations of such plasmas (Shalaby Reference Shalaby2024).

Acknowledgements

Thanks are due to A. Velikovich and F. Fraschetti for valuable inputs.

Editor Luís O. Silva thanks the referees for their advice in evaluating this article.

Funding

A.B. acknowledges support by the Ministerio de Economía y Competitividad of Spain (Grant No. PID2021-125550OB-I00). R.N. acknowledges support from the NSF Grant No. AST-1816420. C.C.H. was partially supported by NSF FDSS grant AGS-1936393 as well as NASA grants 80NSSC20K1273 and 80NSSC23K0099.

Simulations were performed on computational resources provided by TACC's Stampede2 and Purdue's ANVIL through ACCESS (formally XSEDE) allocation TG-AST180008. Some of the work was supported by the Geospace Environment Modeling Focus Group “Particle Heating and Thermalization in Collisionless Shocks in the Magnetospheric MultiScale Mission (MMS) Era”.

R.N. and A.B. thank the Black Hole Initiative (BHI) at Harvard University for support. The BHI is funded by grants from the John Templeton Foundation and the Gordon and Betty Moore Foundation.

Declaration of interest

The authors report no conflict of interest.

Appendix A. Isotropic MHD conservation equations for an oblique shock

The MHD conservation equations for an oblique shock and a fluid of adiabatic index $\gamma =5/3$ read (Kulsrud Reference Kulsrud2005)

(A1)\begin{gather} n_2 v_2 \cos \xi _2 = n_1 v_1, \end{gather}
(A2)\begin{gather}B_2 \cos \theta_2 = B_1 \cos \theta_1 , \end{gather}
(A3)\begin{gather}B_2 v_2 \sin \theta_2 \cos \xi_2 - B_2 v_2 \cos \theta_2 \sin \xi_2 = B_1 v_1 \sin \theta_1 , \end{gather}
(A4)\begin{gather}\frac{B_2^2 \sin ^2\theta_2}{8 {\rm \pi}}+n_2 k_B T_2 + m n_2 v_2^2 \cos ^2\xi_2 = \frac{B_1^2 \sin ^2\theta_1}{8 {\rm \pi}}+n_2 k_B T_1 + m n_1 v_1^2 , \end{gather}
(A5)\begin{gather}- \frac{B_2^2}{4 {\rm \pi}}\sin \theta_2 \cos \theta_2 + m n_2 v_2^2 \sin \xi_2 \cos \xi_2 ={-} \frac{B_1^2 \sin \theta_1 \cos \theta_1}{4 {\rm \pi}} , \end{gather}
(A6)\begin{gather}\mathcal{A} v_2 \sin \xi_2 + \mathcal{C} = m n_1 v_1 \left(\frac{5 k_B T_1 }{2 m}+\frac{B_1^2 \sin ^2\theta_1}{4 {\rm \pi}m n_1}+\frac{v_1^2}{2}\right), \end{gather}

where $m$ is the mass of the particles, $k_B$ the Boltzmann constant and

(A7)\begin{gather} \mathcal{A} ={-} \frac{B_2^2}{4 {\rm \pi}}\sin \theta_2 \cos \theta_2, \end{gather}
(A8)\begin{gather}\mathcal{C} = m n_2 v_2 \cos \xi_2 \left( 5\frac{k_BT_2}{2 m} +\frac{B_2^2 \sin ^2\theta_2}{4 {\rm \pi}m n_2}+\frac{v_2^2}{2}\right). \end{gather}

After some algebra, $\bar {T}_2 \equiv \tan \theta _2$ is found as the solution ofFootnote 10

(A9)\begin{equation} \sum_{k=0}^4 a_k \bar{T}_2^k = 0, \end{equation}

with

(A10)\begin{equation} \left.\begin{gathered} a_0 = 2 M_{A1}^2 \left(\sin 2\theta_1-2 M_{A1}^2 \tan \theta_1\right)^2, \\ a_1={-}\dfrac{1}{8} M_{A1}^2 \tan \theta_1 [{-}4 M_{A1}^2 (31 \cos 2\theta_1+21)+80 M_{A1}^4+4 \cos ^2\theta_1 (21 \cos 2\theta_1+11)], \\ a_2 =\dfrac{1}{2} M_{A1}^2 \cos ^2\theta_1 (15 \cos 2\theta_1+1)-M_{A1}^4 (7 \cos 2\theta_1+3)+2 M_{A1}^6, \\ a_3={-}\dfrac{1}{4} \sin 2\theta_1 M_{A1}^2 ({-}2 M_{A1}^2+\cos 2\theta_1+1), \\ a_4= 3 \cos ^4\theta_1 M_{A1}^2 . \end{gathered}\right\} \end{equation}

Appendix B. Anisotropic MHD conservation equations for an oblique shock

The conservation equations for anisotropic temperatures were derived in Hudson (Reference Hudson1970) and Erkaev et al. (Reference Erkaev, Vogl and Biernat2000). They have been re-derived in Bret & Narayan (Reference Bret and Narayan2022a) with the present notations. They are formally valid even for anisotropic upstream temperatures, with $T_{\parallel 1} \neq T_{\perp 1}$. Writing them for $T_{\parallel 1}=T_{\perp 1}\equiv T_1$, they read

(B1)\begin{gather} n_2 v_2 \cos \xi _2 = n_1 v_1, \end{gather}
(B2)\begin{gather}B_2 \cos \theta_2 = B_1 \cos \theta_1 , \end{gather}
(B3)\begin{gather}B_2 v_2 \sin \theta_2 \cos \xi_2 - B_2 v_2 \cos \theta_2 \sin \xi_2 = B_1 v_1 \sin \theta_1 , \end{gather}
(B4)\begin{gather}\frac{B_2^2 \sin ^2\theta_2}{8 {\rm \pi}}+n_2 k_B ( T_{{\parallel} 2} \cos ^2\theta_2 + T_{{\perp} 2} \sin ^2\theta_2) + m n_2 v_2^2 \cos ^2\xi_2 = \frac{B_1^2 \sin ^2\theta_1}{8 {\rm \pi}}+n_2 k_B T_1 + m n_1 v_1^2 , \end{gather}
(B5)\begin{gather}\mathcal{A} + m n_2 v_2^2 \sin \xi_2 \cos \xi_2 ={-} \frac{B_1^2 \sin \theta_1 \cos \theta_1}{4 {\rm \pi}} , \end{gather}
(B6)\begin{gather}\mathcal{A} v_2 \sin \xi_2 + \mathcal{B} + \mathcal{C} = m n_1 v_1 \left(\frac{5 k_B T_1 }{2 m}+\frac{B_1^2 \sin ^2\theta_1}{4 {\rm \pi}m n_1}+\frac{v_1^2}{2}\right), \end{gather}

where

(B7)\begin{gather} \mathcal{A} = \sin \theta_2 \cos \theta_2 n_2 k_B (T_{{\parallel} 2}-T_{{\perp} 2}) - \frac{B_2^2}{4 {\rm \pi}}\sin \theta_2 \cos \theta_2, \end{gather}
(B8)\begin{gather}\mathcal{B} = v_2 \cos ^2\theta_2 \cos \xi_2 n_2 k_B ( T_{{\parallel} 2} - T_{{\perp} 2} ), \end{gather}
(B9)\begin{gather}\mathcal{C} = m n_2 v_2 \cos \xi_2 \left( \frac{k_B}{2 m}(T_{{\parallel} 2} + 4 T_{{\perp} 2} ) +\frac{B_2^2 \sin ^2\theta_2}{4 {\rm \pi}m n_2}+\frac{v_2^2}{2}\right). \end{gather}

It can be checked that setting $T_{\parallel 2} = T_{\perp 2} = T_2$ gives back the MHD equations, ((A1)–(A6)).

Appendix C. Main quantities for stages 1 and 2

Implementing the algorithm described in figure 7 requires computing the density ratio $r$, the anisotropy $A_2$ and the downstream Alfvén Mach number $M_{A2x}$, for stage 1, stage 2-firehose and stage 2-mirror. The results are presented below.

C.1 Stage 1

Solving the system of ((B1)–(B6)) with prescription (2.3) allows us to derive a polynomial for the quantify $T_2=\tan \theta _2$ defined in (3.1). It has been derived in Paper 1.Footnote 11 Using now the expression (4.1) for the Alfvén velocity, the Alfvén Mach number for stage 1 reads

(C1)\begin{equation} M_{A2x}^2 = \frac{4 (T_2^4+2) \sec^2\theta_1 M_{A1}^2}{(T_2^2-2) \sec^2\theta_1 (4 M_{A1}^2 (r-1)-r \cos 2 \theta_1+r)+2 r (T_2^4+2 T_2^2+4)}. \end{equation}

From (2.3), the anisotropy of stage 1 is simply

(C2)\begin{equation} A_2 = \tfrac{1}{2}\tan^2\theta_2. \end{equation}

The density ratio $r$ is the same as in Paper 1, namely,

(C3)\begin{equation} r=\frac{4 \mathcal{M}_{A1}^2 T_2^3 (1+T_2^2)}{\displaystyle \sum_{k=0}^5b_k T2^k}, \end{equation}

where

(C4)\begin{equation} \left.\begin{gathered} b_0 = 8 \mathcal{M}_{A1}^2 \tan \theta_1-4 \sin 2 \theta_1, \\ b_1 ={-}8 \mathcal{M}_{A1}^2+6 \cos 2 \theta_1+2,\\ b_2 = 0, \\ b_3 = 4 \mathcal{M}_{A1}^2+\cos 2 \theta_1+3, \\ b_4 = 4 \mathcal{M}_{A1}^2 \tan \theta_1-2 \sin 2 \theta_1, \\ b_5 = 2 \cos^2\theta_1. \end{gathered}\right\} \end{equation}

C.2 Stage 2 firehose

If stage 1 has $A_2 < 1-2/\beta _{\parallel 2}$, then stage 2-firehose is the end state. Since the stability criterion differs from that used in Paper 1 by the factor 2, the properties of stage 2 firehose change with respect to Paper 1. A polynomial equation can still be derived for $T_2=\tan \theta _2$ as

(C5)\begin{equation} \sum_{k=0}^3 a_k T_2^k = 0, \end{equation}

with

(C6)\begin{gather} a_0 ={-}8 (\sin 2 \theta_1-2 M_{A1}^2 \tan \theta_1)^2 , \end{gather}
(C7)\begin{gather}a_1 = ({-}2 M_{A1}^2+\cos 2 \theta_1+1) ({-}20 M_{A1}^2+9 \cos 2 \theta_1-1) \tan \theta_1, \end{gather}
(C8)\begin{gather}a_2 = 8 \cos 2 \theta_1 M_{A1}^2-8 (M_{A1}^4+M_{A1}^2)-\cos 4 \theta_1+1, \end{gather}
(C9)\begin{gather}a_3 = (2 M_{A1}^2-\cos 2 \theta_1-1) \sin 2 \theta_1 . \end{gather}

The density jump is then given by

(C10)\begin{equation} r=\frac{M_{A1}^2 T_2}{M_{A1}^2 \tan \theta_1-\sin \theta_1 \cos \theta_1}. \end{equation}

As for the Alfvén Mach number, we found in § 4 that the Alfvén speed vanishes on the firehose instability threshold. Therefore, stage 2-firehose has

(C11)\begin{equation} M_{A2x} ={+}\infty. \end{equation}

As a consequence, when the system eventually settles in stage 2-firehose with $\sigma _x < 1$ (i.e. $M_{A1x}>1$), it is evolutionary.

The anisotropy of stage 2-firehose is also modified with respect to Paper 1, due to the modified stability criterion. It now reads

(C12)\begin{equation} A_2 = 1-\frac{2 (T_2^3+T_2) \cos ^2\theta_1}{2 M_{A1}^2 (T_2-\tan \theta_1)+T_2^3 \cos ^2\theta_1+T_2 \sin ^2\theta_1+\sin 2 \theta_1}. \end{equation}

C.3 Stage 2 mirror

Since the stability threshold for the mirror instability is the same as in Paper 1, the polynomial for $T_2$, the anisotropy $A_2$ and the density ratio are also the same. The Alfvén Mach number reads here

(C13)\begin{equation} M_{A2x}^2 = \frac{2}{3}\frac{\sec^2\theta_1}{r}M_{A1}^2 . \end{equation}

Footnotes

1 For a parallel shock, the MHD conservation equations are identical to the fluid ones (Kulsrud Reference Kulsrud2005). With anisotropic temperature, they are obtained setting $\theta _1=\theta _2=\xi _2=0$ in ((B1)–(B6)) of Appendix B.

2 The nature of this instability will be specified shortly.

3 The jump of a strong fluid shock with adiabatic index $\gamma =5/3$ is 4. And for the parallel case, the MHD and fluid conservation equations are identical (Kulsrud Reference Kulsrud2005).

4 This expression already includes the $\cos \theta$ factor of (2.6). There is no need to multiply it by an additional $\cos \theta$ when computing the Alfvén Mach number (2.6).

5 This parameter is just twice the $S_\parallel$ parameter defined by (4.2). We could use only $\beta _\parallel$. But for clarity, we present the result of Abraham-Shrauner (Reference Abraham-Shrauner1967) with the notations of Abraham-Shrauner (Reference Abraham-Shrauner1967).

6 Figure 5(b) of Bret & Narayan (Reference Bret and Narayan2022a) shows the switch-on solutions of our model with $M_{A2x} \neq 1$. This is because correction (4.1) to the Alfvén speed was not considered.

7 Simulation parameters identical to those of Haggerty et al. (Reference Haggerty, Bret and Caprioli2022).

8 Switch-on shocks have been produced in the laboratory (Craig & Paul Reference Craig and Paul1973), but within a collisional environment, where MHD rules. We here deal with collisionless plasmas.

9 See § 6.3 in Falle & Komissarov (Reference Falle and Komissarov2001).

10 This quantity is written here with a bar to avoid confusion with the downstream temperature $T_2$ in ((A4), (A6)). Such a confusion is excluded in the rest of the paper since the downstream temperature splits into 2 different quantities, namely $T_{\parallel 2}$ and $T_{\perp 2}$.

11 See (4.2)–(4.5) in Paper 1.

References

Abraham-Shrauner, B. 1967 Propagation of hydromagnetic waves through an anisotropic plasma. J. Plasma Phys. 1 (3), 361378.CrossRefGoogle Scholar
Axford, W.I., Leer, E. & Skadron, G. 1977 The acceleration of cosmic rays by shock waves. In International Cosmic Ray Conference, vol. 11, p. 132.Google Scholar
Bale, S.D., Kasper, J.C., Howes, G.G., Quataert, E., Salem, C. & Sundkvist, D. 2009 Magnetic fluctuation power near proton temperature anisotropy instability thresholds in the solar wind. Phys. Rev. Lett. 103, 211101.CrossRefGoogle ScholarPubMed
Bale, S.D., Mozer, F.S. & Horbury, T.S. 2003 Density-transition scale at quasiperpendicular collisionless shocks. Phys. Rev. Lett. 91, 265004.CrossRefGoogle ScholarPubMed
Balogh, A. & Treumann, R.A. 2013 Physics of Collisionless Shocks: Space Plasma Shock Waves. Springer.CrossRefGoogle Scholar
Bell, A.R. 1978 a The acceleration of cosmic rays in shock fronts. I. Mon. Not. R. Astron. Soc. 182, 147.CrossRefGoogle Scholar
Bell, A.R. 1978 b The acceleration of cosmic rays in shock fronts. II. Mon. Not. R. Astron. Soc. 182, 443.CrossRefGoogle Scholar
Blandford, R & Eichler, D 1987 Particle acceleration at astrophysical shocks: a theory of cosmic ray origin. Phys. Rep. 154, 1.CrossRefGoogle Scholar
Blandford, R.D. & Ostriker, J.P. 1978 Particle acceleration by astrophysical shocks. Astrophys. J. Lett. 221, L29L32.CrossRefGoogle Scholar
Bret, A. 2020 Can we trust MHD jump conditions for collisionless shocks? Astrophys. J. 900 (2), 111.CrossRefGoogle Scholar
Bret, A. 2021 Modified jump conditions for parallel collisionless shocks. Phys. Plasmas 28 (8), 082107.CrossRefGoogle Scholar
Bret, A. 2023 a Density jump as a function of magnetic field strength for parallel collisionless shocks with anisotropic upstream pressure. Mon. Not. R. Astron. Soc. 520 (4), 60836090.CrossRefGoogle Scholar
Bret, A. 2023 b Density jump as a function of magnetic field strength for perpendicular collisionless shocks with anisotropic upstream pressure. Mon. Not. R. Astron. Soc. 524 (3), 44984505.CrossRefGoogle Scholar
Bret, A. & Narayan, R. 2018 Density jump as a function of magnetic field strength for parallel collisionless shocks in pair plasmas. J. Plasma Phys. 84 (6), 905840604.CrossRefGoogle Scholar
Bret, A. & Narayan, R. 2019 Density jump as a function of magnetic field for collisionless shocks in pair plasmas: the perpendicular case. Phys. Plasmas 26 (6), 062108.CrossRefGoogle Scholar
Bret, A. & Narayan, R. 2022 a Density jump as a function of magnetic field for switch-on collisionless shocks in pair plasmas. J. Plasma Phys. 88 (3), 905880320.CrossRefGoogle Scholar
Bret, A. & Narayan, R. 2022 b Density jump for oblique collisionless shocks in pair plasmas: allowed solutions. J. Plasma Phys. 88 (6), 905880615.CrossRefGoogle Scholar
Bret, A., Pe'er, A., Sironi, L., Sa̧dowski, A. & Narayan, R. 2017 Kinetic inhibition of magnetohydrodynamics shocks in the vicinity of a parallel magnetic field. J. Plasma Phys. 83 (2), 715830201.CrossRefGoogle Scholar
Brown, C.R., Juno, J., Howes, G.G., Haggerty, C.C. & Constantinou, S. 2023 Isolation and phase-space energization analysis of the instabilities in collisionless shocks. J. Plasma Phys. 89 (3), 905890308.CrossRefGoogle Scholar
Buneman, O. 1993 Computer Space Plasma Physics (ed. H. Matsumoto & Y. Omura), p. 67. Terra Scientific.Google Scholar
Caprioli, D., Haggerty, C.C. & Blasi, P. 2020 Kinetic simulations of cosmic-ray-modified shocks. II. Particle spectra. Astrophys. J. 905 (1), 2.CrossRefGoogle Scholar
Craig, A.D. & Paul, J.W.M. 1973 Observation of ‘switch-on’ shocks in a magnetized plasma. J. Plasma Phys. 9 (2), 161186.CrossRefGoogle Scholar
Drury, L.O.C. 1983 An introduction to the theory of diffusive shock acceleration of energetic particles in tenuous plasmas. Rep. Prog. Phys. 46, 973.CrossRefGoogle Scholar
Erkaev, N.V., Vogl, D.F. & Biernat, H.K. 2000 Solution for jump conditions at fast shocks in an anisotropic magnetized plasma. J. Plasma Phys. 64, 561578.CrossRefGoogle Scholar
Falle, S.A.E.G. & Komissarov, S.S. 1997 On the existence of intermediate shocks. In Computational Astrophysics; 12th Kingston Meeting on Theoretical Astrophysics (ed. D.A. Clarke & M.J. West), Astronomical Society of the Pacific Conference Series, vol. 12, p. 66.Google Scholar
Falle, S.A.E.G. & Komissarov, S.S. 2001 On the inadmissibility of non-evolutionary shocks. J. Plasma Phys. 65 (1), 2958.CrossRefGoogle Scholar
Farris, M.H., Russell, C.T., Fitzenreiter, R.J. & Ogilvie, K.W. 1994 The subcritical, quasi-parallel, switch-on shock. Geophys. Res. Lett. 21 (9), 837840.CrossRefGoogle Scholar
Feng, H.Q., Lin, C.C., Chao, J.K., Wu, D.J., Lyu, L.H. & Lee, L.C. 2009 Observations of an interplanetary switch-on shock driven by a magnetic cloud. Geophys. Res. Lett. 36.CrossRefGoogle Scholar
Gary, S.P. 1993 Theory of Space Plasma Microinstabilities. Cambridge University Press.CrossRefGoogle Scholar
Gerbig, D. & Schlickeiser, R. 2011 Jump conditions for relativistic magnetohydrodynamic shocks in a gyrotropic plasma. Astrophys. J. 733 (1), 32.CrossRefGoogle Scholar
Goedbloed, J.P., Keppens, R. & Poedts, S. 2010 Advanced Magnetohydrodynamics: With Applications to Laboratory and Astrophysical Plasmas. Cambridge University Press.CrossRefGoogle Scholar
Guo, X., Sironi, L. & Narayan, R. 2017 Electron heating in low-Mach-number perpendicular shocks. I. Heating mechanism. Astrophys. J. 851, 134.CrossRefGoogle Scholar
Guo, X., Sironi, L. & Narayan, R. 2018 Electron heating in low Mach number perpendicular shocks. II. Dependence on the pre-shock conditions. Astrophys. J. 858, 95.CrossRefGoogle Scholar
Gurnett, D.A. & Bhattacharjee, A. 2005 Introduction to Plasma Physics: With Space and Laboratory Applications. Cambridge University Press.CrossRefGoogle Scholar
Haggerty, C.C., Bret, A. & Caprioli, D. 2022 Kinetic simulations of strongly magnetized parallel shocks: deviations from MHD jump conditions. Mon. Not. R. Astron. Soc. 509 (2), 20842090.CrossRefGoogle Scholar
Haggerty, C.C. & Caprioli, D. 2020 Kinetic simulations of cosmic-ray-modified shocks. I. Hydrodynamics. Astrophys. J. 905 (1), 1.CrossRefGoogle Scholar
Hasegawa, A. 1975 Plasma Instabilities and Nonlinear Effects, Springer Series on Physics Chemistry Space, vol. 8. Springer.CrossRefGoogle Scholar
Hellinger, P., Trávníček, P., Kasper, J.C. & Lazarus, A. J. 2006 Solar wind proton temperature anisotropy: linear theory and WIND/SWE observations. Geophys. Res. Lett. 33 (9), L09101.CrossRefGoogle Scholar
Hudson, P.D. 1970 Discontinuities in an anisotropic plasma and their identification in the solar wind. Planet. Space Sci. 18 (11), 16111622.CrossRefGoogle Scholar
Karimabadi, H., Krauss-Varban, D. & Omidi, N. 1995 Temperature anisotropy effects and the generation of anomalous slow shocks. Geophys. Res. Lett. 22 (20), 26892692.CrossRefGoogle Scholar
Kennel, C.F., Blandford, R.D. & Wu, C.C. 1990 Structure and evolution of small-amplitude intermediate shock waves. Phys. Fluids B 2 (2), 253269.CrossRefGoogle Scholar
Kulsrud, R. M 2005 Plasma Physics for Astrophysics. Princeton University Press.CrossRefGoogle Scholar
Landau, L.D. & Lifshitz, E.M. 1981 Course of Theoretical Physics, Physical Kinetics, vol. 10. Elsevier.Google Scholar
Maruca, B.A., Kasper, J.C. & Bale, S.D. 2011 What are the relative roles of heating and cooling in generating solar wind temperature anisotropies? Phys. Rev. Lett. 107, 201101.CrossRefGoogle ScholarPubMed
Niemiec, J., Pohl, M., Bret, A. & Wieland, V. 2012 Nonrelativistic parallel shocks in unmagnetized and weakly magnetized plasmas. Astrophys. J. 759 (1), 73.CrossRefGoogle Scholar
Russell, C.T. & Farris, M.H. 1995 Ultra low frequency waves at the earth's bow shock. Adv. Space Res. 15 (8–9), 285296.CrossRefGoogle Scholar
Sa̧dowski, A., Narayan, R., McKinney, J.C. & Tchekhovskoy, A. 2014 Numerical simulations of super-critical black hole accretion flows in general relativity. Mon. Not. R. Astron. Soc. 439 (1), 503520.CrossRefGoogle Scholar
Sa̧dowski, A., Narayan, R., Tchekhovskoy, A. & Zhu, Y. 2013 Semi-implicit scheme for treating radiation under M1 closure in general relativistic conservative fluid dynamics codes. Mon. Not. R. Astron. Soc. 429 (4), 35333550.CrossRefGoogle Scholar
Sagdeev, R.Z. 1966 Cooperative phenomena and shock waves in collisionless plasmas. Rev. Plasma Phys. 4, 23.Google Scholar
Salas, M.D. 2007 The curious events leading to the theory of shock waves. Shock Waves 16 (6), 477487.CrossRefGoogle Scholar
Schwartz, S.J., Henley, E., Mitchell, J. & Krasnoselskikh, V. 2011 Electron temperature gradient scale at collisionless shocks. Phys. Rev. Lett. 107, 215002.CrossRefGoogle ScholarPubMed
Shalaby, M. 2024 In preparation.Google Scholar
Silva, T., Afeyan, B. & Silva, L.O. 2021 Weibel instability beyond bi-Maxwellian anisotropy. Phys. Rev. E 104 (3), 035201.CrossRefGoogle ScholarPubMed
Sironi, L., Spitkovsky, A. & Arons, J. 2013 The maximum energy of accelerated particles in relativistic collisionless shocks. Astrophys. J. 771 (1), 54.CrossRefGoogle Scholar
Spitkovsky, A. 2005 Simulations of relativistic collisionless shocks: shock structure and particle acceleration. In Astrophysical Sources of High Energy Particles and Radiation (ed. T. Bulik, B. Rudak & G. Madejski), American Institute of Physics Conference Series, vol. 801, pp. 345–350.Google Scholar
Thorne, K.S. & Blandford, R.D. 2017 Modern Classical Physics: Optics, Fluids, Plasmas, Elasticity, Relativity, and Statistical Physics. Princeton University Press.Google Scholar
Weibel, E.S. 1959 Spontaneously growing transverse waves in a plasma due to an anisotropic velocity distribution. Phys. Rev. Lett. 2, 83.CrossRefGoogle Scholar
Zeković, V. 2019 Resonant micro-instabilities at quasi-parallel collisionless shocks: cause or consequence of shock (re)formation. Phys. Plasmas 26 (3), 032106.CrossRefGoogle Scholar
Zel'dovich, Y.B. & Raizer, Y.P. 2002 Physics of Shock Waves and High-Temperature Hydrodynamic Phenomena. Dover Publications.Google Scholar
Figure 0

Figure 1. System under scrutiny. It is identical to that of Paper 1. Although the conservation equations used are valid for any upstream temperatures, we here, as in Paper 1, only treat the strong sonic case $T_{\perp 1}=T_{\parallel 1}=0$. We work in the frame of reference where $\boldsymbol {v}_1$ is normal to the front.

Figure 1

Figure 2. Upper row: density jump $r$ for different angles $\theta _1$ in terms of $M_{A1x}=M_{A1}/\cos \theta _1$, defined by (2.6). Sometimes there is more than 1 solution. Second row: angle $\theta _2$ of the field $\textbf {B}_2$ with the $x$ axis (see figure 1). Third row: evolutionarity criterion. Some branches, namely the ones crossing the shaded areas, are to be excluded. Lower row: same as upper row, but without the branches excluded by the evolutionarity criterion. There is now but 1 solution for a given $M_{A1x}$, or none. The black circles for $\theta _1=0.5{\rm \pi} /2$ are the results of our MHD simulations (see § 3.1).

Figure 2

Figure 3. Magnetohydrodynamics simulation in which the evolutionarity criterion results in the formation of 2 sub-shocks instead of one. The dashed blue lines show the initial fluid configuration at $t=0$, when two fluid columns collide at $x=0$. The solid red lines show the configuration at $t=1.5$. Each fluid column develops two shocks, a fast shock at $|x|=1.58$ and a slow shock at $|x|=0.34$. The little dip in density $n$ at $x=0$ for $t=1.5$ is a numerical artefact.

Figure 3

Figure 4. Stability diagram. The solid lines picture the criterion used here. The dashed lines pertain to other criteria discussed in § 4.

Figure 4

Figure 5. (a) All solutions offered by our model for $\theta _1=0$. (b) All the solutions fulfil the evolutionarity criterion. Note that this is different from the isotropic MHD problem, where some solutions for $\theta _1=0$ enter the shaded forbidden zone (figure 2, row 3). Note also that, with definition (4.1) of the Alfvén speed, $M_{A2x}=+\infty$ on the firehose threshold. This is why the green branch on the left is sent to $M_{A2x}=+\infty$ on the centre plot. Notably, $M_{A2x}=1$ for the switch-on shocks of our model (red segment on centre plot). For this analysis, the horizontal axis must be $\sigma _x=\sigma \cos ^2\theta _1$, which for $\theta _1=0$ makes no difference. (c) End result once the ‘$\theta _2$ closest to $\theta _1$’ filter has been applied. This is the result found in Bret & Narayan (2018) and Haggerty et al. (2022).

Figure 5

Figure 6. (a) Expected value of $B_{\perp 2}/B_1$, in red for the switch-on solutions with $\theta _2\neq 0$, and in black for the solution with $\theta _2=\theta _1=0$. The circles show the values of $\sigma$ which have been simulated. (b) Results of the PIC simulations, with the shock density profile on top, and the ratio $B_{\perp 2}/B_1$ below. The $\sigma$ value for a curve is given by the circle of the same colour on the left plot. The same results are obtained for $\theta _1=2^\circ$. The flow is along the $x$ axis and $d_e=c/\omega _p$ is the electron inertial length.

Figure 6

Figure 7. Flowchart of the resolution of our model. In case stage 1 is mirror unstable for a pair $(\sigma,\theta _1)$, there can be a degeneracy in the choice of stage 2-mirror for the same pair $(\sigma,\theta _1)$. In such a case, we choose the stage 2-mirror with the $\theta _2$ closest to the $\theta _2$ of the stage 1 it comes from, as was already done in Bret & Narayan (2022b). Such a situation never happens with stage 2-firehose.

Figure 7

Figure 8. Detailed progression of the application of the filtering process described by the flow chart in figure 7, for $\theta _1=0.3{\rm \pi} /2$. In (d), the green colour indicates that these stage 1 solutions are firehose unstable. Then in (g), the green colour shows the density jump of stage 2-firehose. The red branch in ( f) with $A_2 \sim -1$ discards stage 2 firehose solutions with $r >4$, not shown in (e).

Figure 8

Figure 9. Result of the filtering process elaborated in this work for various angles $\theta _1 \in [0, {\rm \pi}/2]$. All of them but $\theta _1=0.3{\rm \pi} /2$ are the ones considered in figure 8 of Bret & Narayan (2022b). The grey dashed lines picture the MHD result, evolutionary filtered. The green solutions pertain to (stage 2-firehose, green), and the red ones to (stage 2-mirror, red).