2.1 Simulations

We investigate the Earth's foreshock region using Vlasiator (Palmroth et al., 2018), a hybrid-Vlasov simulation capable of describing ion kinetics whilst encompassing global scales. Vlasiator solves the Vlasov equation for grid-discretized particle distribution functions, with closure provided by Ohm's law augmented by the Hall term. Electrons are considered a charge-neutralizing massless fluid with no electron pressure gradient term. We investigate the foreshock using a 2D–3V simulation, with 3D moments and velocity distribution functions but a 2D spatial domain. The geocentric solar ecliptic (GSE) simulation domain is $X\in (-\mathrm{48.66}{R}_{\mathrm{E}};\mathrm{64.35}{R}_{\mathrm{E}})$ and $Z\in (-\mathrm{59.65}{R}_{\mathrm{E}};\mathrm{39.24}{R}_{\mathrm{E}})$ and a single cell width in the Y direction. The spatial resolution is 300 km (1.3 times the solar wind ion inertial length) and the velocity–space resolution for both ion species (protons and alpha particles) is 30 km s⁻¹.

Our simulation is set to have solar wind values of β=0.7, M_ms=5.6, and M_A=6.9. We initialize the simulation with a solar wind of $n_{\mathrm{p}}=\mathrm{1}{\mathrm{cm}}^{-\mathrm{3}}$ and $n_{\mathit{\alpha }}={\mathrm{10}}^{-\mathrm{2}}{\mathrm{cm}}^{-\mathrm{3}}$. Due to the mass ratio, we set T_p=0.5 MK, and T_α=1.0 MK. The solar wind speed is set to $u_{\text{sw}}=\mathrm{750}{\mathrm{km}\mathrm{s}}^{-\mathrm{1}}$ in the $-{\widehat{e}}_x$ direction, simulating fast solar wind conditions and ensuring efficient simulation initialization. Despite the use of fast solar wind conditions, the Alfénic Mach number (7) and plasma beta (0.7) are typical for a variety of solar wind conditions, ensuring the validity of the initial conditions.

We set an interplanetary magnetic field (IMF) of B_x=3.54 nT and $B_z=-\mathrm{3.54}\mathrm{nT}$, resulting in a 45^∘ cone angle. Although the simulation plane is meridional, the foreshock dynamics are comparable with an ecliptical Parker spiral setup. The Earth's magnetic dipole is a ${\widehat{e}}_z$-aligned line dipole, resulting in a realistic magnetopause standoff distance (Daldorff et al., 2014). The simulation has an inner boundary at $\mathrm{3}\times {\mathrm{10}}^{\mathrm{4}}\mathrm{km}\approx \mathrm{4.7}{R}_{\mathrm{E}}$, modelled as a perfectly conducting ionosphere. Thus, the run is nearly identical to the simulation presented in Blanco-Cano et al. (2018), with the addition of alpha particles as an independent, self-consistent species. Our alpha particle density is set to only 1 % of the solar wind, so mass loading and effects on ULF wave properties are expected to be small. In order to constrain memory usage, we set the minimum stored phase–space densities (as explained in von Alfthan et al., 2014) to $f_{\text{min,p}}={\mathrm{10}}^{-\mathrm{15}}{\mathrm{s}}^{\mathrm{3}}{\mathrm{m}}^{-\mathrm{6}}$ and $f_{\text{min,}\mathit{\alpha }}={\mathrm{10}}^{-\mathrm{17}}{\mathrm{s}}^{\mathrm{3}}{\mathrm{m}}^{-\mathrm{6}}$.

Additionally, in subsection 3.1 we compare our main simulation to an equatorial Vlasiator simulation with a spatial resolution of 228 km and solar wind values of β=2.3, M_ms=5.9, M_A=10, $n_{\mathrm{p}}=\mathrm{3.3}{\mathrm{cm}}^{-\mathrm{3}}$, $n_{\mathit{\alpha }}=\mathrm{3.3}\times {\mathrm{10}}^{-\mathrm{2}}{\mathrm{cm}}^{-\mathrm{3}}$, and $u_{\text{sw,}x}=-\mathrm{600}\mathrm{km}{\mathrm{s}}^{-\mathrm{1}}$. The IMF is set to 5 nT, with a 5^∘ cone angle.

2.2 Observations

