# Parametric study of non-relativistic electrostatic shocks and the structure of their transition layer

###### Abstract

Nonrelativistic electrostatic unmagnetized shocks are frequently observed in laboratory plasmas and they are likely to exist in astrophysical plasmas. Their maximum speed, expressed in units of the ion acoustic speed far upstream of the shock, depends only on the electron-to-ion temperature ratio if binary collisions are absent. The formation and evolution of such shocks is examined here for a wide range of shock speeds with particle-in-cell (PIC) simulations. The initial temperatures of the electrons and the 400 times heavier ions are equal. Shocks form on electron time scales at Mach numbers between 1.7 and 2.2. Shocks with Mach numbers up to 2.5 form after tens of inverse ion plasma frequencies. The density of the shock-reflected ion beam increases and the number of ions crossing the shock thus decreases with an increasing Mach number, causing a slower expansion of the downstream region in its rest frame. The interval occupied by this ion beam is on a positive potential relative to the far upstream. This potential pre-heats the electrons ahead of the shock even in the absence of beam instabilities and decouples the electron temperature in the foreshock ahead of the shock from the one in the far upstream plasma. The effective Mach number of the shock is reduced by this electron heating. This effect can potentially stabilize nonrelativistic electrostatic shocks moving as fast as supernova remnant (SNR) shocks.

###### pacs:

52.65.Rr, 52.35.Tc, 52.35.Qz, 98.38.Mz^{†}

^{†}thanks: Electronic mail: Mark.E.Dieckmann@itn.liu.se

## I Introduction

A supernova explosion accelerates a significant fraction of the material of the progenitor star up to a few percent of the light speed c. The initial density of the blast shell plasma is many orders of magnitude larger than that of the surrounding plasma of the interstellar medium (ISM); the blast shell expands freely. The reduction of the plasma density at the front of the radially expanding blast shell and its piling up of the ISM plasma imply that the latter eventually starts to affect the expansion. A forward shock between the upstream ISM plasma and the downstream plasma and a reverse shock between the downstream plasma and the (upstream) blast shell plasma form close to the blast shell front. The downstream region between both shocks expands in time. The shock formation initiates the Sedov-Taylor phase of the supernova remnant (SNR) blast shell’s expansion Truelove (), in which it can propagate over astronomical distances before it is stopped by its interaction with the ISM.

The low particle collision frequency in the ISM implies that the forward shock is collision-less. Its dynamics is dominated by collective plasma processes, which imposes constraints on the range of parameters for which stable shocks can form. The temperature of the cool ionized component of the ISM is about 1 eV and the sound speed m/s. The sound speed is defined by for electrons and ions with temperatures and (ion mass and Boltzmann constant ), assuming that the adiabatic constant is equal for both species. A fast SNR shock, like the north-eastern (NE) shock of the SNR RCW86 Helder (), would have a Mach number in the ISM plasma. The weak magnetization of the ISM Ferriere () suggests that SNR shocks are essentially unmagnetized.

Electrostatic unmagnetized plasma shocks can be generated in laboratory plasmas Honzawa (); Romagnani () and in simulation plasmas ForslundA (); ForslundB () by the collision of two plasma clouds at an appropriate speed. Collisions between identical plasma clouds can result in two types of unmagnetized electrostatic shocks ForslundB (). Sub-critical shocks can convert the entire kinetic energy of the inflowing upstream plasma into heat. The shock-reflected ion beam of super-critical shocks provides an additional energy dissipation mechanism and such shocks are stable at larger Mach numbers than the sub-critical shocks. Analytic estimates of the maximum Mach number exist, at which we find stable electrostatic shocks. The largest estimates of the maximum Mach number of shocks, which develop out of the collision of two identical plasma clouds, are derived under the assumption that the electron temperature exceeds by far the ion temperature and that the particle velocity distributions are Maxwellian’s far upstream of the shock ForslundA (); ForslundB (). The peak Mach number of subcritical shocks is in this case , while that of super-critical shocks is . The peak Mach number of shocks decreases with a decreasing ratio Bardotti (). No stable electrostatic shocks exist under these approximations above a value , which is orders of magnitude below those of SNR shocks.

An asymmetry between the colliding plasma clouds in terms of the electron temperature and density yields a double layer Langmuir () that can raise the maximum value of Sorasio (). Double layers are unstable, because the thermal pressure of the dense hot downstream plasma is not necessarily balanced by the deceleration or reflection of the upstream plasma. Experimental observations indicate that most double layers are rarefactive and some are compressive Raadu () in their rest frame. They may thus not result in stable shocks, but they can trigger the formation of stable electrostatic shocks SarriA (), even if the initial Mach number of the flow speed is 100 Dieckmann ().

