Articles | Volume 38, issue 6
Regular paper
24 Nov 2020
Regular paper |  | 24 Nov 2020

A deep insight into the ion foreshock with the help of test particle two-dimensional simulations

Philippe Savoini and Bertrand Lembège

Two-dimensional (2D) test particle simulations based on shock profiles issued from 2D full particle-in-cell (PIC) simulations are used in order to analyze the formation processes of ions back streaming within the upstream region after interacting with a quasi-perpendicular curved shock front. Two different types of simulations have been performed based on (i) a fully consistent expansion (FCE) model, which includes all self-consistent shock profiles at different times, and (ii) a homothetic expansion (HE) model in which shock profiles are fixed at certain times and artificially expanded in space. The comparison of both configurations allows one to analyze the impact of the front nonstationarity on the back-streaming population. Moreover, the role of the space charge electric field El is analyzed by either including or canceling the El component in the simulations. A detailed comparison of these last two different configurations allows one to show that this El component plays a key role in the ion reflection process within the whole quasi-perpendicular propagation range. Simulations provide evidence that the different field-aligned beam (FAB) and gyro-phase bunched (GPB) populations observed in situ are essentially formed by a Et×B drift in the velocity space involving the convective electric field Et. Simultaneously, the study emphasizes (i) the essential action of the magnetic field component on the GPB population (i.e., mirror reflection) and (ii) the leading role of the convective field Et in the FAB energy gain. In addition, the electrostatic field component El is essential for reflecting ions at high θBn angles and, in particular, at the edge of the ion foreshock around 70. Moreover, the HE model shows that the rate BI% of back-streaming ions is strongly dependent on the shock front profile, which varies because of the shock front nonstationarity. In particular, reflected ions appear to escape periodically from the shock front as bursts with an occurrence time period associated to the self-reformation of the shock front.

1 Introduction

While upstream ions of the incoming solar wind interact with the curved terrestrial bow shock, a certain percentage is reinjected back into the solar wind and propagates along the interplanetary magnetic field (IMF); they form the so-called ion foreshock. This population has been extensively studied both with the help of experimental data (Tsurutani and Rodriguez1981; Paschmann et al.1981; Bonifazi and Moreno1981a, b; Fuselier1995; Eastwood et al.2005; Oka et al.2005; Kucharek2008; Hartinger et al.2013) and numerical simulations (Blanco-Cano et al.2009; Lembège et al.2004; Savoini et al.2013; Kempf et al.2015; Savoini and Lembège2015; Otsuka et al.2018).

Even if we restrict ourselves to the quasi-perpendicular region (i.e., for 45θBn90, where θBn is the angle between the local shock normal and the IMF), different types of back-streaming ions are identified, and (a) some are characterized by a gyrotropic velocity distribution and form the field-aligned ion beam population (hereafter FAB); conversely, (b) others exhibit a nongyrotropic velocity distribution and form the gyro-phase bunched ion population (hereafter GPB). None of these populations yet have a well-established origin, and different mechanisms have been proposed for years (Möbius et al.2001; Kucharek et al.2004), including (i) scenarios based on the specular reflection (Sonnerup1969; Paschmann et al.1980; Schwartz et al.1983; Schwartz and Burgess1984; Gosling et al.1982) with or without the conservation of the magnetic moment and (ii) scenarios which invoke the leakage of some magnetosheath ions producing a low-energy FAB population (Edmiston et al.1982; Tanaka et al.1983; Thomsen et al.1983). Nevertheless, the origin of FAB ions could be due to (iii) the diffusion of some reflected ions (called gyrating ions; these ions are reflected by the supercritical shock front but do not manage to escape into the upstream region and go into the downstream region after their initial gyration (Schwartz et al.1983)). The diffusion can either be generated by upstream magnetic fluctuations (Giacalone et al.1994) or more directly by the shock ramp itself (with a pitch angle scattering during the reflection process) (Kucharek et al.2004; Bale et al.2005). All scenarios have some drawbacks and are not able to clearly explain the origin of both populations. On the other hand, GPB are preferentially observed at some distance from the curved shock front (Thomsen et al.1985; Fuselier et al.1986a), and their synchronized nongyrotropic distribution is part of a low-frequency monochromatic waves trapping (Mazelle et al.2003; Hamza et al.2006), or of beam plasma instabilities (Hoshino and Terasawa1985). As a conclusion, it is quite difficult to discriminate between these different scenarios which can be present simultaneously or separately in time.