In this study, we also analyse observations from the MMS mission (Burch et al., 2016) in the Earth's foreshock. We use data from three different instruments, namely the fluxgate magnetometer (FGM; Russell et al., 2016), the Hot Plasma Composition Analyzer (HPCA; Young et al., 2016), and the dual ion spectrometers (DIS; Pollock et al., 2016), which is part of the fast plasma investigation (FPI) instrument. For all instruments, we use survey-mode data only. The FGM produces magnetic field measurements with 16 s⁻¹ time resolution in fast survey mode.

Energy–time spectrograms from HPCA measurements are used to identify whether the spacecraft are located in the foreshock or in the solar wind. We also calculate partial densities for the different ion species, following the procedure described in the HPCA science algorithms and user manual available on the MMS Science Data Center web page (https://lasp.colorado.edu/mms/sdc/public/datasets/hpca/, last access: 1 October 2020). HPCA data are made available in survey or burst modes. Data are collected in 0.5 spin (10 s) time resolution at 64 energies, 16 azimuth angles, and 16 elevation angles for five different species (H⁺, He⁺, He²⁺, O⁺, and O²⁺). In the survey-mode data used here, the data are reduced in energy and angles to 16 energies, eight azimuths, and eight elevation angles. The total energy range of the data is between $\sim \mathrm{1}\mathrm{eV}/e$ and 40 keV∕e, with a mass resolution of $M/\mathrm{\Delta }M\sim \mathrm{8}$.

Finally, we calculate the ion partial density (without species separation), using measurements from a FPI–DIS, as a means of comparison to the partial densities obtained from HPCA. The FPI produces burst sky maps which consist of ion count arrays of 32 energies × 32 azimuth angles × 16 polar angles that are accumulated every 150 ms. Then, 30 consecutive DIS burst sky maps are summed in order to produce the survey-mode sky maps with a time resolution of 4.5 s. The energy range of the DIS is between 10 eV∕e and 30 keV∕e.

3.1 Foreshock edge

As shown in Fig. 1, the ion and ULF foreshocks are not identical in extent. The ULF foreshock edge visible in both panels connects to the bow shock at $Z\sim -\mathrm{5}{R}_{\mathrm{E}}$, whereas the ion foreshock edge intersects the bow shock already at $Z\sim +\mathrm{5}{R}_{\mathrm{E}}$ (Fig. 1b). At the bottom of the figure, at $Z=-\mathrm{50}{R}_{\mathrm{E}}$, we see the ion foreshock extending up to X=40 R_E, whereas the ULF foreshock only extends to a position ∼10 R_E further downstream. We also see that, in Fig. 1b, the ratio of suprathermal alphas to protons in the foreshock shows significant deviations from the incoming solar wind ratio of 1 ∕ 100. Throughout most of the deep foreshock, the $n_{\mathit{\alpha }\text{,st}}{n}_{\text{p,st}}^{-\mathrm{1}}\times {\mathrm{10}}^{\mathrm{2}}$ ratio, shown in Fig. 1b, tends to values ≳2, whereas at the foreshock edge, it falls below 0.2. This abundance of He²⁺ within the deep foreshock is likely a computational artefact, with H ions being efficiently scattered into the diffuse population, which is not tracked, whereas He²⁺ remains in non-gyrotropic partial rings and clumps. Right at the edge of the foreshock we see spatially periodic structures with, again, large alpha-to-proton ratios, likely caused by bursty reflection at the quasi-perpendicular shock and the alpha particles having larger gyroradii, thus gyrating further into the upstream. The ULF and ion dynamics of the foreshock edge are further examined in Fig. 2.

Figure 2(a) Magnetic field out-of-plane component B_y indicating the ULF foreshock extent. Contours indicate proton (black) and helium (green) suprathermal densities at 0.1 % (dotted contours) and 0.5 % (solid contours) of solar wind values. Four thick lines indicate the cut-throughs across the foreshock boundary. (b–e) Profiles across the foreshock at the positions shown in (a). Profiles are shown for suprathermal proton density (black solid line), suprathermal helium density (green solid line), and magnetic field magnitude (blue solid line). The horizontal dashed lines indicate 0.5 % (dashed line) of the solar wind density.

Download

Figure 2a shows an excerpt from Fig. 1a, featuring a portion of the foreshock edge. The colour map is again the out-of-plane magnetic field component, overlaid with contours of proton (black) and helium (green) suprathermal densities. We focus on the contours upstream of the ULF foreshock. The solid contours are drawn where the suprathermal density is 0.5 % of the species' solar wind density and the dotted contours at 0.1 %, respectively. The thick grey lines indicate four cut-throughs across the foreshock edge. The magnetic field strength and the suprathermal proton and helium densities along these cuts are plotted in Fig. 2b–e. As shown in these plots, there is variation in the profiles of ions across the foreshock edge, moving into the foreshock from right to left. Close to the bow shock (Fig. 2b–c), there is a somewhat rapid increase in ion densities near 10 R_E, followed by a gradual increase over the following 1–2 R_E and finally plateauing values. Further out (Fig. 2d and e), we see a more gradual increase in suprathermal ion densities over several R_E and a less clear plateau. The suprathermal ion density threshold at 0.5 % of the species solar wind density is shown as a horizontal dashed line.

We note that non-thermal particles are found several R_E upstream of the foreshock waves, consistent with previous works showing that no ULF wave activity is observed in conjunction with the region closest to the foreshock edge where field-aligned beams are expected (and found). More importantly, we find significant amounts of suprathermal helium throughout the field-aligned beam region and into the ULF foreshock. As shown by the black and green contours in Fig. 2a, the helium foreshock edge is located equal to or slightly downstream of the proton foreshock edge. This shift is more pronounced in the dotted contours and increases when moving further away from the bow shock. At about 30 R_E from the shock, the difference is of the order of 1 R_E. The shift of the foreshock edge is also noticeable in the line profiles (Fig. 2b–e). This suggests that the proton foreshock is slightly more extended than the helium foreshock, and measurements made at the very edge of the foreshock can suggest significantly lower helium fractions despite the helium abundances rising to proton-comparable levels a few R_E deeper within the foreshock. We also point out that the plateau of alpha particles inside the foreshock edge has more fluctuations to it than the plateau of protons.

We also note that the suprathermal ion density contours in Fig. 2 and the time-averaged ratios seen in Fig. 1 are not smooth, having a wavy shape instead. Supplementary Video A shows the wavy shape is due to intermittent reflection of ions at the bow shock, with bursts of ions propagating away from the shock along the field lines. This periodic and intermittent enhanced reflection of particles may be due to mesoscale reformation of the bow shock (Battarbee et al., 2020e) and can cause further discrepancies and variation in field-aligned beam densities.

Figure 3MMS3 crossing from the solar wind into the foreshock region on 30 December 2018. (a) The MMS position relative to the bow shock just before the foreshock crossing in a plane that contains B. The solid blue line shows the bow shock position in the plane containing MMS and the dashed line in the plane containing Earth. The grey lines are the magnetic field lines. (b) The survey-mode measurements of the magnetic field magnitude and components, showing the spacecraft entering the ULF foreshock. (c) The number densities for H⁺ and He²⁺ suprathermal ions measured by HPCA with spacecraft frame energy cutoffs indicated in the legend, with the helium density multiplied by 10 to ease density evaluation. The FPI ion number density is shown for comparison. (d, e) The HPCA ion energy spectrograms for protons and helium, with the dashed line indicating the respective suprathermal ion cutoff energies.

Download

Figure 4MMS3 crossing from the solar wind into the foreshock on 18 November 2018. The format is the same as in Fig. 3 but with He²⁺ multiplied by 25.

Download

We now compare our numerical results with observations from the MMS spacecraft at the foreshock edge. Figures 3 and 4 show two time intervals during which the MMS3 spacecraft crossed from the solar wind into the foreshock region. Figures 3a and 4a show the spacecraft position relative to the bow shock model, which is scaled with the solar wind dynamic pressure (Farris et al., 1991; Schwartz, 1998). The plane is chosen so that B is in the $X_{\text{GSE}}-Y^{\prime }$ plane and that Y^′ is in the positive Y_GSE direction. Figures 3b–e and 4b–e, from top to bottom, show the IMF magnitude and components, the suprathermal densities of H⁺ (black; HPCA), He²⁺ (green; HPCA), and all ions (red; FPI) and proton and He²⁺ energy–time spectrograms from the HPCA instrument. The He²⁺ suprathermal densities have been multiplied by 10 or 25 to ease comparison of suprathermal ion density gradients. The lowest energy used in the calculation of the suprathermal (partial) densities for each species is also stated on the panel. This energy is about 4 times higher for alphas compared to the protons, due to larger mass of He²⁺. We note that the energy–time spectrograms show energy per charge $E_i{Z}_i^{-\mathrm{1}}$, with the lowest energy included in the suprathermal population shown with the dashed line set to the values indicated in Figs. 3c and 4c. Also, only those suprathermal ions measured by HPCA (up to 40 keV∕e) are included in the densities, but we estimate the lost contribution due to higher energy ions to be small due to very small phase–space densities.

During both time intervals the transition between the foreshock and the undisturbed solar wind was not due to any IMF rotational discontinuity, as can be seen from the MMS magnetic field components. We also checked the solar wind measurements propagated to the bow shock nose from the OMNI database (King and Papitashvili, 2005) and found no sharp IMF change during those intervals. This means that the foreshock was undergoing gradual motion due to slow IMF rotations, so the spacecraft did not observe travelling foreshocks (Kajdič et al., 2017) or foreshock cavities (Schwartz et al., 2006; Billingham et al., 2008, 2011). This gradual motion of the foreshock region over the spacecraft location allows for magnetic field and suprathermal density profiles to be compared with those obtained from our simulations (Fig. 2).

Figures 3c and 4c show that the proton and He²⁺ suprathermal densities do not always behave in the same manner. In Fig. 3c, we can clearly see that the proton and He²⁺ suprathermal densities are well correlated across the foreshock edge, although H⁺ has a slightly stronger initial beam before 22:06 Coordinated Universal Time (UTC). Both increase with a similar relative gradient as the spacecraft crosses from the solar wind into the foreshock. In Fig. 4c, the suprathermal proton density starts to increase ∼3 min before the suprathermal He²⁺ density and remains low well into the foreshock. This suggests that the proton foreshock extends further outward than the He²⁺ foreshock. In other words, the proton foreshock extends to field lines connected to larger θ_Bn values than the He²⁺ foreshock.

In order to better understand the different foreshock edge crossing behaviour seen in Figs. 3 and 4, we compare the main Vlasiator simulation used in this study to an equatorial plane simulation with a quasi-radial IMF, as detailed in Sect. 2.1. In Fig. 5, we show an enlarged view of the foreshock edge of both runs, with colour maps indicating the out-of-plane component of the magnetic field and black and green contours indicating proton and helium suprathermal densities at 0.5 % (solid) and 0.1 % (dotted) inflow densities, averaged over 2 min. In comparing the Vlasiator plots with MMS observations, it is important to keep in mind that the Vlasiator definition for suprathermals includes reflected particles with simulation frame energies comparable with the solar wind bulk, whereas the MMS suprathermal density is defined with a markedly higher minimum energy threshold. Comparison of Fig. 5a–b indicates that simulations can reproduce the two different foreshock edge behaviours, with Fig. 5a (our main run) presenting a relatively rapid drop-off and Fig. 5b (the quasi-radial IMF comparison run) showing a much more gradual fall-off of suprathermal ion densities and a greater difference between the two ion species.

The IMF prior to the 30 December 2018 foreshock edge crossing had the IMF pointing in the dawn, anti-sunward direction with its cone angle $\sim \mathrm{55}{}^{\circ }$. Consequently, this was also the orientation of the foreshock. At the time, MMS3 was located near the nose of the Sun–Earth line, at approximately $(\mathrm{12.5},\mathrm{0.8},\mathrm{5.1})R_{\mathrm{E}}$ in GSE coordinates. This location, together with large IMF cone angle, means that the foreshock edge was crossed in a way that is qualitatively similar to that shown in Fig. 5a.

In the case of the 18 November 2018 event, the IMF cone angle was $\sim \mathrm{35}{}^{\circ }$ and pointing in the northern, anti-sunward direction prior to the foreshock encounter, and consequently, this was also the foreshock orientation. The GSE coordinates of MMS3 were $(\mathrm{9.8},\mathrm{11.3},\mathrm{5.7})R_{\mathrm{E}}$, i.e. far from the Sun–Earth line. Thus, this foreshock edge crossing is comparable with the one shown in Fig. 5b.

Figure 5An enlarged view of the foreshock edge, with the magnetic field out-of-plane component B_oop indicating the ULF foreshock extent. Contours indicate suprathermal ion densities, as in Fig. 2a, but averaged over 2 min. (a) The meridional plane 45^∘ IMF simulation used in the majority of this paper. (b) An equatorial plane 5^∘ IMF run for comparison. The IMF orientation appears to affect the suprathermal ion density gradient at the foreshock edge.

Download

3.2 Velocity distribution functions and their properties

We now examine the properties of ion velocity distribution functions (VDFs) at the edge of and within the foreshock. Figure 6 shows 2D projections of proton and helium velocity distribution functions in the solar wind frame at three positions close to the foreshock edge, extracted from our Vlasiator simulation at time 1000 s. In Fig. 6 (panels 1–2 and 4, subpanels $(\parallel ⟂)$), we show projections that have been generated by averaging the instantaneous VDFs over the v_B×V direction, whereas subpanels labelled $(⟂⟂)$ have been averaged over the v_B direction. We use two different colour scales to differentiate the phase–space densities of the two ion species, with ranges selected to account for the input solar wind abundance ratio. Figure 6 (panels 1–2 and 4) shows VDFs from virtual spacecraft locations (A, B, and C, respectively). The Fig. 6 $(\parallel ⟂)$ subpanels also show an ellipsoid located at the position in velocity space where particles would end up if they were specularly reflected from the solar wind population at the closest point of the bow shock. The shock location was determined according to a plasma compression criterion of n_p>2n_p,sw, and the shock normal direction was estimated to be equal to a vector pointing from the closest shock location to the virtual spacecraft location. We emphasize that this estimate is a rough one, to be improved upon in future studies, and does not account for ion propagation times, drifts, or the existence of an electron pressure gradient cross-shock potential. Figure 6 (panel 3) shows an overview of the foreshock region, indicating the locations of the spacecraft on top of a colour map depicting the temperature anisotropy calculated from the whole proton VDF. The dark orange region at the upstream edge of the foreshock, with parallel temperatures in excess of perpendicular temperatures, is indicative of the FAB region of the proton foreshock. Magnetic fields lines are depicted with black curves.

In panel 1 of Fig. 6, at the position of virtual spacecraft A located very close to the foreshock edge, we see very clear proton and helium FABs, though the helium beam appears to have more structure. Parallel velocities of both beams are between 1500 and 2000 km s⁻¹ in the solar wind frame or $\mathit{\gtrsim }{u}_{\text{sw}}=\mathrm{750}\mathrm{km}{\mathrm{s}}^{-\mathrm{1}}$ in the simulation frame, which translates to roughly 3–5 keV nuc⁻¹. We note that this is larger than the energy at which specularly reflected particles would be found, which suggests that these particles have experienced shock drift acceleration (SDA) at the quasi-perpendicular shock front. This indicates that the source of the FABs in panel 3 of Fig. 6 is located at roughly $X=\mathrm{16}{R}_{\mathrm{E}};Z=\mathrm{0}{R}_{\mathrm{E}}$. Panel 2 of Fig. 6 depicts VDFs at position B, at the boundary between the FAB region and the ULF foreshock. At this location, the ion beam seems to be transitioning into a more gyrating distribution, with a gyrophase-bunched signature (visible in the $(⟂⟂)$ panel at $v_{B\times V}>\mathrm{500}\mathrm{km}{\mathrm{s}}^{-\mathrm{1}}$) for both species but particularly for helium at $v_{B\times V}>\mathrm{1000}\mathrm{km}{\mathrm{s}}^{-\mathrm{1}}$. Parallel velocities have begun to decrease, extending down to 1000 km s⁻¹ in the solar wind frame of reference. A similar upstream ion velocity decrease was also seen in Battarbee et al. (2020e). Still, at this position, both helium and protons show a qualitatively similar shape to their distribution functions. Conversely, in panel 4 of Fig. 6, at position C, the proton VDF resembles a low-energy field-aligned beam, extending from 1000 to 1800 km s⁻¹ in v_B, whereas the helium distribution has very little trace of a beam, instead consisting of a broken up gyrating ion population. An animated version of Fig. 6 can be found as Supplementary Video B. We also note that, early in this video, at virtual spacecraft A, the helium VDF in particular shows what appears to be a gyrophase-bunched population but which remains stationary in velocity space. We deduct that it is in fact spatial sampling of the helium foreshock edge, visible due to the large gyroradius of alpha particles. We also point out that the proton beam energy found, in particular, in panel 1 of Fig. 6 does not extend to the tens of keV∕e found in some spacecraft observations, possibly a result of the sparse velocity space implementation in Vlasiator.

Figure 6Velocity distribution functions (VDFs) for protons and alpha particles and their locations in the foreshock. (3) Map of the foreshock, with three VDF positions indicated with capital letters. The colour indicates the proton temperature anisotropy, with magnetic field lines in black. Panels (1), (2), and (4) show sets of four projections of ion VDFs at virtual spacecraft locations in the solar wind frame. Coordinates are chosen to represent two positions at the foreshock boundary (A and B) and one just within the ULF foreshock (C). In each set, the subpanels are labelled as v_B versus $v_{B\times (B\times V)}$ ($\parallel ⟂$; left columns) or v_B×V versus $v_{B\times (B\times V)}$ ($⟂⟂$; right columns) and protons (top rows) or helium (bottom rows). The black and white circles estimate the areas of specular reflection.

Download

Figure 7 is similar to Fig. 6 but with virtual spacecraft locations chosen to represent regions further within the foreshock and also closer to the quasi-parallel bow shock. Panels 1 and 2 of Fig. 7 show spacecraft D and E, which are located close to the bow shock and depict mostly gyrating and intermediate ion populations. The helium populations of gyrating ions appear to have more structure to them. It is also noteworthy that both proton and helium $(\parallel ⟂)$ subpanels of panel 1 in Fig. 7 have what looks like a FAB propagating towards the shock at a parallel velocity greater than the solar wind speed, with velocities v_B<0. Although the estimate of the velocity space location of specularly reflected particles does not account for particle travel times or shock shape evolution, we find that the near-circular ellipsoids indicating potential specular reflection are usually found in VDF regions where there are enhancements, suggesting that specular reflection indeed plays a role in the generation of these populations. As the bow shock reforms as a mesoscale process with distinctly non-planar features, the direction of specular reflection varies, leading to the intermittent and gyrating partial rings seen in these VDFs. This mapping of specular reflection can also be seen in Supplementary Video C, which is an animated version of Fig. 7. Panel 4 of Fig. 7 shows proton and helium VDFs further within the deep foreshock and away from the bow shock. The proton population appears to be a gyrating ion population, and the helium population is a low-energy beam population, but viewing the VDFs at different times (see Supplementary Video C) shows that both ions usually resemble gyrating ion populations.

Figure 7Velocity distribution functions for protons and alpha particles and their locations in the foreshock. (3) Map of the foreshock with three VDF positions indicated with capital letters. The colour indicates proton temperature anisotropy, with magnetic field lines in black. Panels (1), (2), and (4) show sets of four projections of ion VDFs at virtual spacecraft locations in the solar wind frame. Coordinates are chosen to represent two positions close to the quasi-parallel bow shock (D and E) and one deep within the outer foreshock (F). In each set, subpanels are labelled as v_B versus $v_{B\times (B\times V)}$ ($\parallel ⟂$; left columns) or v_B×V versus $v_{B\times (B\times V)}$ ($⟂⟂$; right columns) and protons (top rows) or helium (bottom rows). The black and white circles estimate the areas of specular reflection.

Download

In light of the complex simulated VDF shapes, and in order to obtain a better understanding of global foreshock ion characteristics, we show, in Fig. 8, plots of global per-species temperature anisotropies and a measure of per-species non-gyrotropy. We apply the temperature anisotropy to the whole VDF, allowing us to identify regions where VDFs have FAB-like features showing up as anisotropy values much smaller than 1. The agyrotropy measure (Swisdak, 2016) is as follows:

where $P_{\parallel }=P_{\mathrm{11}}$ is the parallel pressure, $P_{⟂}=\mathrm{0.5}(P_{\mathrm{22}}+P_{\mathrm{33})}$ is the perpendicular pressure, and P₁₂, P₁₃, and P₂₃ are off-diagonal pressure tensor components that can be used to evaluate the complexity of the VDF and the role of gyrating ions. For a completely gyrotropic distribution, Q_ag=0 and Q_ag=1 would signify a maximal deviation from gyrotropy. Figure 8a–b show agyrotropies for protons and helium, respectively, and Fig. 8c–d show temperature anisotropies. The panels in Fig. 8 also show out-of-plane magnetic field contours at a level of ±0.2 nT, indicating the extent and structure of the ULF foreshock. We find that the FAB region of protons is clearly visible in Fig. 8c at the edge of the ion foreshock as a dark orange band ($T_{⟂,\mathrm{p}}{T}_{\parallel ,\mathrm{p}}^{-\mathrm{1}}\approx \mathrm{0.2}$). In Fig. 8a we see the same structure at the foreshock edge. It is visible as a band of medium green blobs (Q_ag,p∼0.01) close to the shock nose, transitioning to paler green thin streaks (Q_ag,p∼0.001) away from the shock.

Figure 8Foreshock characteristics for proton and helium at 1100 s. (a, c) Protons and (b, d) helium are shown. (a, b) Agyrotropy Q_ag (see Eq. 1), where a value of 0 corresponds with perfect gyrotropy. (c, d) Temperature anisotropy $T_{⟂}{T}_{\parallel }^{-\mathrm{1}}$ calculated for the total distribution function, including both solar wind and suprathermal parts. Black contours show the out-of-plane magnetic field B_y fluctuations at a level of ±0.2 nT, indicating the extent of the ULF foreshock region.

Download

For helium, in Fig. 8d, this FAB region at the edge of the foreshock is also clearly visible, with very low anisotropy values. Interestingly, helium also shows a parallel pressure signature (dark orange low anisotropy values) deeper in the foreshock. This region coincides with a weakening of the ULF foreshock, visible in the destructuring of wave fronts (black contours) in the vicinity of ($X=\mathrm{10}{R}_{\mathrm{E}};Z=-\mathrm{40}{R}_{\mathrm{E}}$). Referring to Fig. 1b, we see that this band of low-energy, FAB-type alpha particles coincides with an increase in the measured helium fraction, which continues downstream from that point.

Figure 8a–b show, in particular for protons, an enhancement in agyrotropy at the boundary between the FAB region and the ULF foreshock. This dark green region is a signature of gyrating and gyrophase-bunched ions at this boundary, in agreement with Meziane et al. (2004), Mazelle et al. (2007), and Andrés et al. (2015). This effect can also be seen at virtual spacecraft location B in Fig. 6 (panel 2; ($p⟂⟂$) and ($\mathit{\alpha }⟂⟂$) subpanels), as a gyrophase-bunched extension of the ion distribution. We note that this increase in agyrotropy and, thus, gyrating ions matches the ULF foreshock and previous studies well, down to about $Z=-\mathrm{25}{R}_{\mathrm{E}}$, but the boundary becomes less well defined further out, away from the shock. We also note that Fig. 8b shows darkened bands of increased agyrotropy on both the upstream and downstream edges of the inner-foreshock-heightened parallel pressure alpha particle band, which suggests there might be similar gyrophase sampling taking place as for the foreshock FAB beam proper.

For both protons and helium, we see striped enhancements of agyrotropy right at the outer ion foreshock boundary, indicative of spatial sampling of gyrating ion beams right at the outermost foreshock edge. The gyration of these ion beams, accelerated at the quasi-perpendicular shock front and made non-uniform by the rippling of the shock front, is particularly visible in Supplementary Video D, which is an animated version of Fig. 8. Finally, we note that large regions of the foreshock close to the quasi-parallel bow shock show signatures of temperature anisotropies ≳1 and enhanced agyrotropies, which are likely a signature of specular reflection of ions at the quasi-parallel bow shock.

3.3 Foreshock waves

Figure 9Time series and wavelets at virtual spacecraft positions C (a–h) and D (i–p) at the edge of the ULF foreshock. Panels (a–d) and (i–l) show proton (helium) quantities in black (green). (a, i) Total number densities. (b, j) Suprathermal number densities on a logarithmic scale. (a, b, i, and j) The helium number densities have been multiplied by 100. (c, k) Temperature anisotropies. (d, l) Agyrotropy Q_ag measured as a 5 s running average. (e, m) Magnetic field magnitudes $\left|B\right|$. (f, n) Magnetic field out-of-plane components B_y. (g, o) Wavelet power spectra of $\left|B\right|$ fluctuations. (h, p) Wavelet power spectra of B_y fluctuations. The black contours show the 95 % confidence level, and the cross-hatched regions denote the cone of influence. The white dashed lines show the expected frequency of B_y fluctuations (see text).

Download

Figure 10Time series and wavelets at virtual spacecraft positions E (a–h) and F (i–p) in the ULF foreshock. Panels (a–d) and (i–l) show proton (helium) quantities in black (green). (a, i) Total number densities. (b, j) Suprathermal number densities on a logarithmic scale. (a, b, i, and j) The helium number densities have been multiplied by 100. (c, k) Temperature anisotropies. (d, l) Agyrotropy Q_ag measured as a 5 s running average. (e, m) Magnetic field magnitudes $\left|B\right|$. (f, n) Magnetic field out-of-plane components B_y. (g, o) Wavelet power spectra of $\left|B\right|$ fluctuations. (h, p) Wavelet power spectra of B_y fluctuations. The black contours show the 95 % confidence level, and the cross-hatched regions denote the cone of influence. The white dashed lines show the expected frequency of B_y fluctuations (see text).

Download

Figures 9 and 10 show measurements from virtual spacecraft placed at locations C, D, E, and F in the foreshock (see Figs. 6 and 7). The top row in each plot displays the total number density n_tot for protons and helium, and the second row displays their suprathermal number densities n_st. In both rows, the helium number densities have been multiplied by a factor of 100 to ease the comparison with proton number densities. Note that the suprathermal number densities are shown on a logarithmic scale, since the values vary significantly. In the third row, temperature anisotropy $T_{\perp }{T}_{\parallel }^{-\mathrm{1}}$ is displayed, and the fourth row displays the agyrotropy measure Q_ag. The time series of Q_ag have been smoothed out with a 5 s running average to help readability due to a large number of high-frequency fluctuations. The fifth and sixth rows display the total magnetic field $\left|B\right|$ and its out-of-plane component B_y, respectively.

The fluctuations of H⁺ and He²⁺ total densities follow each other quite closely at all locations, though the amplitude of the He²⁺ density oscillations is larger than that of protons at point D and E. These larger He²⁺ density variations seem to be well correlated with the fluctuations of the He²⁺ suprathermal density at position D (Fig. 9j) but not so much at point E (Fig. 10b). At all locations analysed with these virtual spacecrafts, the suprathermal ion ratio $n_{\mathit{\alpha }\text{,st}}{n}_{\text{p,st}}^{-\mathrm{1}}$ is larger than the solar wind ion ratio $n_{\mathit{\alpha }\text{,sw}}{n}_{\text{p,sw}}^{-\mathrm{1}}$, as shown already by Fig. 1b. The agyrotropy is also more pronounced for He²⁺ than for H⁺, as illustrated in Fig. 8a–b, and with the VDF discussed in the previous section.

The two bottom rows of Figs. 9 and 10 (panels g, h, o and p) display the wavelet power spectra of $\left|B\right|$ and B_y, in order to investigate foreshock wave activity. The wavelet transform (Torrence and Compo, 1998) is calculated using the Morlet wavelet. The black contours in the power spectra show the 95 % confidence level, and the cross-hatched regions bordering the spectra depict the cone of influence, where edge effects arising due to the time series' endpoints become important. The horizontal dashed line on the B_y wavelet power spectra is drawn at 50 s, which is close to the expected spacecraft frame period of foreshock fast magnetosonic waves for these upstream conditions, according to empirical models (Le and Russell, 1996; Takahashi et al., 1984). At all four positions, an enhancement in the B_y component power spectra can be seen in the vicinity of this period, which is in agreement with previous works showing that fast magnetosonic waves permeate the foreshock in Vlasiator simulations as in spacecraft observations (Palmroth et al., 2015; Turc et al., 2018). We note that none of the virtual spacecraft selected here display the typical quasi-monochromatic foreshock waves shown in these previous studies, due to their relative proximity to the bow shock (≲7 R_E). In the vicinity of the shock, the wave activity is more complex due to nuanced interactions between the waves and the diffuse ion population (Greenstadt et al., 1995; Turc et al., 2018).

At point F we find weaker wave activity, despite its location deep within the foreshock. This is most likely due to the low density of suprathermal gyrating or beam particles in this part of the foreshock (see Fig. 10j). This results in a lower wave growth rate as this parameter depends on the beam density for the beam–beam instabilities at play in the foreshock (Gary, 1993). This weaker wave activity is accompanied by a lower temperature anisotropy for the He²⁺ ions, due to the suprathermal He²⁺ population being in the form of low-energy, field-aligned beams in this region (see Fig. 7, panel 4). Comparing the positions of our virtual spacecraft with Fig. 8 and the contours indicating well-structured or more broken up ULF wave fronts, we see that the weakest wave power (point F) is seen at the most broken up position, and the strongest power (point D) is found at a position of well-structured wave fronts.

According to previous works, the plasma rest frame frequency of foreshock fast magnetosonic waves generated by proton beams is of the order of 10 % of the proton gyrofrequency, probably due to the cyclotron resonance which gives rise to the waves (Hoppe and Russell, 1982; Eastwood et al., 2005; Wilson, 2016). He²⁺ ions have a gyrofrequency that is half that of protons. He²⁺ ion beams could, thus, generate fast magnetosonic waves at a period of ∼100 s with the upstream parameters used in our simulation. The wavelet power spectra in Figs. 9 and 10 show enhanced wave power around ∼100 s, but this part of the spectra is mostly inside the cone of influence of the wavelet transform, making it difficult to draw firm conclusions. Moreover, the analysis of the foreshock wave properties in an identical Vlasiator run without a helium population (not shown) reveals similar enhancements of the wave power at ∼100 s. It is therefore unlikely that these fluctuations are due to helium beam instabilities. Thus, despite a 1 % solar wind helium content providing non-negligible mass loading, we find the helium component to not have a significant impact on foreshock wave populations.