Shocks can be stabilized for by a magnetic field. An example is the Earth’s quasi-perpendicular bow shock Bale (). It has a Mach number 10-12 with respect to the ion acoustic speed, if we take an electron temperature of 10 eV of the solar wind and a shock speed of 500 km/s in the rest frame of the solar wind. This Mach number exceeds the stability limit of electrostatic shocks and the magnetic field must thus provide an additional stabilization. Indeed the Earth bow shock is a fast MHD mode shock Sckopke (). Shocks that propagate into an upstream medium with a uniform magnetic field have been examined in the context of SNR shocks with particle-in-cell (PIC) simulations LembegeA (); Shimada (); LembegeB (); Lee (); Scholer (); Chapman (); Umeda (); Amano (); Matsukiyo (); Pohl (). These studies have addressed their stability, their efficiency in accelerating electrons and they have examined wave instabilities in the foreshock BretA (). Usually a relatively strong magnetic field is used due to computational constraints. It is important to determine if a magnetic field is essential for a shock formation and for its stable propagation, given that the magnetic field of the ISM is relatively weak.

Unmagnetized electrostatic shocks with a Mach number may not exist in the ISM. This Mach number is however a least upper bound. It is based on the electron temperature of the ISM and neglects any pre-heating of the ISM by the shock. Unless the shock is quasi-perpendicular, shock-accelerated electrons or hot downstream electrons can escape upstream of the shock. They increase the ion acoustic speed, which in turn decreases the Mach number. An electron foreshock is observed, for example, at the Earth’s bow shock Eastwood (). Hot electrons are detected ahead of this shock even though the temperature difference between the electrons in the (downstream) magnetosheath and the upstream solar wind is well below that between the downstream region of SNR shocks and the ISM. Hotter electrons ahead of an SNR shock result in a higher and therefore in a lower value of . This pre-heating mechanism is independent of the heating of the foreshock by the cosmic ray precursor VolkA (); Drury (); Vink () and the wave precursors Ghavamian ().

The possible existence of electrostatic unmagnetized SNR shocks motivates our numerical study. The shock is generated by letting an unmagnetized plasma collide head-on with a reflecting wall at the speed . This corresponds to a collision of two equal plasma clouds at 2, because the particles are reflected elastically. We perform 7 separate simulations in which the collision speed between the plasma and the wall is varied from 1.06 to 2.48. The temperatures of the co-moving and spatially uniform electrons and ions are equal. The downstream plasma and its boundary, the shock, expand away from the wall and into the upstream region. The shock speed is thus always larger than the collision speed.

Our results are as follows. All shocks are super-critical and move at a speed higher than 1.6 in the upstream reference frame. The simulations thus demonstrate that stable shocks can form even if the ions are as hot as the electrons, which is not possible according to some analytic models Bardotti (). A higher speed of the inflowing plasma results in a sharper shock discontinuity and in a larger density of the shock-reflected ion beam ForslundB (). The fraction of ions that can cross the shock decreases and the expansion of the downstream plasma is slowed down; doubling from 1.06 to 2.12 can only raise the shock’s Mach number from 1.69 to 2.31. Shocks form instantly only for 1.9 . An increase of beyond this value delays the shock formation by tens of inverse ion plasma frequencies. A collision at results in a shock with and no shock forms for during the resolved timescale, which is comparable to that achieved in laboratory experiments. It may thus not always be possible to observe fast electrostatic shocks in laser-plasma experiments, while they can develop in the essentially unbounded astrophysical plasmas.

The density of the shock-reflected ion beam, which increases with the collision speed, raises the overall ion number density ahead of the shock. This foreshock region goes on a positive potential relative to the upstream plasma, which is maintained by the ambipolar electrostatic field of the plasma density gradient. Upstream electrons are accelerated towards the shock by this electric field. The accelerated electrons interact with the electrons, which leak from the downstream region into the foreshock. Their mixing raises the mean thermal energy of the electrons in the foreshock to about 1/3 of the downstream one for all considered cases and the electron velocity distribution becomes a flat-top one rather than a Maxwellian one. The simulations suggest that this electron temperature depends on the relative speed between the pre-accelerated upstream electrons and the leaking downstream electrons and not on the electron temperature far upstream of the shock. The implication is that once an electrostatic shock is present it is no longer appropriate to determine its stability properties using the Mach number it would have in the far upstream, because it propagates through a much hotter foreshock plasma.