Our previous papers (Savoini et al.2013; Savoini and Lembège2015) were focused on the origin of these two populations. A large-scale, two-dimensional particle-in-cell (PIC) simulation of a curved shock has been used, where full curvature and time-of-flight effects for both electrons and ions are self-consistently included. Our simulations have shown that both FAB and GPB populations and their typical associated pitch angle distributions observed experimentally (Fuselier et al.1986b; Meziane2005) have been retrieved not far from the front (up to 2–3RE, where RE is the Earth's radius). Moreover, results have shown that these two populations can be generated directly by the macroscopic electric E and magnetic B fields present at the shock front itself. In other words, the differences observed between FAB and GPB populations are not the result of distinct reflection processes but are the consequence of the time history of ions interacting with the shock front. The FAB population loses their initial phase coherency by suffering several bounces along the front, which is in contrast with the GPB population which suffers mainly one bounce (i.e., mirror reflection process). This important result was not expected and greatly simplifies the question on each population origin (Savoini and Lembège2015)

Nevertheless, some further questions, which are difficult to investigate with full PIC simulations (because of the self-consistency), still need to be answered in order to analyze several aspects of the reflection process. For this reason, we use complementary test particle simulations herein to clarify the respective impact of the shock curvature and the time variation in the macroscopic fields at the shock front on the back-streaming ion reflection process. The main questions presently addressed are summarized as follows:

  1. Is the reflection process noncontinuous in time (burst-type reflection process) or not? In this case, how is it linked to the θBn angle variation (i.e., space dependence) and/or to the shock profile variation (i.e., time dependence)?

  2. What is the impact of the space charge electric field localized within the shock ramp on the reflection process?

  3. What kind of reflection mechanisms can be identified in present simulations?

The paper is organized as follows. Section 2 briefly summarizes the conditions of the previous 2D PIC simulations (Savoini and Lembège2015) and of present particle test simulations. In Sects. 3 and 4, results of test particles are presented, and the ion reflection processes are investigated. The discussion and conclusions will be presented in Sects. 5 and 6, respectively.

2 Numerical simulation conditions

The numerical conditions concerned in the present paper are similar to those described in Savoini et al. (2013) and Savoini and Lembège (2015). In short, we used a 2D, fully electromagnetic, relativistic particle code based on a standard finite-sized particle technique (similar to Lembège and Savoini1992, 2002 for planar shocks).

2.1 Self-consistent full PIC simulations

The code solves Maxwell and Poisson's equations in the Fourier space (so-called pseudospectral code), which allows one to separate the electric field contribution in two distinct parts, namely (i) a longitudinal or electrostatic component, hereafter denoted by subscript l (built up by the space charge effects El=ρ/ϵo), and (ii) a transverse or induced component, hereafter denoted by subscript t (coming from the temporal variations in the magnetic field ×Et=-B/t). The longitudinal component is essentially built up within the shock front due to the different dynamics of ions and electrons, whereas the induced component is mainly generated by the propagating shock front itself (see Fig. 1, panel 2a) through the convective term Et=-U×B (where U corresponds to the bulk shock front velocity since we are in the solar wind frame). In addition, the subscripts and stand for parallel and perpendicular directions to the local magnetic field, respectively. In Fig. 1 and in following, the XY reference frame is the solar wind frame, with the third direction along Z pointing backward into the plot. Then, Et has the direction of the increasing Y, and B has the same direction as the present U vector.

Figure 1Panel (1) plots the simulation plane geometry; the B̃o magnetostatic field is mainly outside and directed downward from the plane. Panels (2a, b) illustrate the evolution of the magnetic field B̃tz in the fully consistent expansion (FCE) model (time dependent), respectively, at t̃init=2.4τ̃ci and t̃simul=5.4τ̃ci. In panel (2a), used as a reference, the vector velocity Ũ=ṽshock has been superimposed to illustrate the shock front propagation (white thick arrows); the arrow length is not at the right scale. In addition, the projection of the B̃o magnetic field lines has been reported (oblique white thin lines). Panels (3a, b) illustrate one example of the curved magnetic field B̃tz in the homothetic expansion (HE) model (time independent), where the shock profile is fixed in time but expands in space via an expanding factor proportional to the shock front velocity vshock×t; this shock profile has been chosen at the time t̃init=2.4τ̃ci of the self-consistent simulation. From this time, the shock front dilates by a factor of 2.6 compared to its initial shape.


In our configuration, the magnetostatic field is partially lying outside the simulation plane (see Savoini and Lembège2001 for more details). Then, the simulation is limited to the whole quasi-perpendicular shock (i.e., for 45θBn90). We use a magnetic piston for which the geometry is adapted to initiate a shock front with a curvature radius large enough to be compared with the upstream ion Larmor gyroradius ρ̃ci (all normalized quantities are indicated with a tilde ̃); this is the same normalization which is used in the previous self-consistent PIC simulations (Savoini and Lembège2001; Savoini et al.2013; Savoini and Lembège2015). The curvature increases during the simulation. This configuration has two consequences; (i) first, as the time increases and the shock front expands, its velocity slightly decreases and so does the Alfvén Mach Number MA from ≈5 to ≈3, where the velocity is measured at θBn=90 used as a reference angle; (ii) the time-of-flight effects are self-consistently included. Indeed, this ballistic process is observed when the upstream magnetic field lines are convected by the incoming solar wind. In present simulations (based on upstream rest frame), the curved shock front expands and scans different θBn values. As a result, back-streaming particles, collected at a given upstream location, come from different parts of the curved shock front, depending on their respective velocity.

Initial plasma conditions are summarized as follows: light velocity c̃=3 and temperature ratio between ion and electron population Tel/Tio=1.58. A mass ratio mi/me=84 is used in order to save CPU time, and the Alfvén velocity is ṽA=0.16. The number of cells in the simulation plane along each axis is NCX = NCY =8192150ρ̃ci, with the size of a grid cell Δx=Δy1ρ̃ce. The shock is in the supercritical regime. with a time-averaged Alfvén Mach number MA≈4 measured at θBn=90. In order to observe the early stage of the ion foreshock formation, the end time of the simulation is t̃simul=5.4τ̃ci (where τ̃ci is the upstream ion gyro-period), which is large enough to investigate the interaction of incoming ions with the shock front and the further formation of back-streaming ions.

2.2 Test particle simulations

In the present paper, we use all field components issued from the same previous PIC simulation as in Savoini and Lembège (2015); all components have been saved every ΔT̃=τ̃ci/20. Test particle simulations appear to be a straightforward way to evaluate the action of different field components on the ion dynamics. Indeed, the feedback effects of particles on electromagnetic fields are excluded in test particle simulations, and one can modify or cancel some field components independently from each other. This allows one to identify their specific actions on particles and on the resulting ion reflection processes.

Figure 1 plots an example of the two configurations used hereafter in this paper. Panels 2a–b show the fully consistent expansion model (hereafter named FCE), which corresponds to the results where test ions interact with the E and B fields issued from the self-consistent simulation and where both spatial inhomogeneities and nonstationarities are fully included. If this configuration is easy to understand, the so-called particular approach named the homothetic expansion model (hereafter named HE) shown in panels 3a–b is complementary. In this case, particles interact with a propagating fixed front profile, i.e., all-time profile variations are excluded; only spatial inhomogeneities of the shock front profile chosen at the selected time are included, as detailed in Sect. 4.

In the two configurations of FCE and HE, we inject test particles distributed within 10 individual sampling boxes located along the curved shock front (Fig. 2). This procedure allows us to analyze the impact of the front curvature (local θBn) on the formation of back-streaming ions. We follow a total of 1 million test particles. Then, each box has the same number of particles N=100 000, and they are initialized as a Maxwellian distribution with a thermal velocity vthi, which is the same as in the self-consistent simulation (Savoini and Lembège2015).

Figure 2Initial location of the 10 sampling boxes (labeled from zero to nine) which map the upstream ion foreshock region. All ions belonging to a given box are represented by the same color for statistical analysis only (Sects. 3 and 4). The θBn propagation angle, where each box is initially centered at time t=0, is reported above the corresponding identification number of the box, but these colors will not be used anymore in this paper.


Let us point out that the use of finite-sized sampling boxes at different initial θBn angles only estimates the location where the particles hit the shock but does not provide an exact value for the local θBn seen by the particles when these interact with the expanding shock front. Nevertheless, it is quite helpful when classifying the different types of particle interactions with the curved front. The sizes of all the identical boxes are chosen so that (i) along the curved shock front each box has an angular extension of 4, which is small enough to scan the different orientations of θBn and large enough for statistical constraints, and (ii) along the local shock normal, where each box has a length large enough (Lsize=2000Δx) to ensure that most particles interact with the shock front during a noticeable time range (i.e., DT3τ̃ci).

Then, Sect. 3 will present the results obtained in the FCE model, which is the usual configuration representing the time evolution of test particles with the self-consistent shock front profile. Section 4 will then introduce the more unusual HE expansion model (i.e., homothetic simulation approach).

3 Numerical results: the fully consistent expansion model (FCE)

To clarify the presentation, we will split this section into two parts in which we analyze, respectively, (i) the dynamics of the back-streaming ions (alias BI) in the different boxes and their main features, in terms of spatial and time evolution, and (ii) their behavior when some field components at the shock front are included/artificially excluded.

Figure 3FCE configuration: spatial distribution of the back-streaming particle density within the simulation XY plane at the end of the simulation time (t̃end5.4τ̃ci, where τ̃ci is the upstream ion cyclotron period). All boxes are plotted from NBox=0 to 9, defined initially at θBn90 and 45, respectively; the bottom panel labeled NBox=10 shows aggregate boxes in which we have reported the edge of the ion foreshock (dashed line) and the angle θBn=45 (dotted line) for reference. Considering the small number of ions involved in the reflection process, we have used a Gaussian interpolation which gives the relative density weight of each ion. Then, the color coding (vertical bar) gives only an indication of the relative density amplitude. The location of the curved shock front is defined at the middle of the front ramp (thick black line) at the last time t̃end. Moreover, we have reported the space-integrated percentage value BI% of back-streaming ions within each corresponding box. In order to exclude the gyrating ions present near the front from the back-streaming population, we have eliminated ions found within a small area ≈2–3ρci upstream of the shock front. For this reason, a very thin white area is visible along the curved front where no particles are present. For NBox=0–3, the arrows point to the two spots (see the text for explanations). The edge of the foreshock is represented by a mixed line at θBn=70. The dashed line (defined for θBn=45) is indicated for reference.


3.1 General features of the back-streaming ions

Figure 3 plots the spatial distribution of the back-streaming ions' density for the different sampling boxes (defined in Fig. 2) at the end time of the simulation for the FCE model. Different information can be summarized as follows:

  • i.

    The percentage of the back-streaming ions BI% is obtained by computing the ratio of the back-streaming ions over the number of ions which have interacted with the shock front. This number increases when moving further into the foreshock (i.e., for decreasing θBn) from BI%=0.1 for NBox=0 to BI%≈14 for NBox=9. This θBn dependence is in agreement with previous experimental observations (Ipavich et al.1981; Eastwood et al.2005; Mazelle et al.2005; Turner et al.2014) and numerical simulations (Savoini and Lembège2015; Kempf et al.2015).

  • ii.

    The upstream edge of the ion foreshock (mixed line in Fig. 3 for NBox=10) is not parallel to the IMF but is the result of the time-of-flight effect included in our simulations (Savoini et al.2013). At the end of the simulation, this edge starts from the shock at the same critical angle, the so-called θio,fore70, as found in our previous self-consistent simulations (Savoini et al.2013; Savoini and Lembège2015).

  • iii.

    The back-streaming ion density is not uniform along the shock normal but exhibits different maxima, according to our observations. Not only is the spatial distribution not the same for all boxes but is even not uniform within a given same box, i.e., back-streaming ions do not escape uniformly away from or along the shock front. For instance, boxes NBox=0–3 evidence two distinct spots near the shock front, as indicated by black arrows. As θBn decreases (i.e., NBox=5–9), the right-hand spot disappears and the back-streaming population increasingly aligns along the upstream magnetic field Bo. Accordingly, the width of the reflection area (i.e., the angular extension of the ion foreshock defined very near the shock front) shrinks from ≈50ρci (for NBox=0) to ≈17ρci (for NBox=9).

These two distinct spots may be explained by the different time histories of the back-streaming ions within the shock front, as reported in Savoini and Lembège (2015). Actually, the interaction time strongly differs from one ion to another, depending on its gyrating feature when it hits the shock front for the first time. Short and long interactions over time can be defined, depending on whether the reflection process is respectively associated to a short or long displacement of the ion along the shock front before escaping upstream. If the individual trajectory of the reflected ions has been already evidenced in Savoini and Lembège (2015), present test particle simulations allow one to generalize the results via a statistical approach versus their initial angular locations (i.e., the NBox number). Figure 4 plots the distribution of θBn angles seen by the particles when these hit for the first time the shock front (hereafter named θBnhit in red) and when they finally exit the shock front to escape upstream (hereafter named θBnexit in blue). These statistical results are obtained by computing these angles for each particle. As a consequence, angle values are computed neither at the same time nor at the same location along the curved front (even if they are initially located in the same box). In other words, each particle sees different local shock front profiles in terms of the spatial inhomogeneity and time nonstationarity of the shock front.

Figure 4FCE configuration: plots of the ion distribution functions for each box NBox=0–9 versus the local θBn angle computed when ions hit the shock front for the first time (red distribution function of the so-called θBnhit) and when these leave and escape upstream (blue distribution function versus the so-called θBnexit). The angles θBnhit and θBnexit have been reported in NBox=2 for reference.


First, let us note that the averaged value of θBnhit corresponds mainly to the initial location of the box, and therefore, the most important feature is not the angle values themselves but rather the difference between the averaged values of θBnhit and θBnexit when ions hit and leave the shock front, respectively. For this reason, we will use the angular range of the particles' interaction with the front as defined by ΔintθBn=θBnexit-θBnhit. Obviously, θBnhit decreases as NBox increases until it approaches the limit of the quasi-perpendicular domain of propagation, i.e., θBn=45 for NBox=9.

Second, the distribution functions of θBnexit strongly differ according to the respective box. For NBox=0–2, two (blue) peaks occur, namely one for high θBnexit (for which ΔintθBn≈45) and the other for low θBnexit (ΔintθBn15). In terms of the time trajectory, the presence of these two peaks suggests that some ions have spent different interaction times (subscript int) within the shock front. Some escape after a short interaction time (i.e., small ΔintθBn), while others escape after a long interaction time (i.e., large ΔintθBn), where the terms short and long refer to a small and large drift along the shock front, as already analyzed in Savoini and Lembège (2015). In other words, a small drift refers to one bounce whereas a large drift refers to multiple bounces along the shock front.

Moreover, as NBox increases (i.e., NBox≥3), the lower θBnexit distribution (i.e., correspondingly the largest ΔθBn) decreases rapidly in amplitude and disappears from NBox=6 (i.e., θBnhit56). Simultaneously, the other peak (i.e., correspondingly the smaller ΔintθBn) becomes dominant for all higher order boxes, meaning that fewer and fewer ions are associated to large drifts along the shock front.

Figure 5FCE configuration: plots of the ion distribution function (for each box NBox=0-9) versus the interaction time range ΔT̃int spent by each particle within the shock front. As shown, this interaction time range is not continuous, but it evidences distinct bursts of reflected ions (hereafter named f1 and f2), respectively, defined at ΔT̃int2500.25τ̃ci and ΔT̃int9501τ̃ci, where τ̃ci is the upstream cyclotronic period. A third burst f3 of reflected ions can be identified around ΔT̃int15001.5τ̃ci for NBox=0-2 but can be neglected as compared with the others.


Third, in order to complete information deduced from Fig. 4, Fig. 5 plots the number of reflected ions versus the time spent within the shock front. This interaction time ΔT̃int is defined as the time difference between the time associated to θBnexit and to θBnhit. Different main maxima of back-streaming ions' density are evidenced, namely f1, f2; a third maximum f3 can be also observed for boxes NBox=0–4, but its amplitude is too weak to be relevant in this discussion. One important feature is that f1 and f2 appear in all boxes and are independent of the box number. More precisely, f1 appears about ΔT̃int0.25τ̃ciτ̃cishock while f2 is observed at ΔT̃int1τ̃ci4τ̃cishock, where τ̃cishock is the local gyroperiod estimated within the shock front (at the middle of the ramp). This indicates that the reflection process is not uniform in time but leads to the formation of ion bursts associated to the shock dynamics, even if the number of ions which spend several gyroperiods τ̃cishock (i.e., ≈4 bounces) within the shock front is rapidly negligible. In addition, for NBox=0–2, f1 and f2 have a similar amplitude, which is not the case for NBox=4–9. In fact, a close look at f1 and f2 shows that f2 does not decrease in magnitude but rather the amplitude of f1 drastically increases from 10 (NBox=0) to 2500 (NBox=9). Then, the f2 population is always present but becomes negligible for lower θBn as compared with f1; this explains why we do not observe two distinct spots for NBox=4–9.

One helpful aspect of the test particle approach is including or excluding some electromagnetic field components in order to analyze their impact on the particles dynamics. Indeed, it is clear that some electric field component (i.e., El×B drift) and strong magnetic gradients drift (i.e., -B×B drift) can be a prerequisite for a large drift along the front (i.e., ΔintθBn15), whereas it could be unnecessary for the other case ΔintθBn≈45. Unfortunately, the shock front magnetic gradient cannot be canceled without the shock itself; then, we will focus our study by including or excluding the electric components which will shed new light on the origin of back-streaming ions filling the foreshock.

3.2 Impact of electric field components

Savoini and Lembège (2015) have analyzed the impact of E×B drift velocity on the dynamics of back-streaming ions (Gurgiolo et al.1983) and, more particularly, as a source of FAB and/or GPB populations. This study has shown that the origin of both populations can be easily explained in terms of E×B drift associated or not to a diffusion in the velocity space, but it was not able to explain the details of the reflection mechanism itself. Then, herein, we will focus on the role of the electrostatic field component Ẽl built up within the shock front (i.e., space charge effects). This longitudinal component, defined along the normal to the shock front, can be associated to the electrostatic potential wall responsible for some reflected ions. In the case of a constant shock profile in time with a planar geometry, this reflection does conserve the energy since the potential is the same before and after the reflection, and the total work of the electric force is canceled. Nevertheless, in more realistic conditions, this scenario is not valid anymore for ions which drift along the shock front and suffer both time and space electric and magnetic field variations. Then, in the following sections, we will preferentially use the field El rather than the potential Φ.

Figure 6FCE configuration: characteristics of the ion foreshock with and without the electrostatic field component Ẽl (i.e., Ẽlx=Ẽly=0). Panel (a) shows the percentage of the back-streaming ions BI% versus the box number. Black and red lines are defined for Ẽl0 and Ẽl=0, respectively. Panel (b) shows the density of the back-streaming particles in the same format as Fig. 3 but when Ẽl=0. Only the view which aggregates all boxes (i.e., NBox=10) is shown in order to evidence the location of the edge of the ion foreshock (dotted line) in each case.


Figure 6a shows the percentage BI% of back-streaming ions versus the box number where the black and red curves are defined for Ẽl0 and Ẽl=0 respectively. The impact of the Ẽl field (i.e., the potential wall) on the reflection process is clearly apparent for all θBn values (namely for each NBox number). In particular, the percentage BI% strongly decreases as Ẽl=0, which illustrates the dominant role of Ẽl field whatever the box is. This is especially true for the lower box number NBox=0–3 (i.e., high θBn approaches 90) where very few back-streaming ions are observed. This point is not surprising if one remembers that this electrostatic field decelerates the incoming ions (i.e., accelerates ions in the present reference frame) and contributes to the reflection process. In other words, for NBox=0 (i.e., largest θBn), escaping ions have to be accelerated to higher parallel velocity, as reviewed in Burgess et al. (2012). Let us stress that Fig. 6a exhibits a clear change in the slope of BI% increase at the box NBox=2 centered around θBn=70. Herein, we will consider this value as the reference angle identifying the starting location of the ion foreshock edge attached to the shock front. This value is in reasonable agreement with the value (θio,fore66) found approximately in the previous self-consistent PIC simulations Savoini and Lembège (2015).

Another consequence is illustrated in Fig. 6b, which shows that the edge of the ion foreshock is shifted due to the lack of reflected ions, and it starts around θBn55. Clearly, the contribution of the electric field is important for ions which populate the edge of the foreshock and need to escape at high θBn.

Another way to observe the strong impact of Ẽl on the dynamics of reflected ions is illustrated in Fig. 7, which shows the two θBnhit and θBnexit distributions in the same format as that of Fig. 4. The number of back-streaming ions decreases drastically for all boxes and is zero for NBox=0. Furthermore, the density is not uniform for all boxes and appears to be much more important for NBox=5–9 than for NBox=1–4. The θBnexit distribution is strongly modified and a comparison between Figs. 4 and 7 can be summarized as follows:

  1. The boxes NBox=1–2 evidence a total absence of reflected ions with a small range ΔintθBn, and only the θBnexit distribution around 60 persists. This result shows that Ẽl field plays a key role, more specifically on back-streaming ions suffering a one bounce reflection near the edge of the ion foreshock (i.e., ΔintθBn≈45).

  2. The θBnexit distribution in both cases (Ẽl0 in Fig. 4 and Ẽl=0 in Fig. 7) is roughly similar in corresponding boxes for all NBox≥6. This can be interpreted as either that the ions have been accelerated enough during their reflection at the shock front or that they need a lower parallel velocity to escape upstream. As a consequence, the El component is not mandatory anymore, and the mirror magnetic mechanism at the shock front can be invoked as the only reflection process.

In summary, the comparison between Figs. 4 and 7 evidences that Ẽl components are essential in the ion reflection for a high θBn angle (>56, i.e., NBox=6–9) where ions need a strong acceleration process but play a less important role at lower θBn angles. Conversely, for NBox=6–9, the reflection process takes place with a very small ΔintθBn, with or without the electric field. In other words, the large shock drift invoked for θBn56 seems to be mainly supported by the convective electric field Ẽt components present at the shock front.

Similarly, the one bounce reflection always occurs, even in the absence of Ẽl field (i.e., in the absence of the shock front potential wall), and can then be associated to a magnetic reflection which seems to be very efficient, especially at lower θBn. This one bounce reflection (i.e., f1) corresponds essentially to a short interaction time, as illustrated in Fig. 5 (ΔT̃int1τ̃cishock). Then, the ion energy gain is essentially due to a Fermi-type acceleration.

Figure 7FCE configuration: same plots as for Fig. 4 but when electrostatic field components Elx and Ely are artificially excluded (El=0). The angles θBnhit and θBnexit have been reported in NBox=2 for reference purposes.


3.3 Impact of the shock front nonstationarity

Previous studies have largely evidenced that a quasi-perpendicular shock front can be intrinsically nonstationary due to different mechanisms (for a review, see Lembège et al.2004; Marcowith et al.2016). Then, it is important to analyze the impact of such nonstationarity on the temporal ion foreshock dynamics. As a first step, we plot in Figure 8 the time evolution of the back-streaming ions' percentage BI% as these leave the front and escape into the upstream region, where BI% is the instantaneous rate computed during a short time range ΔT̃=τ̃ci/20. The time T̃init=1248 is the initial time when test particles are launched into the time-dependent simulation.

Figure 8FCE configuration: time history of the reflected ion percentage BI%; each value is computed during a short integrated time interval ΔT̃=τ̃ci/20 for the different boxes. Within this interval, only the newly reflected ions are memorized. As for Fig. 6, black and red lines correspond to the case where the electric field is included (Ẽl0) or excluded (Ẽl=0), respectively. We have used the same normalization in order to compare both cases. The blue dotted straight lines are not obtained from a data linear approximation but only indicate the mean values of the time modulations.


  • i.

    Results of Fig. 8 are obtained as the Ẽl field components are included (black curve) and artificially excluded (red curve). One concludes that the percentage BI% strongly decreases as Ẽl components are excluded, and that the impact of Ẽl field is emphasized for lower NBox. In other words, the back-streaming ions mainly appear for higher NBox>5 (i.e., for lower θBn), even in absence of Ẽl field.

  • ii.

    The different results may be classified into two groups, where (i) the first one concerns boxes NBox=0–4 showing a slow increase (almost monotonic) in the reflection rate and (ii) a second group NBox=5–9 evidences a steep increase followed by a flat-top shape around BI%≈1, even if it increases slightly with NBox. At the end of the simulation, the strong decrease in BI% observed for all boxes corresponds to the time when all ions of the different boxes have been swapped by the propagating shock front and then no more ions are back streaming.

    For the first group NBox=0–4, a delay is observed in the formation of back-streaming ions between the different boxes, although test particles are initially evenly distributed in the whole boxes. For NBox=0, back-streaming ions appear around T̃5248, i.e., well after the initial release time T̃init=1248 (i.e., ΔT̃=5248-1248, which covers 2.6τ̃ci), whereas this time delay decreases to T̃2018 (i.e., ΔT̃=2018-1248=770, which covers 0.5τ̃ci) as NBox increases. This illustrates the larger time delay of ions having interacted with the front to escape upstream at high θBn. For NBox=0–2, ions have to stay longer within the shock front to finally escape at lower θBnexit (Fig. 7), which is illustrated by the increase of BI% as the time evolves (Fig. 8). The second group (NBox=5–9) concerns boxes which are already at a lower angle θBn with easier escaping conditions. In this case, ions are reflected continuously with some time variation in BI%.

  • iii.

    Another interesting feature is the presence of different modulations which are superimposed on a time-averaged reflection ion rate (blue dashed line), especially for the boxes NBox=1–5 and the boxes NBox=7–8. Nearly all boxes exhibit these modulations which represent about 40 % of the averaged BI%, although the amplitude of these modulations varies versus time and NBox. These evidence a nonstationary ion escaping rate. These modulations almost disappear in the case of Ẽl=0 (red line). This illustrates that the presence of the Ẽl field component is a key ingredient in the formation of these modulations since this electric field component is also involved in the shock front self-reformation, as described in previous works (Lembège and Savoini1992; Scholer et al.2003; Matsukiyo and Scholer2006).

    This result confirms the importance of the electrostatic field component at the shock front in the reflection processes of the back-streaming ions and, most importantly, on the ion foreshock nonstationarity behavior as described in a previous paper (Savoini and Lembège2015). Nevertheless, it is quite difficult to establish a one-to-one correspondence between these modulations and the nonstationarity of the shock front because, during the sampling time interval ΔT̃=τ̃ci/20, the front nonstationarity and the time-of-flight effects have mixed ions coming from either different times and/or different θBnexit regions (even if they are in the same box). So this self-consistent approach is not totally adapted to resolve this question. A complementary approach is necessary, based on simplified simulations with fixed shock front profiles in expansion (nonstationary effects are excluded). This motivates the homothetic expansion model (HE) described in the next section.

4 Numerical results: the homothetic expansion (HE) model

4.1 Descriptions of the HE model

Let remember that all simulations are made in the solar wind reference frame (i.e., the curved shock expands into the upstream region where the solar wind is at rest). As a consequence, if one follows test particles within this configuration, we have to mimic this behavior. In order to proceed, we apply a homothetic transformation (homogeneous dilatation in all directions) with an expansion factor deduced from the shock front velocity determined at selected times, as illustrated in panels 3a–b of Fig. 1. Special attention has been taken in the determination of this homothetic factor λ=vshock×t, since the shock front velocity vshock at a given time must fit with the corresponding value issued from the PIC simulations (Savoini and Lembège2015). With this information, we are able to expand the shock front through a cubic interpolation as it propagates with vshock in an expanding simulation plane (i.e., the grid cells stay constant Δ̃x=Δ̃y=cte but the number of the grid cells increases accordingly). In other words, all points of electromagnetic fields at the shock front follow the relation OM(vshock×t)OM, where vshock is the value of the shock velocity, as computed from our self-consistent 2D PIC simulations at the selected time, and OM is the vector between the initial location of the shock front (i.e., the point O) and any point of the field array (i.e., the point M). At this stage, we have to point out that the velocity vshock remains artificially constant during the whole simulation, which is not the case for the FCE model where vshock slightly decreases. Then, the same procedure is repeated for another selected time, so that one can analyze the impact of different shock front inhomogeneities and curvature on ion dynamics; let us note that time-of-flight effects are always included. Each front profile is analyzed within a similar simulation time range 3τ̃ci. In summary, similar simulations are performed for 174 different times in order to simulate all the different shock profiles provided by the 2D PIC self-consistent simulation from t̃init=1.2τ̃ci to t̃simul=5.4τ̃ci.

4.2 General features of the back-streaming ions

Figure 9 has been achieved in the same format as Fig. 8 by performing 100 independent simulations (i.e., we take only the first 100 simulations so that all test particles hit the propagating shock front). During this range, the shock front is forced to expand (see Sect. 4.1). For example, the chosen time T̃=44564.2τ̃ci on the abscissa axis corresponds to one particular shock front profile (i.e., including all electromagnetic field components issued from the previous self-consistent PIC simulations) from which we have measured the instantaneous shock front velocity and that we follow during the time range covering 3τ̃ci.

The purpose is to determine (i) whether some shock profiles are more appropriate than others for the formation of back-streaming ions and (ii) if, yes, whether better reflection takes place for some particular angular range of θBn.

Figure 9HE configuration: percentage BI% of back-streaming ions measured at the end of each simulation, where each time corresponds to a given fixed shock front. For each shock profile in homothetic expansion, the simulation covers 3τ̃ci, allowing one to obtain a well-developed ion foreshock, and BI% represents the ratio of the back-streaming ions over the total number of upstream ions which are released at the beginning of the simulation within a given box. As in Fig. 7, black and red lines correspond to the case when the El field is included and artificially excluded, respectively. The concerned shock profiles are chosen only at late times of the full PIC simulation (from T̃=3456 to 5474), where the curvature radius of the shock front is relatively large (R̃>70ρ̃ci).


The comparison of Figs. 8 and 9 provides the following information. The maximum BI% value is much higher for the HE model than for the previous FCE model for each corresponding box. For example, for NBox=9, BI%max1.2 in the FCE simulations as compared with BI%max15 in the HE configuration. In fact, one has to remember that for the FCE model the BI%max value represents an instantaneous reflection rate versus the shock front evolution, whereas this rate is a time-integrated value for the HE model. Indeed, in this model, the shock front profile stays the same during the whole simulation, and then, if this profile allows the reflection of some incident ions, they will be reflected continuously, leading to a high BI% number.

Obviously, the main information is not the BI% value itself but rather its evolution versus time and for different shock profiles. Other main results issued from Fig. 9 may be summarized as follows:

  1. The NBox=0 box evidences almost no reflection for the majority of the shock front profiles, which indicates that nonstationarity effects present in the FCE configuration (i.e., Fig. 3) are needed for feeding back-streaming ions along the edge of the foreshock.

  2. Boxes NBox=1–8 clearly show some strong modulations in the percentage BI% versus the shock profile of concern, which corresponds to a quasi-periodic bursty emission of back-streaming ions. The BI% rate periodically reaches a maximum value followed by a minimum around zero. The corresponding time period ΔT̃max is about 460=0.5τ̃ci (between two successive maxima). The temporal width of each maximum is about ΔT̃range256=0.25τ̃ci. These modulations mean that conditions for the formation of back-streaming particles are not continuous but correspond to some specific shock front profiles. In addition, these modulations appear synchronized in time for the different boxes 1–7, which implies that the local reflection conditions are not strongly dependent on the θBn angle but rather depend on the shock profile at certain times.

  3. In contrast, the boxes NBox=8–9 also evidence the same kind of modulations but with greatly reduced amplitude; these are even nonexistent between T̃3456 and 4456, which indicates a low sensitivity to the shock front profile when approaching θBn=45.

  4. Similarly, the maximum values BI%max (black curve) change drastically with box numbers, namely from small amplitudes, BI%≤6 for NBox=1–2, to very high values, BI%≈30 for NBox=3–5, before decreasing again for NBox=6–9. These variations may be understood by taking into account the reflection processes present at these different θBn and, more specifically, in regards to the electric field components. For boxes NBox=0–2, reflection is almost impossible without the Ẽl component (i.e., electric potential wall). As θBn decreases (for NBox=3–5), the reflection becomes easier and both magnetic and electric fields contribute to the percentage of reflected ions. Finally, for the last boxes NBox=6–9, the peak amplitude decreases but the contribution of Ẽl becomes less important in the reflection process, as evidenced by comparing both black (Ẽl0) and red (Ẽl=0) curves. Instead, another process, essentially driven by the magnetic field (mirror reflection) contributes more since the peak amplitude of the red curve increases progressively as NBox increases from 6 to 9.

  5. Only the last peak around T̃5456 has a different behavior in comparison with others. In particular, we observe that, for these different times or profiles, the presence of the Ẽl leads to higher BI%. This behavior can be understood because vshock is lower for these times, and the Ẽl component is necessary for decelerating ions and reflecting them. Without this component, only the magnetic reflection is present, and BI% has the same amplitude as the previous maximum.

  6. Finally, Fig. 9 confirms the key role of El field in back-streaming ion formation, except when approaching θBn=45 (NBox=8–9) while another reflection process is also at work (in absence of El). This represents an indirect way of stressing the noticeable impact of the magnetic field in this angular range. This magnetic reflection process is more evidenced at lower θBn since ions need less parallel velocity to be reflected back into the upstream region. This statement can be quantified more precisely, as explained in Sect. 5.

This result is an illustration of the impact of the electrostatic component at the shock front. As is well known, this component works to decelerate incoming ions (i.e., accelerate in our solar wind frame) and to accelerate electrons (i.e., decelerate in our solar wind frame) to the downstream region (Savoini and Lembège1994; Bale et al.2005). As a consequence, this electrostatic component is an essential ingredient in the formation of back-streaming ions, especially at higher θBn.

5 Discussions

A previous paper (Savoini and Lembège2015) demonstrated that all reflected ions suffered the same E×B drift in the velocity space, which can account for the pitch angle distributions observed at the shock front. In fact, the key point is the time spent by particles within the front shock, which finally decides whether ions will escape to form the FAB (with a pitch angle α0) or GPB populations (with a pitch angle α0), where α is the angle between the velocity vector and the magnetic field. In other words, the FAB population may be associated to a large drift along the shock front (and/or long interaction time) during which particles see a time-varying shock front and lose their phase coherency; this case corresponds to a large angular range ΔintθBn mentioned in Sect. 3.1. In contrast, the ions of the GPB population have a shorter interaction time with the shock front associated to a small angular range ΔintθBn. Present test particle simulations allow one to have a deeper insight to the spatial origin of the observed FAB and GPB populations. Then, we have to analyze the ion velocity distribution more carefully.

Figure 10 plots the local perpendicular velocity distribution functions f(v1,v2) in both FCE and HE approaches (where v⟂1 and v⟂2 refer the ion perpendicular velocity components defined with respect to the local magnetic field). All plots are obtained at the end of the simulations and take into account the different populations observed in Figs. 4 and 7 (i.e., both distinct peaks of θBnexit angle for lower NBox are included); results issued from FCE (Fig. 10, left panel) and HE (Fig. 10, right panel) configurations are considered. For the HE configuration, we chose an initial time T̃=4848 (see Fig. 9 for reference) which corresponds to a maximum of BI% in order to have enough reflected ions in the velocity space. Results from the three different boxes, NBox=1,4 and 8, are represented in order to give an overview of the whole ion foreshock components.

Figure 10Local perpendicular ion velocity space (v⟂1, v⟂2) of all back-streaming ions computed at the end time of the simulations (i.e., after 3τci for all simulations) for boxes NBox=1,4 and 8 in the FCE approach when the Ẽl is included (left panel); the case Ẽl=0 is not plotted since the percentage BI% is very weak (see Fig. 8). The right panel shows similar results for the same boxes in the HE configuration (corresponding to T̃=4848) in both cases where Ẽl is included and artificially excluded; statistical results where BI% is too weak are not shown. Red and blue hold for maximum and minimum density value in the velocity space.


Results of the FCE configuration (left panels) can be analyzed in Figs. 4, 7 and 10. Plots of the El≠0 case (Fig. 10) show that the NBox=1 has a low number of upstream reflected ions, which leads to poor statistics and a noisy f(v1,v2) distribution. Nevertheless, it evidences approximately a distribution with a maximum slightly noncentered at v=0. Then, this distribution can be viewed as a mixing of GPB and FAB populations, even if the GPB population with a pitch angle different zero is the largest one. When moving further into the foreshock region (i.e., lower θBn angle with NBox=4), the number of reflected ions drastically increases and we observe more clearly the characteristic partial ring of the GPB population (as in Savoini and Lembège2015). In addition, the center of the ring is also partially filled in because of partial diffusion due to particles having a large ΔintθBn range, corresponding to the θBnexit peak around 60 in Fig. 4, and/or by the intrinsic time fluctuations of the front which tends to blur out the velocity distribution both in perpendicular and parallel directions. At last, in agreement with the associated small θBnexit of Fig. 4, NBox=8 (in Fig. 10) also evidences a non-Maxwellian-like distribution (α0). When we look at the El=0 case, we observe roughly the same behavior for all boxes, even if the decrease in the back-streaming ions number in box NBox=8 makes the comparison difficult.

A further analysis requires a similar approach with the HE configuration in which we follow a succession of independent expanding shock profiles in order to exclude the impact of the time fluctuations on the velocity distribution f(v,v); then, no ion diffusion associated to these fluctuations is allowed. Results of the HE configuration (right panel, Fig. 10) show reflected ions for NBox=1 when Ẽl0, but once again, no reflected ions can be seen when Ẽl=0. This evidences the importance of the electrostatic potential wall in order to reflect upstream ions for high θBn. On the other hand, if the amplitude of the perpendicular velocity is roughly the same for the two different configurations (FCE and HE), the HE case shows a very well-formed ring, in contrast with the FCE case, which exhibits a diffuse velocity space. This illustrates that field time variations (FCE configuration) are much more efficient at diffusing particles than the fields spatial variations (HE configuration). Similarly, NBox=4 exhibits a clear ring, which is a feature of the GPB population. The center of the ring is not partially filled in since time velocity diffusion is excluded. These results demonstrate that the formation of the FAB-like population is also mainly due to ion velocity diffusion related to the time fluctuations of the shock front which can have different origins, as described in previous works (Kucharek et al.2004; Bale et al.2005). Finally, for NBox=8, the velocity distribution drastically changes from a ring (Ẽl0) to a localized bump (Ẽl=0) roughly similar to the FCE case. It is clear that the number of reflected ions decreases drastically as Ẽl field components are artificially suppressed. But, more important is that the formation of a nongyrotropic distribution does not depend strongly on the Ẽl component and is mainly controlled by the convective electric field through the Et×B drift in the velocity space, as already described in Savoini and Lembège (2015).

Let us remember that each distribution results from a combination of all particles originating from one given box (no matter where they end up spatially). This differs from the more common strategy based on measurements of local ion distributions, as performed in Savoini and Lembège (2015), but which did not establish, at that time, which part of the curved shock front the FAB and GPB ions are issued from. A further analysis is needed and is left for later work.

Figure 11Sketch of the ion foreshock in the quasi-perpendicular shock region, illustrating the angular areas along the curved front where the four main identified processes contributing to and/or impacting on the formation of back-streaming ions apply (namely the longitudinal electric field El, the magnetic field, the shock front nonstationary and the convective electric field Et). Each process is illustrated by different thick band along the curved front which is shifted, one from the other, in order to avoid overwhelming the sketch. One color (red, blue, green and black, respectively) is associated to each process. The varying intensity of the color indicates where the process is strong (full color) or weak (white). This allows one to identify, at a glance, the angular areas where the different processes are complementary or accumulating. The upstream interplanetary magnetic field (IMF) is reported in gray (dotted and dashed lines), as is the shock front itself. Typical directions θBn=90 and 70 (blue and red dashed lines) defined between the local shock normal n and the IMF indicate the location where the electron and ion foreshock edge initiates from the curved shock front, respectively. The electron foreshock edge is indicated as a reference.


6 Conclusions

Present test particle simulations reinforce the scenario described in Savoini and Lembège (2015) and have allowed us to investigate the formation of the ion foreshock more deeply. In summary, the impact of the three quantities or effects has been identified as follows: (i) the electric field (separately the Ẽl and the Ẽt components), (ii) the magnetic field and (iii) the shock front nonstationarity. The synoptic of Fig. 11 summarizes the importance of each impact versus the shock front curvature from the edge of the ion foreshock (θBn70) to θBn=45. In this sketch, the colors vary from strong (full color) to weak (white) intensities as a function of their respective influence on the ion dynamics. These different effects are the following:

  1. Impact of the El field on the ion reflection process. As is well known, the built-up potential wall at the shock front (i.e., the electric field El along the shock normal n) is mainly responsible for the deceleration (i.e., acceleration in our reference frame) of the incoming upstream ions by the shock front. The El component has essentially two distinct impacts where, (i) without this electric component, no reflected ions are observed for θBn>62, whereas in presence of this electric component reflected ions (having suffered even one bounce) can be observed at the edge of the ion foreshock (θBn70); (ii) at lower angles (θBn50) many ions are reflected without the help of the Ẽl component and can be associated to a magnetic mirror reflection. Then, in Fig. 11, Ẽl is only reported as being strong around the edge of the ion foreshock to emphasize its mandatory action for high θBn angles.

  2. Impact of the Et field on the ion reflection process. Figure 10 evidences that the convective electric component Et is always present in our simulation (we are in the solar wind reference frame and then Et≠0 within the curved propagating shock front). Our previous work (Savoini and Lembège2015) was only able to show that the E×B drift scenario in the velocity space foreseen by Gurgiolo et al. (1983) was at the origin of two distinct GPB (i.e., one bounce) and FAB (i.e., multibounces) populations only separated from the particle time history within the shock front.

    But this scenario was not able to distinguish the relative importance between the two electric field components of El and Et, respectively. This question has been clarified in the present paper since FAB and GPB populations formed by the E×B drift in the velocity space are evidenced with or without the El component (Fig. 10). Then, the El field seems to be less dominant than the Et field. Finally, Fig. 11 illustrates the E×B drift impact by a dark color that is almost uniform within the whole quasi-perpendicular region.

  3. Impact of the B field on the ion reflection process. The magnetic field is important for several reasons: (i) its increase at the shock front builds up the El component (space charge effect), which reflects back incoming low-energy ions, and (ii) more importantly, it is also directly responsible for the reflection of some ions (i.e., through the magnetic mirror reflection) and for the drift along the shock front of the multibounce ions (i.e., FAB population). Then, this population gains energy as ions propagate in the Et direction along the shock front. Nevertheless, as θBn decreases from 90 to 45, the ion reflection becomes easier since the parallel guiding center velocity needed to overcome the shock front velocity decreases (Paschmann et al.1980). This behavior is clearly illustrated herein by the increase in the percentage of reflected ions BI% as θBn decreases (i.e., NBox increases). This behavior persists even in absence of El, where BI% is only reduced by a factor of 2.5, as illustrated by Fig. 6. Then, if the magnetic field is important in the whole quasi-perpendicular region, we emphasize, in Fig. 11, its stronger impact near θBn45 where it is mainly responsible for the reflection (i.e., magnetic mirror) and acceleration (i.e., Fermi-type) of ions (Webb et al.1983).

  4. Impact of the shock front nonstationarity. Present simulations show that the reflection process is not continuous in both time and space but strongly depends on the local shock front profile met by incoming ions at their hitting time. This behavior is difficult to identify in experimental measurements since the particles coming from different shock locations and at different times are mixed; in contrast, this can be easily evidenced in our HE test particles configuration. This configuration evidences that particular shock profiles are more suitable for the formation of back-streaming ions than other ones. Indeed, we observe modulations of the BI% percentages in Fig. 9 which are much more pronounced than in our FCE configuration (Fig. 8). These modulations are so strong that BI% drops to zero periodically, which means that for certain shock front profiles no ion can escape into the upstream region. This behavior is observed for all NBox at the same time (i.e., same shock profile), which implies that the ion reflection does not depend on the location along the shock front but essentially on the global profile of the shock at a given time. In the present simulations, we can identify four distinct and noticeable bursts (i.e., maxima BI% values) with an average cyclic occurrence period of 1τ̃cishock (where τ̃cishock is defined at the shock ramp). Surprisingly, NBox=0 (Fig. 9) does not show the same bursts as the others, which suggests that the time variations in the shock front (and associated particle diffusion, as suggested by Kucharek et al.2004) are mandatory for obtaining back-streaming ions around the edge of the ion foreshock (i.e., high θBn). This point will require a further investigation.

In summary, present results show that the formation of the ion foreshock is not a continuous process but must be considered as time dependent, which leads to bursty emissions of back-streaming ions. Three different contributions have been evidenced, namely (i) the El component in the ion global reflection process, in particular for high θBn, (ii) the magnetic field B essentially observed when El=0 for lower θBn, such as the magnetic mirror reflection, and (iii) finally, the Et×B drift in the velocity space mainly sustained by the convective electric field which is necessary to generate both FAB and GPB populations as described in Savoini and Lembège (2015).

Unfortunately, the impact of the shock front nonstationarity on the ion foreshock is difficult to analyze (see, for example, Fig. 7) for two different reasons, namely (i) the time-of-flight effects mix reflected ions coming from different shock profiles, and (ii) even if some shock profiles are more efficient than others at reflecting ions, their respective impacts disappear rapidly since they are being blurred out by the impact of less efficient profiles on particles as time evolves.

Data availability

Test particle data have been deposited in the archive of the LPP laboratory at (Savoini2020). A “Readme” file has been added, which includes the creators, title, and date of last access as well as the description of the HDF5 data structure produced and a Python program that can be used to read the data directly and display the particle trajectories.

Author contributions

PS and BL contributed to the design and implementation of the research, the analysis of the results and the writing of the paper. The data were produced by PS.

Competing interests

The authors declare that they have no conflict of interest.


Numerical simulations were performed at the TGCC high-performance computer center located in Bruyères-le-Châtel (near Paris), and we thank the staff for their support (DARI project; grant no. A0050400295). ISSI (Bern, Switzerland) is thanked for supporting the collaboration network for “Resolving the Microphysics of Collisionless Shock Waves”. The authors wish to thank Yann Pfau-Kempf and the anonymous referee for their helpful comments.

Financial support

This research has been supported by Labex PLAS@PAR (grant no. ANR-11-IDEX-0004-02) and by the Centre National d'Études Spatiales (CNES; grant nos. APR-W-EEXP/10-01-01-05 and APR-Z-ETP-E-0010/01-01-05).

Review statement

This paper was edited by Minna Palmroth and reviewed by Yann Pfau-Kempf and one anonymous referee.


Bale, S. D., Balikhin, M. A., Horbury, T. S., Krasnoselskikh, V. V., Kucharek, H., Möbius, E., Walker, S. N., Balogh, A., Burgess, D., Lembège, B., Lucek, E. A., Scholer, M., Schwartz, S. J., and Thomsen, M. F.: Quasi-Perpendicular Shock Structure and Processes, Space Sci. Rev., 118, 161–203,, 2005. a, b, c

Blanco-Cano, X., Omidi, N., and Russell, C. T.: Global Hybrid Simulations: Foreshock Waves and Cavitons under Radial Interplanetary Magnetic Field Geometry, J. Geophys. Res., 114, 01216,, 2009. a

Bonifazi, C. and Moreno, G.: Reflected and Diffuse Ions Backstreaming from the Earth's Bow Shock: 2. Origin, J. Geophys. Res., 86, 4405–4413, 1981a. a

Bonifazi, C. and Moreno, G.: Reflected and Diffuse Ions Backstreaming from the Earth's Bow Shock. I Basic Properties, J. Geophys. Res., 86, 4397–4404, 1981b. a

Burgess, D., Möbius, E., and Scholer, M.: Ion Acceleration at the Earth's Bow Shock, Space Sci. Rev., 173, 5–47,, 2012. a

Eastwood, J., Lucek, E., Mazelle, C., and Meziane, K.: The Foreshock, Space Sci. Rev., 118, 41–94,, 2005. a, b

Edmiston, J. P., Kennel, C. F., and Eichler, D.: Escape of Heated Ions Upstream of Quasi-Parallel Shocks, Geophys. Res. Lett., 9, 531–534,, 1982. a

Fuselier, S. A.: Ion Distributions in the Earth's Foreshock Upstream from: The Bow Shock, Adv. Space Res., 15, 43–52, 1995. a

Fuselier, S. A., Thomsen, M., Gary, S., and Bame, S.: The Phase Relationship between Gyrophase-Bunched Ions and MHD-like Waves, Geophys. Res. Lett., 13, 60–63, 1986a. a

Fuselier, S. A., Thomsen, M. F., Gosling, J. T., Bame, S. J., and Russell, C. T.: Gyrating and Intermediate Ion Distributions Upstream from the Earth's Bow Shock, J. Geophys. Res., 91, 91–99,, 1986b. a

Giacalone, J., Jokipii, J. R., and Kota, J.: Ion Injection and Acceleration at Quasi-Perpendicular Shocks, J. Geophys. Res., 99, 19351–19358, 1994. a

Gosling, J. T., Thomsenn, M. F., Bame, S. J., Feldman, W. C., Pashmann, G., and Sckopke, N.: Evidence for Specularly Reflected Ions Upstream from the Quasi-Parallel Bow Shock, Geophys. Res. Lett., 9, 1333–1336, 1982. a

Gurgiolo, C., Parks, G. K., and Mauk, B. H.: Upstream Gyrophase Bunched Ions: A Mechanism for Creation at the Bow Shook and the Growth of Velocity Space Structure through Gyrophase Mixing, J. Geophys. Res., 88, 9093–9100,, 1983. a, b

Hamza, A. M., Meziane, K., and Mazelle, C.: Oblique Propagation and Nonlinear Wave Particle Processes, J. Geophys. Res., 111, 04104,, 2006. a

Hartinger, M. D., Turner, D. L., Plaschke, F., Angelopoulos, V., and Singer, H.: The Role of Transient Ion Foreshock Phenomena in Driving Pc5 ULF Wave Activity, J. Geophys. Res., 118, 299–312,, 2013. a

Hoshino, M. and Terasawa, T.: Numerical Study of the Upstream Wave Excitation Mechanism. I-Nonlinear Phase Bunching of Beam Ions, J. Geophys. Res., 90, 57–64, 1985. a

Ipavich, F. M., Galvin, A. B., Gloeckler, G., Scholer, M., and Hovestadt, D.: A Statistical Survey of Ions Observed Upstream of the Earth's Bow Shock – Energy Spectra, Composition, and Spatial Variation, (Upstream Wave and Particle Workshop, 86, 4337–4342,, 1981. a

Kempf, Y., Pokhotelov, D., Gutynska, O., Wilson, III, L. B., Walsh, B. M., von Alfthan, S., Hannuksela, O., Sibeck, D. G., and Palmroth, M.: Ion Distributions in the Earth's Foreshock: Hybrid-Vlasov Simulation and THEMIS Observations, J. Geophys. Res., 120, 3684–3701,, 2015. a, b

Kucharek, H.: On the Physics of Collisionless Shocks: Cluster Investigations and Simulations, J. Atmos. Sol.-Terr. Phys., 70, 316–324,, 2008. a

Kucharek, H., Möbius, E., Scholer, M., Mouikis, C., Kistler, L. M., Horbury, T., Balogh, A., Réme, H., and Bosqued, J. M.: On the origin of field-aligned beams at the quasi-perpendicular bow shock: multi-spacecraft observations by Cluster, Ann. Geophys., 22, 2301–2308,, 2004. a, b, c, d

Lembège, B. and Savoini, P.: Nonstationarity of a Two-dimensional Quasiperpendicular Supercritical Collisionless Shock by Self-Reformation, Phys. Fluids, 4, 3533–3548, 1992. a, b

Lembège, B. and Savoini, P.: Formation of Reflected Electron Bursts by the Nonstationarity and Nonuniformity of a Collisionless Shock Front, J. Geophys. Res., 107, 1037,, 2002. a

Lembège, B., Giacalone, J., Scholer, M., Hada, T., Hoshino, M., Krasnosel'skikh, V., Kucharek, H., Savoini, P., and Terasawa, T.: Selected Problems in Collisionless-Shock Physics, Space Sci. Rev., 110, 161–226,, 2004. a, b

Marcowith, A., Bret, A., Bykov, A., Dieckman, M. E., Drury, L. O., Lembège, B., Lemoine, M., Morlino, G., Murphy, G., Pelletier, G., Plotnikov, I., Reville, B., Riquelme, M., Sironi, L., and Novo, A. S.: The Microphysics of Collisionless Shock Waves, Rep. Prog. Phys., 79, 046901,, 2016. a

Matsukiyo, S. and Scholer, M.: On Reformation of Quasi-Perpendicular Collisionless Shocks, Adv. Space Res., 38, 57–63,, 2006. a

Mazelle, C., Meziane, K., Lequéau, D., Wilber, M., Eastwood, J. P., Rème, H., Sauvaud, J. A., Bosqued, J. M., Dandouras, I., McCarthy, M., Kistler, L. M., Klecker, B., Korth, A., Bavassano-Cattaneo, M. B., Pallocchia, G., Lundin, R., and Balogh, A.: Production of Gyrating Ions from Nonlinear Wave-Particle Interaction Upstream from the Earth's Bow Shock: A Case Study from Cluster-CIS, Planet. Space Sci., 51, 785–795,, 2003. a

Mazelle, C., Meziane, K., and Wilber, M.: Field-aligned and gyrating ion beams in a planetary foreshock, AIP Conf. Proc., 781, 89–94,, 2005. a

Meziane, K.: AIP Conference Proceedings, in: THE PHYSICS OF COLLISIONLESS SHOCKS: 4th Annual IGPP International Astrophysics Conference, AIP, 116–122,, 2005. a

Möbius, E., Kucharek, H., Mouikis, C., Georgescu, E., Kistler, L. M., Popecki, M. A., Scholer, M., Bosqued, J. M., Rème, H., Carlson, C. W., Klecker, B., Korth, A., Parks, G. K., Sauvaud, J. C., Balsiger, H., Bavassano-Cattaneo, M.-B., Dandouras, I., DiLellis, A. M., Eliasson, L., Formisano, V., Horbury, T., Lennartsson, W., Lundin, R., McCarthy, M., McFadden, J. P., and Paschmann, G.: Observations of the spatial and temporal structure of field-aligned beam and gyrating ring distributions at the quasi-perpendicular bow shock with Cluster CIS, Ann. Geophys., 19, 1411–1420,, 2001. a

Oka, M., Terasawa, T., Saito, Y., and Mukai, T.: Field-Aligned Beam Observations at the Quasi-Perpendicular Bow Shock: Generation and Shock Angle Dependence, J. Geophys. Res., 110, 05101,, 2005. a

Otsuka, F., Matsukiyo, S., Kis, A., Nakanishi, K., and Hada, T.: Effect of Upstream ULF Waves on the Energetic Ion Diffusion at the Earth's Foreshock. I. Theory and Simulation, Astrophys. J., 853, 117,, 2018. a

Paschmann, G., Sckopke, N., and Asbridge, J.: Energization of Solar Wind Ions by Reflection from the Earth's Bow Shock, J. Geophys. Res., 85, 4689–4693, 1980. a, b

Paschmann, G., Sckopke, N., Papamastorakis, I., Asbridge, J. R., Bame, S. J., and Gosling, J. T.: Characteristics of Reflected and Diffuse Ions Upstream from the Earth's Bow Shock, J. Geophys. Res., 86, 4355–4364,, 1981. a

Savoini, P.: Data for self consistent test particles, available at:, last access: 13 November 2020. a

Savoini, P. and Lembège, B.: Electron Dynamics in Two and One Dimensional Oblique Supercritical Collisionless Magnetosonic Shocks, J. Geophys. Res., 99, 6609–6635, 1994. a

Savoini, P. and Lembège, B.: Two-Dimensional Simulations of a Curved Shock: Self-Consistent Formation of the Electron Foreshock, J. Geophys. Res., 106, 12975,, 2001. a, b

Savoini, P. and Lembège, B.: Production of Nongyrotropic and Gyrotropic Backstreaming Ion Distributions in the Quasi-Perpendicular Ion Foreshock Region, J. Geophys. Res., 120, 7154–7171,, 2015. a, b, c, d, e, f, g, h, i, j, k, l, m, n, o, p, q, r, s, t, u, v, w, x

Savoini, P., Lembège, B., and Stienlet, J.: On the Origin of the Quasi-Perpendicular Ion Foreshock: Full-Particle Simulations, J. Geophys. Res., 118, 1132–1145,, 2013. a, b, c, d, e, f

Scholer, M., Shinohara, I., and Matsukiyo, S.: Quasi-Perpendicular Shocks: Length Scale of the Cross-Shock Potential, Shock Reformation, and Implication for Shock Surfing, J. Geophys. Res., 108, 1014,, 2003. a

Schwartz, S. J. and Burgess, D.: On the Theoretical/Observational Comparison of Field-Aligned Ion Beams in the Earth's Foreshock, J. Geophys. Res., 89, 2381–2384,, 1984. a

Schwartz, S. J., Thomsen, M. F., and Gosling, J. T.: Ions Upstream of the Earth's Bow Shock – A Theoretical Comparison of Alternative Source Populations, J. Geophys. Res., 88, 2039–2047,, 1983. a, b

Sonnerup, B. U. Ö.: Acceleration of Particles Reflected at a Shock Front, J. Geophys. Res., 74, 1301–1304,, 1969. a

Tanaka, M., Goodrich, C. C., Winske, D., and Papadopoulos, K.: A Source of the Backstreaming Ion Beams in the Foreshock Region, J. Geophys. Res., 88, 3046–3054,, 1983. a

Thomsen, M., Gosling, J., and Bame, S.: Gyrating Ions and Large-Amplitude Monochromatic MHD Waves Upstream of the Earth's Bow Shock, Jo. Geophys. Res., 90, 267–273, 1985. a

Thomsen, M. F., Gosling, J. T., and Schwartz, S. J.: Observational Evidence on the Origin of Ions Upstream of the Earth's Bow Shock, J. Geophys. Res., 88, 7843–7852,, 1983. a

Tsurutani, B. T. and Rodriguez, P.: Upstream Waves and Particles – An Overview of ISEE Results, J. Geophys. Res., 86, 4317–4324, 1981. a

Turner, D., Angelopoulos, V., Wilson, L., Hietala, H., Omidi, N., and Masters, A.: Particle Acceleration during Interactions between Transient Ion Foreshock Phenomena and Earth's Bow Shock, Geophys. Res. Abstr., EGU2014-2276, EGU General Assembly 2014, Vienna, Austria, 2014.  a

Webb, G., Axford, W., and Terasawa, T.: On the Drift Mechanism for Energetic Charged Particles at Shocks, Astrophys. J., 270, 537–553, 1983. a

Short summary
Numerical simulations have been used to investigate some acceleration mechanisms in order to explain the origin of the energized back-streaming ions observed at the Earth's bow shock. This paper used test particles in two different configurations with self-consistent and fixed shock front profiles. The comparison of these two configurations allows us to analyze, in detail, the impact of the shock front nonstationarity and the role of the built-up electric field in the acceleration process.