We outline the particle-in-cell simulation scheme and our initial conditions in section 2, we present our results in section 3 and we discuss them in section 4.

## Ii The simulation code and the initial conditions

### ii.1 The model

The mechanism how a shock forms can be understood with the qualitative model shown in Fig. 1. A plasma cloud is reflected at the position =0. The cloud is spatially uniform and unmagnetized. It consists of electrons and ions, each with the number density , that move at the speed modulus towards x=0. We assume in this qualitative model that the ions are cold, so that the ion front does not spread out in time. The reflected ions and the incoming ions interpenetrate and a plasma sheet with an ion density develops close to = 0. It expands at the speed to larger values of . We call this sheet the overlap layer. The higher electron mobility implies that some electrons stream out of this layer. They form a thin negatively charged sheet just outside of it. An ambipolar electrostatic field is built up and the overlap layer goes onto a positive potential with respect to the inflowing upstream plasma. The electric field decelerates the ions that enter the overlap layer and this plasma compression raises the density of this layer above the value . This field also confines the electrons and accelerates the surface ions of the overlap layer SarriA (); Sack ().

If the electric field is strong enough to equalize the speed of all ions in the overlap layer, then the ions of both plasma slabs will mix and form a single population. This downstream region is bounded by a forward shock and it expands due to ion accumulation. The upstream ions and electrons are heated up and compressed as they enter the downstream region. If all incoming ions enter the downstream region and mix with its plasma, then the shock is sub-critical. A fraction of the ions is reflected by a super-critical shock. They outrun the shock and form together with the incoming plasma the foreshock plasma. The upstream region is the one that has not yet been reached by shock-reflected ions.

We employ a 1D simulation geometry, which separates the shock physics from that of secondary instabilities and simplifies the interpretation of the results. Instabilities can be triggered by the counterstreaming ion beams in the foreshock with a relative speed that is larger than the ion acoustic speed. However, these instabilities can only develop if wave vectors are resolved that are oriented obliquely to the difference vector between the mean speeds of both ion beams ForslundC (); Karimabadi (); Kato (). Such vectors are excluded here by the alignment of the collision velocity vector with the simulation direction. The maximum shock speed is a few times the ion acoustic speed and thus much smaller than the electron thermal speed, which suppresses the Buneman instability that can drive electrostatic waves with a beam-aligned wave vector Buneman ().

### ii.2 The solved equations and the initial conditions

A particle-in-cell (PIC) simulation code Dawson () solves Maxwell’s equations together with the relativistic Lorentz force equation for an ensemble of computational particles (CPs). The collective charge distribution and current distribution on the grid is obtained from the interpolation of the charge and current of all CPs onto the grid. The collective charge and current distributions evolve in time the electromagnetic fields through the discretized forms of the Ampére’s and Faraday’s laws

(1) | |||

(2) |

They fulfill Gauss’ law and either as constraints or through correction steps. The electric and magnetic fields update in turn the momentum of each CP through the relativistic Lorentz force equation.

(3) |

This equation updates the momentum of the particle of species . The position is updated as . We use the EPOCH PIC code Brady ().

Our initial conditions are as follows. A plasma consisting of electrons and ions with mass ratio is introduced into a one-dimensional simulation box with length . We use perfectly reflecting boundary conditions at and open boundary conditions at . The spatially uniform plasma fills up the entire simulation box at . The particle number density of each species is and the initial temperature of their Maxwellian velocity distributions is . The electron thermal speed m/s and m/s. The Debye length, which includes the contribution from the ions, is , where is the electron plasma frequency. The ion plasma frequency is . The box length is resolved by 14400 cells of size . We represent the electron species and the ion species by 2000 particles per cell, respectively.

We perform 7 separate simulations with different collision speeds : (Run 1), (Run 2), (Run 3), (Run 4), (Run 5), (Run 6) and (Run 7).

## Iii Simulation results

Non-relativistic electrostatic shocks are in their most basic form one-dimensional plasma structures. We assume that they are planar. Such shocks are fully described by the distribution of the electrostatic potential along the shock normal direction and by the phase space distributions and of the electrons and ions, respectively. The shock’s electrostatic potential has to be strong enough to reflect the incoming upstream ions, so we normalize it as . The potential difference between the far upstream region and the downstream region should be , because the shock speed exceeds . Our reference potential is at the reflecting boundary. We subtract from the offset , which is the spatio-temporal average of the fully developed downstream potential. The offset differs between the simulation runs, because we obtain in some cases electrostatic ion phase space holes IonHole () at the wall, which put the wall and the downstream plasma on different potentials.

The potential reveals if and when the potential difference between the plasma overlap layer and the incoming plasma cloud becomes large enough to reflect the incoming upstream ions. In what follows, we discuss separately the results of the simulation runs 1-7 followed by a comparison of the shocks in terms of plasma compression and electron heating.

Run 1: The plasma collides at the speed with the wall and a shock should form, which is confirmed by the potential in Fig. 2(a).

A positive potential forms at the wall on electron time scales. The front of this structure expands in the simulation frame at the constant speed or for the ion-to-electron mass ratio that we use. The front of the potential thus expands at the speed in the reference frame of the incoming (upstream) plasma. Typical values of the upstream potential reach . The potential falls off smoothly in space and the width of the transition layer is of the order of .

The ion distribution in Fig. 2(b) shows that this transition layer is actually wider. The plasma has reached a phase space distribution close to , which is stable over the considered time scale. The velocity distribution is not a single Maxwellian one that is centred at . We find instead a local minimum of the phase space density at and with a value, which is about 20% below the peak one. The shock can thus not fully convert the flow energy of the incoming upstream ions into heat. The width of the transition layer of the electrons in Fig. 2(c) is even wider. It expands in time because the electrons react to the positive potential of the shock-reflected ion beam. The Buneman instability is suppressed by thermal effects and can not account for the observed electron temperature rise in the foreshock region. The electrons are heated instead by a turbulent mixing of electrons, which are accelerated by the foreshock potential towards the shock, with the electrons that leak from the downstream into the foreshock region. The downstream electrons (in the interval at ) do not have a Maxwellian distribution. The distribution shows a flat top centered at and a steep fall-off of the phase space density for . A flat-top distribution is also observed in the case of the foreshock electrons.

The non-Maxwellian electron distribution alters the shock stability condition Sharmin () and the ions do not get fully thermalized when they cross the shock. These two aspects may explain why a shock can form at the high ion temperature and Mach number we consider here, which exceeds significantly the limit determined by the analytic model in Ref. Bardotti ().

Run 2: Raising the collision speed to leaves unchanged the qualitative structure of the shock. The potential distribution in Fig. 3(a) shows that the shock forms again on electron time scales and that it expands to larger at a constant speed , which corresponds to . This shock thus expands at a lower speed in the simulation frame than the one in Run 1. However, the larger implies that its speed in the upstream frame of reference is higher with .

The ion distribution in Fig. 3(b) shows some qualitative differences compared to the one in Fig. 2(b). The shock-reflected ion beam appears to be denser, which we confirm below. We observe in Fig. 3(b) a pronounced velocity change of the incoming upstream ions at , which we do not find in this form in Fig. 2(b). This stronger velocity change should lead according to the continuity equation to an increased plasma compression. The electron phase space distribution in Fig. 3(c) shows again a flat-top velocity distribution downstream of the shock and a broad transition layer towards the upstream region.

Figure 3(a) shows a weak depletion of the downstream potential that accelerates in time. This depletion crosses the position at the time and it reaches at . This electrostatic field structure is connected to a depletion of the ion phase space density at and (Fig. 3(b)). We discuss ion phase space holes in more detail below.

Run 3: The collision speed between the plasma and the wall is set to . The maximum of the positive electrostatic potential in Fig. 4(a) still develops close to the wall. The potential jump between the overlap layer and the upstream region does however not reach the value necessary for a shock formation on electron time scales as before. The potential jump reaches the value , which is required for a shock formation, at the time . The front of the potential moves at the speed or to larger . Its Mach number in the upstream reference frame is . Figure 4(a) shows a potential structure that outruns the shock. This secondary structure reaches the position at . Its speed is comparable to , which corresponds to the speed at which the reflected ions move to larger . The potential jump from to at the shock front after , which separates the downstream region and the shock’s foot (foreshock), is not enough to reflect ions that move at the speed towards the shock. However, the potential difference between the downstream region and the far upstream region is still large enough to reflect ions.

Figure 4(b) reveals that a shock is present at . This structure reflects the incoming ions, which have the lowest speed in the reference frame of the shock, and it accelerates those ions from the downstream region, which have the largest speed in the direction of the shock. The downstream population of the ions shows again a double-peaked velocity distribution that is getting more pronounced as we move away from the wall. This separation is upheld by the potential, which changes at from the value at to at . The shock transition layer has an extension comparable to a few hundred . The electron distribution in Fig. 4(c) shows again the subdivision into downstream electrons with , foreshock electrons in the interval and upstream electrons at even larger .

The potential distribution in Fig. 4(a) reveals a localized depletion of the potential in the downstream region that is similar to the one in Fig. 3(a) but stronger. Its curved trajectory implies that it accelerates in time. It reaches the shock position at the time . This potential depression is caused by an ion phase space hole, which is demonstrated by Fig. 5.

Ion phase space holes are structures, in which ions gyrate in a localized depletion of the electrostatic potential and form a phase space vortex. This negative potential is upheld by a local depletion of the ion charge density IonHole (). These vortices are tied to the ion distribution and move at speeds comparable to the ion thermal speed . Its apparent acceleration in Fig. 4(a) arises from a change of the bulk speed of the ion beam that carries it.

Run 4: A collision speed results in a shock formation that is similar to that in Run 3. The potential in Fig. 6(a) demonstrates that initially only a potential of is reached close to the wall. A shock forms at , which propagates at the speed or to the right, yielding a Mach number in the upstream frame of reference. A shock foot is visible that expands at a speed to larger values of .

The ion distribution in Fig. 6(b) reveals a phase space vortex at , a broad spatial interval with an ion distribution that is practically uniform along and incoming and reflected ions at . The electron distribution in Fig. 6(c) is qualitatively similar to those in the previous runs.

Run 5: A simulation with shows a qualitatively different shock formation stage compared to those in the previous cases. A weak potential forms first. This structure detaches from the wall and its peak potential is located at at . This potential starts to grow and reaches at . The potential decreases to at . We observe a cyclic reformation of the potential with a period of . Despite this intermittent behaviour, the front of the potential propagates at a uniform speed or into the upstream plasma. The front thus moves upstream at .

The particle distributions at in Figs. 7(b,c) agree with those observed in the previous two runs, implying that the shock structure does not strongly depend on the details of its formation. An ion phase space hole has formed at the wall and the electrons have a flat top distribution downstream of the shock and in the foreshock.

The cause of the initial oscillations of the potential is a non-stationary ion phase space hole at the reflecting boundary. Figure 8 compares the ion phase space distributions close to the wall at the times and .

The ion phase space hole and thus the minimum of the potential depression are centered at at . This implies that the plasma overlap layer is on a higher potential than the wall. The ion phase space hole’s center and the minimum of the electrostatic potential have moved to at the later time. The difference of the electrostatic potential between the plasma overlap layer and the wall has been reduced.

Run 6: A collision speed is the largest one that results in a shock during a simulation time . The shock does not form instantly, in line with the results of the simulation runs 3-5. It forms close to the wall and the potential distribution in Fig. 9(a) shows no intermittent behaviour. The latter is thus not a consequence of large values of . A strong and steady foreshock potential is observed, which is again tied to the shock-reflected ion beam. We notice that a straight line fit of the shock front does not intersect at and the shock has thus not formed at the front of the plasma cloud.

The ion distribution in Fig. 9(b) agrees qualitatively with those in the simulation runs 3-5 with respect to the subdivision into the ion phase space hole at the wall, a downstream region with a double-peaked velocity distribution and the foreshock structure consisting of incoming and reflected ion beams.

An increase of the collision speed to does not result in a shock formation during , which is demonstrated by the corresponding displayed in Fig. 10.

The potential difference between the plasma overlap layer and the far upstream region is not sufficient to yield a shock; no localized strong jump of the electrostatic potential develops that could be associated with a shock. The ion phase space distribution (not shown) does not show a shock either.

A comparison of the ion distributions in the simulation runs 1-6 shows two trends. The first trend is a shrinking of the ionic shock transition layer that connects the downstream region with its spatially uniform double-peaked velocity distribution and the upstream region that consists of incoming and reflected ions. The width of this interval, which is characterized by maxima in the ion velocity distribution that diverge with increasing values of , is in Fig. 2(b) and it shrinks with increasing to reach a few tens of for the highest collision speed. The second trend, which has already been reported in Ref. ForslundB (), is that the density of the shock-reflected ion beam increases. The velocity distribution of the incoming ions and the reflected ions is practically symmetric with respect to for in Fig. 9(b). The shock thus reflects practically all incoming ions and the low number density of ions that make it downstream explains the low expansion speed of or in Fig. 9(a). Its Mach number is .

A comparison of the electron phase space distributions obtained from the simulation runs 1-6 also reveals a trend. The distributions show a striking similarity if the velocity axis is scaled by . An important feature is that the thermal spread of the electrons in the foreshock is practically unchanged by the choice of . We attribute this to the way they are heated up. Figures 2(c)-4(c) and Figs. 6(c), 7(c), 9(c) show that the electron’s phase space distribution is bounded by two populations. The electrons with the largest positive values of are supplied by downstream electrons that leak through the shock. We can determine the electron population that can leak into the foreshock as follows. The shock potential in the case of Run 6 can reflect ions that move from the upstream region towards the shock at the speed . The potential can thus confine downstream electrons with a velocity . The fastest downstream electrons in Fig. 9(c) have a speed modulus and they can not be confined by the shock potential. They have lost a significant fraction of their kinetic energy by the time they have travelled far upstream but they are still fast in the foreshock region. The foreshock electrons with the largest negative speeds are supplied by the upstream electrons, which are accelerated towards the shock by the positive potential of the foreshock. Both counter-streaming electron populations thermalize in the foreshock.

In what follows we want to provide further support for the trends outlined above. We analyse for this purpose the normalized ion density and the mean kinetic energy per electron with . We normalize by the mean kinetic energy per electron of a Maxwellian velocity distribution with the temperature .

We compare the ion density distributions in the Figs. 11(a) for run 1-3 and 12(a) for run 4-6. The corresponding mean thermal energies of the electrons are compared in the Figs. 11(b) and 12(b). All curves are measured at the time . We capture the full spatial interval of width affected by the shock-reflected ions.

The ion densities downstream and, thus, the plasma compression increase with an increasing value of . This is a well-known fact that is here confirmed by the simulations. A somewhat surprising observation is that the downstream density in run 1 is only , which corresponds to a mere superposition of the density contributions of the incoming and reflected ions. The electrostatic potential in Fig. 2(a) was not strong enough to slow down significantly the incoming ions (See Fig. 2(b)) and the ions are thus not strongly compressed.

The curves for demonstrate that the width of the shock transition layer shrinks with an increasing value of . The distribution for run 1 changes from the downstream value to the foreshock value over several hundreds of , while the ion density is practically discontinuous on the displayed spatial scale in run 6. The ion density in the shock’s foot increases with . It is immediately ahead of the shock in run 1 and it gradually increases up to in run 6. The latter value implies that most ions are reflected by the shock; the shock acts as a piston ForslundB ().

The curves for the mean thermal energy per electron follow the trend of the ion densities. The mean electron thermal energy is only slightly elevated in the foreshock region in run 1. The foreshock electrons get hotter with increasing and the typical thermal energy of electrons in the foreshock has more than doubled in run 6 compared to the value far upstream. The shock in run 6 thus moves through a plasma with an elevated sound speed since . This temperature rise is not large in our simulations. It can nevertheless be important because it is caused by the mixing of leaking downstream electrons that move upstream and upstream electrons that have been accelerated towards the shock by the foreshock potential. The energy available to the incoming upstream electrons is dominated by the bulk kinetic energy they gain in the foreshock potential and not by their thermal energy. The foreshock’s electron temperature should thus be robust against changes of the electron temperature far upstream. The mean thermal energy of the electrons converges with increasing values of to the initial value and the rise in the electron temperature is thus limited to the foreshock region.

The mean thermal energies of the foreshock electrons and the downstream electrons increase with . We can compute from Figs. 11(b) and 12(b) the ratio of the mean electron thermal energies in the foreshock and downstream. We compare for this purpose the respective values at the positions and in Fig. 11(b) and those at and in Fig. 12(b). The measured ratios are displayed in Fig. 13 as a function of the shock’s Mach number.

We find that the temperature of the foreshock electrons is in all cases about 1/3 of the value in the downstream region.

## Iv Discussion

We have performed here a series of one-dimensional particle-in-cell (PIC) simulations of electrostatic shocks. The purpose has been to gain insight into the conditions under which such shocks form and into the structure of their transition layer for a wide range of shock speeds. These aspects are important in the context of supernova remnant (SNR) shocks. Their large nonrelativistic speeds exceed the ion acoustic speed in the ISM by a factor . This factor is much larger than the limit for stable shocks, which is derived from various theoretical models ForslundA (); ForslundB (); Bardotti (). It is thus essential to know the degree to which the ISM’s temperature is important for SNR shocks. In other words, can a shock, once it is present, generate the electron temperatures it needs to be stable? This parametric study is also relevant for laboratory experiments. Important questions are what kind of shocks can form during the limited time over which observations are possible and what properties these shocks have.

Our study addressed shocks in initially unmagnetized plasma that develop out of the collision of two equal plasma clouds at a speed that is less than 1% of the light speed . Such shocks remain electrostatic because the low flow speeds yield only weak currents and magnetic fields. A speed in excess of results in the growth of magnetic fields Kazimura (); Kato (); Pohl () that are strong enough to affect the shock dynamics and they become dominant at ultrarelativistic speeds even if the plasma was initially unmagnetized Spitkovsky (). The chosen collision speeds are realistic for the late expansion phase of SNR blast shells. Our parametric study has provided the following results.

Shocks form for the considered electron and ion temperatures up to a Mach number 2.5. A shock with this Mach number reflects practically all incoming ions and acts as a piston. Raising the plasma collision speed by another 8% results in a shock-less interpenetration of the plasma clouds. Slower shocks with a Mach number form on electron time scales while the faster ones develop slower, typically during about 10-100 inverse ion plasma frequencies. The instability, which results in the formation of shocks, does not set in immediately after a collision of plasma clouds at a high speed and the formation time of the shock may be unpredictable. We can not exclude an eventual formation of a shock in our simulation that used the fastest collision speed. One possibility is that the ion acoustic instability, which can develop in multidimensional systems, pre-heats the plasma. The Mach number of the flow is reduced by the increasing ion acoustic speed and a stable shock may eventually form. Another shock formation mechanism is the destabilization of nonlinear plasma structures such as ion phase space holes that can trigger the formation of shocks even in very fast flows Dieckmann ().

The shock-reflected ion beam results in electron heating in all simulations. This heating is not caused by instabilities such as the Buneman instability that is excluded here by the low relative flow speed between electrons and ions. It arises from the electron acceleration in the ambipolar electric field of the plasma density gradient. The foreshock electron temperature in our simulations has been about 1/3 of the electron temperature in the downstream plasma and it is probably independent of the upstream temperature. An electrostatic shock moves in this case through a medium with an ion sound speed that depends on the electron temperature in the post-shock plasma and on the density of the shock-reflected ion beam. The shock-reflected ion beam, which is essential for the electron pre-heating, implies that the shock must be super-critical.

The difference between the electron temperature in the foreshock and in the upstream is moderate in our simulations and at solar system shocks. The difference can become significant at astrophysical shocks, where the post-shock electron temperature can be of the order of keV while that in the upstream is about 1 eV. A foreshock temperature of SNR shocks that is a significant fraction of the downstream temperature implies that the Mach number of SNR shocks is in fact reduced by at least an order of magnitude compared to the least upper bound , which estimates the ion acoustic speed using the ISM temperature eV.

The shock does not convert instantly the directed flow energy of the upstream into downstream heat. The shock potential and the ion distribution gradually change over a spatial interval that can be large if the Mach number is low, as our simulation shows. This has potentially important consequences for experimental observations of slow shocks. Their potential jump between the downstream region and the foreshock is moderate and it falls off over hundreds of Debye lengths. The electric field amplitude is thus low and smeared out over at least tens of Debye lengths. Slow shocks may thus not always be detectable in experiments that measure the electric field amplitude.

Finally our simulations have shown that the ions downstream have a bi-Maxwellian velocity distribution and they are not in a thermal equilibrium. The shock-reflected ion beam is thus not the only process by which a supercritical shock can get rid of an excess flow energy of the upstream plasma. This aspect together with the non-Maxwellian electron velocity distributions may help explaining why electrostatic shocks are stable even if electrons and ions have the same temperature, which is not possible according to some analytic models.

Acknowledgements: MED wants to thank Vetenskapsrådet for financial support. MP acknowledges support through grant PO 1508/1-1 of the Deutsche Forschungsgemeinschaft (DFG). LR acknowledges support from the ULIMAC grant, Triangle de la Physique RTRA. GS thanks the Leverhulme foundation (grant ECF-2011-383) for financial support. The computer time and support has been provided by the High Performance Computer Centre North (HPC2N) in Umeå.

## References

- (1) J. K. Truelove and C. F. McKee, Astrophys. J. 120 299 (1999).
- (2) E. A. Helder, J. Vink, and C. G. Bassa, Astrophys. J. 737 85 (2011).
- (3) K. M. Ferriere, Rev. Mod. Phys. 73 1031 (2001).
- (4) T. Honzawa, Plasma Phys. 15 467 (1973).
- (5) L. Romagnani, S. V. Bulanov, M. Borghesi, P. Audebert, J. C. Gauthier, K. Lẅenbrück, A. J. Mackinnon, P. Patel, G. Pretzler, T. Toncian, and O. Willi, Phys. Rev. Lett. 101 025004 (2008).
- (6) D. W. Forslund and C. R. Shonk, Phys. Rev. Lett. 25 1699 (1970).
- (7) D. W. Forslund and J. P. Freidberg, Phys. Rev. Lett. 27 1189 (1971).
- (8) G. Bardotti, and S. E. Segre, Plasma Phys. 12 247 (1970).
- (9) I. Langmuir, Phys. Rev. 33 954 (1929).
- (10) G. Sorasio, M. Marti, R. Fonseca, and L. O. Silva, Phys. Rev. Lett. 96 045005 (2006).
- (11) M. A. Raadu and J. J. Rasmussen, Astrophys. Space Sci. 144 43 (1988).
- (12) G. Sarri, G. C. Murphy, M. E. Dieckmann, A. Bret, K. Quinn, I. Kourakis, M. Borghesi, L. O. C. Drury, and A. Ynnerman, New J. Phys. 13 073023 (2011).
- (13) M. E. Dieckmann and A. Bret, Astrophys. J. 694 154 (2009).
- (14) S. D. Bale, M. A. Balikhin, T. S. Horbury, V. V. Krasnoselskikh, H. Kucharek, E. Mobius, S. N. Walker, A. Balogh, D. Burgess, B. Lembege, E. A. Lucek, M. Scholer, S. J. Schwartz, M. F. Thomsen, Space Sci. Rev. 118 161 (2005).
- (15) N. Sckopke, G. Paschmann, S. J. Bame, J. T. Gosling, C. T. Russell, J. Geophys. Res. 88 6121 (1983).
- (16) B. Lembege, and P. Savoini, Phys. Fluids 4 3533 (1992).
- (17) N. Shimada, and M. Hoshino, Astrophys. J. 543 L67 (2000).
- (18) B. Lembege, J. Giacalone, M. Scholer, T. Hada, M. Hoshino, V. Krasnoselskikh, H. Kucharek, P. Savoini, T. Terasawa, Space Sci. Rev. 110 161 (2004).
- (19) R. E. Lee, S. C. Chapman, and R. O. Dendy, Astrophys. J. 604 187 (2004).
- (20) M. Scholer, and S. Matsukiyo, Ann. Geophys. 22 2345 (2004).
- (21) S. C. Chapman, R. E. Lee, R. O. Dendy, Space Sci. Rev. 121 5 (2005).
- (22) T. Umeda, M. Yamao, and R. Yamazaki, Astrophys. J. 695 574 (2009).
- (23) T. Amano, and M. Hoshino, Astrophys. J. 690 244 (2009).
- (24) S. Matsukiyo, Phys. Plasmas 17 042901 (2010).
- (25) J. Niemiec, M. Pohl, A. Bret, and V. Wieland, Astrophys. J. 759 73 (2012).
- (26) A. Bret, Astrophys. J. 699 990 (2009).
- (27) J. P. Eastwood, E. A. Lucek, C. Mazelle, K. Meziane, Y. Narita, J. Pickett, R. A. Treumann, Space Sci. Rev. 118 41 (2005).
- (28) H. J. Völk, L. O. C. Drury, and J. F. McKenzie, Astron. Astrophys. 130 19 (1984).
- (29) L. O. Drury, F. A. Aharonian, D. Malyshev, and S. Gabici, Astron. Astrophys. 496 1 (2009).
- (30) J. Vink, R. Yamazaki, E. A. Helder, and K. M. Schure, Astrophys. J. 722 1727 (2010).
- (31) P. Ghavamian, J. M. Laming, and C. E. Rakowski, Astrophys. J. 654 L69 (2007).
- (32) C. Sack, and H. Schamel, Phys. Rep. 156 311 (1987).
- (33) D. W. Forslund, and C. R. Shonk, Phys. Rev. Lett. 25 281 (1970).
- (34) H. Karimabadi, N. Omidi, and K. B. Quest, Geophys. Res. Lett. 18 1813 (1991).
- (35) T. N. Kato, and H. Takabe, Phys. Plasmas 17 032114 (2010).
- (36) O. Buneman, Phys. Rev. Lett. 115 503 (1959).
- (37) J. M. Dawson, Rev. Mod. Phys. 55 403 (1983).
- (38) C. S. Brady, A. Lawrence-Douglas, and T. D. Arber, Phys. Plasmas 19 063112 (2012).
- (39) B. Eliasson, and P. K. Shukla, Phys. Rep. 422 225 (2006).
- (40) S. Sultana, G. Sarri, and I. Kourakis, Phys. Plasmas 19 012310 (2012).
- (41) Y. Kazimura, F. Califano, J. Sakai, T. Neubert, F. Pegoraro, and S. V. Bulanov, J. Phys. Soc. J. 67 1079 (1998).
- (42) A. Spitkovsky, Astrophys. J. 682 L5 (2008). 14 012901 (2007).