Plasma Properties in the Earth's Magnetosheath Near the Subsolar Magnetopause: Implications for Geocoronal Density Estimates

Combined in situ ion measurements and remote sensing of energetic neutral atoms are used to determine the geocoronal Hydrogen density at large (∼10 RE) distances from the Earth. This method for determining the geocoronal density requires global magnetospheric modeling. Observations in the Earth's subsolar magnetosheath from the Magnetospheric Multiscale mission are used to determine the accuracy of using global models to predict the geocoronal density. On average, gas dynamic and magnetohydrodynamic (MHD) models and observations are in reasonable agreement, with differences <25%. In addition, the MHD model subsolar magnetopause is about 0.5 RE sunward of the observed location. However, variations around averages are large (up to a factor of 2), indicating that global models introduce relatively large uncertainties in geocoronal density estimates. Finally, the critical ion flux in the Interstellar Boundary Explorer IBEX‐Hi energy range is often minimally affected by fluctuations of a factor of 2 in the density.

• Magnetohydrodynamic modeling tends to overestimate the subsolar magnetopause standoff distance and underestimate the magnetosheath density • The subsolar magnetopause moves constantly, even under quasi-steady solar wind pressure, with displacements at least as large as ∼1 R E • Models for magnetosheath line-of-sight ion fluxes introduce uncertainties of up to a factor of 2 in the geocoronal density estimate at 10 R E

Supporting Information:
Supporting Information may be found in the online version of this article.

Introduction to the Earth's Geocorona at Large Distances From the Planet
The Earth's hydrogen geocorona has been observed as far as the moon (Baliukn et al., 2019).This extensive, tenuous neutral atmosphere has important implications for magnetospheric dynamics.For example, charge exchange with the geocorona is a major contributor to plasma ring current loss in the magnetosphere (e.g., Kistler et al., 1989).
Hydrogen atoms in the geocorona scatter solar Lyman-Alpha and line-of-sight measurement of these Far Ultraviolet (FUV) emissions are used to construct and validate geocoronal models (e.g., Bailey & Gruntman, 2013;Østgaard et al., 2003;Rairden et al., 1986;Zoennchen et al., 2011Zoennchen et al., , 2013) ) in particular, its density out to about 8 R E from the Earth (Zoennchen et al., 2017).At radial distances less than about 8 R E , this technique provides the most accurate determination of the geocoronal density.Beyond 8 R E , Lyman-Alpha emissions are weak, signals and backgrounds are comparable, and the resulting models have large uncertainties.
In contrast, ENA imaging of this region has substantially lower backgrounds.Because ENA emission is the product of charge exchange between plasma ions and the geocorona, geocoronal density is derived from ENA imaging with knowledge (through models and/or in situ measurements) of the source ion population.Equation 1 relates the line-of-sight integrated ion flux (J Ion (E,l)) along the l direction with the ENA flux (J ENA (E)) at a given energy through the energy-dependent charge exchange cross section σ(E) and the geocorona density n H (l).
Determining the geocoronal density from Equation 1 requires knowledge of the ion flux along the line-of-sight as well as how the geocoronal density varies with radial distance from the Earth.Because both FUV and ENA line-of-sight emissions are optically thin, modeling must be applied for proper interpretation of the measurements.For ENA imaging, a single in situ measurement along a ENA line-of-sight measurement provides crucial context for derivation of the geocoronal density.
A modified version of Equation 1 was solved for n H at the subsolar magnetopause (i.e., at a geocentric distance of 10 R E ) using simultaneous ENA and in situ ion observations (Fuselier et al., 2010), yielding n H (r = 10 R E ) = ∼5-10 cm −3 .Recently, Equation 1 was solved directly using simultaneous ENA and in situ observations combined with a gasdynamic model of the magnetosheath density (Fuselier et al., 2020), yielding n H (r = 10 R E ) = 11-17 ± 3 cm −3 .
This recent determination of n H was challenged (Sibeck et al., 2021).Using predictions from a global magnetohydrodynamic (MHD) simulation that incorporated the observed upstream solar wind as a boundary condition, it was argued that the neutral hydrogen geocoronal density was a factor of 4-5 times higher than that reported in Fuselier et al. (2020) (see Table S1 in Supporting Information S1).Half of this discrepancy between the two estimates of n H results from a factor of 2 discrepancy between the observed ion densities at the subsolar magnetopause reported by Fuselier et al. (2020) and the predicted ion densities at that location from an MHD model reported by Sibeck et al. (2021).Since any determination of n H (l) from Equation 1 requires a global model of the magnetosheath, it is important to determine how well these models reproduce subsolar magnetopause observations.This paper investigates some of the observed properties of the subsolar magnetopause that are important for determining n H (l) and compares these observed properties with gasdynamic and MHD model predictions.

Observations
Plasma and magnetic field observations are from the MMS spacecraft (Burch et al., 2016) and remote ENA observations are from the Interstellar Boundary Explorer (IBEX) mission (McComas et al., 2009).Plasma moments and ion fluxes are from the MMS Fast Plasma Investigation (Pollock et al., 2016) and the Hot Plasma Composition Analyzer (Young et al., 2016) and magnetic field observations are from the Magnetometers (Russell et al., 2016).Since properties of the subsolar magnetosheath on scales of 1,000's of km, that is, much larger than the ∼10s of km spacecraft separation, are considered here, data from a single spacecraft, namely MMS3, are sufficient.Upstream conditions are determined from solar wind data convected to the subsolar magnetopause.
The first ∼2 years of MMS observations (September 2015-March 2017) were surveyed for magnetopause encounters within 1 R E of the subsolar point.The MMS orbit was designed to pass near the subsolar point for 10.1029/2023GL105553 3 of 9 these two sweeps of the dayside (see Fuselier et al., 2016).After 2017, the orbit evolved and the spacecraft rarely cross the magnetopause within 1 R E of the subsolar point.
Nineteen events were identified, with multiple, inbound and outbound, complete magnetopause crossings in an event (see Table S2 in Supporting Information S1).Crossings within minutes of each other were not ruled out and event #19 includes the magnetopause encounter in Fuselier et al. (2020) and Sibeck et al. (2021).
Figure 1 shows a pair of crossings on 18 December 2015 (event #9).Top to bottom are the ion and electron omnidirectional fluxes, ion and electron densities, ion velocities, solar wind dynamic pressures convected to the magnetopause from an upstream solar wind monitor, and spacecraft distances to the magnetopause computed from the bulk ion velocity in the magnetosphere (see below).MMS3 was in the subsolar magnetosheath at the beginning of the interval.This region is characterized by high fluxes of several hundred eV to ∼1 keV ions and several tens of eV electrons, high densities, and low bulk flow velocities.The spacecraft crossed the magnetopause at 09:32 UT at a standoff distance of 10.8 R E .It remained in the magnetosphere for almost 8 min and simultaneously observed two populations of ions and electrons at high and low energies, before crossing the magnetopause at ∼09:40 UT and returning to the subsolar magnetosheath for the remainder of the time interval., solar wind dynamic pressures convected to the magnetopause, and distances from MMS to the magnetopause when the spacecraft was in the magnetosphere.The ion density at the inbound magnetopause crossing (09:32 UT, left blue arrow) is almost a factor of 2 lower than that at the outbound crossing (09:39:30 UT, right blue arrow).Both magnetopause crossings probably occurred because of the relatively small changes in the solar wind dynamic pressure.The bottom panel shows that, while the spacecraft was in the magnetosphere, the magnetopause moved almost 1 R E sunward by 0938 UT, before moving rapidly earthward and crossing over the spacecraft at 09:39:30 UT.These large changes in the magnetopause location are not reflected in the changes in the solar wind dynamic pressure.
The magnetopause crossings were likely the result of an ∼10% decrease in the solar wind dynamic pressure before 09:32 UT and a similar increase in pressure before 09:40 UT.In addition to this small change in the dynamic pressure, the IMF clock angle changed from 170° at 09:32 UT to 220° at 09:40 UT with a fairly large change in the cone angle as well (not shown).In contrast to these relatively small changes in upstream conditions, ion densities in the subsolar magnetosheath at the inbound crossing were almost a factor of two lower than densities at the outbound crossing.This large change in magnetosheath density is not reflected in the upstream solar wind conditions convected to the magnetopause.
In the magnetosphere, changes in the energy of the ion population below 100 eV in Figure 1 are related to changes in the bulk convection velocity of the plasma, and these velocity changes are in turn related to magnetopause motion (see, Sauvaud et al., 2001;Song et al., 2019).Integrating the convection velocity component normal to the magnetopause (i.e., Vx in Figure 1) over time yields changes in the distance of the magnetopause from the spacecraft, which are shown in the bottom panel in Figure 1 over the 8 min that the spacecraft was in the magnetosphere.The magnetopause remained relatively close to the spacecraft for the first few minutes, but then moved sunward from 09:35 to 09:38 UT as Vx increased to ∼100 km/s.The magnetopause was nearly 1 R E sunward of the spacecraft at 09:38 UT, before moving earthward rapidly and passing over the spacecraft at 09:40 UT.Thus, the magnetopause position changed by about 1 R E without any significant change in the solar wind dynamic pressure.This change is likely related to the rotation of the IMF because the IMF, in particular the cone angle, influences the shape of the magnetopause (see, e.g., Grygorov et al., 2017;Sibeck et al., 2000).

Comparison With Gas-Dynamic and MHD Models
In the gas-dynamic model of the magnetosheath (Spreiter et al., 1966), the plasma density in the subsolar magnetosheath is 4.43 times the upstream solar wind density (with no uncertainty assumed for this model parameter).
For event #9 in Figure 1, the upstream solar wind density was 4.91 ± 0.05 cm −3 and 5.16 ± 0.07 cm −3 for the inbound and outbound magnetopause crossings, respectively.Therefore, the gas-dynamic model predicts subsolar magnetosheath densities of 21.8 cm −3 and 22.9 cm −3 for the inbound and outbound crossings, respectively.
Averaging the observed magnetosheath density in Figure 1 over 1 min from 09:31 to 09:32 UT for the inbound crossing and from 09:40 to 09:41 UT for the outbound crossing, the observed, average subsolar magnetopause densities were 19.2 ± 2.6 cm −3 and 36 ± 3.8 cm −3 , respectively (the uncertainties are the standard deviation of the mean of the 4.5 s density measurements from FPI).Thus, the gas-dynamic model predicts densities that are approximately the same as the density observed by MMS at the inbound crossing and about 60% less than that observed at the outbound crossing.
For event #9, the subsolar magnetopause density was also predicted from the Open Geospace General Circulation Model (OpenGGCM) (e.g., Raeder et al., 1996), using the following (static) parameters as input: Vsw = −400 km/s, Nsw = 5 cm −3 , IMF(GSM) = (0, 0.5, −4) nT, dipole tilt = −25°.The density profile along the Earth-Sun line is shown in Figure S1 in Supporting Information S1.The subsolar magnetopause location is defined here as the maximum in the second derivative of this density profile.For event #9, the magnetopause location was 10.8 R E and the density was 15.8 cm −3 .Thus, the subsolar magnetopause standoff distance in the MHD simulation and in the observations are in agreement.However, the MHD density at the subsolar magnetopause is lower than the observed density at the inbound crossing and substantially lower than the observed density at the outbound crossing.
Subsolar densities from the gas-dynamic model and densities and locations from MHD simulations were also derived for the other 18 events in Table S1 in Supporting Information S1. Figure 2 compares observed and predicted subsolar densities for the 19 events using the ratio of observed to predicted values.The top panel shows observed/gasdynamic density ratios and the bottom panel shows observed/MHD density ratios.There is considerable scatter around a ratio of 1 with on apparent trend with northward versus southward IMF or with inbound versus outbound crossings.Averaging over all crossings in the 19 events, the gasdynamic model tends to over-predict the subsolar magnetosheath density, with a weighted mean of 0.77 ± 0.01 and no difference between averages of inbound only and outbound only ratios.The MHD model tends to under-predict the subsolar density, with a weighted mean of 1.09 ± 0.02 and weighted means of 1.09 ± 0.03 and 1.14 ± 0.03 for outbound and inbound ratios, respectively.Because of the substantial scatter of the ratios, an important metric is the median values.The median for the observed/gasdynamic ratio = 0.86 and the median for the observed/MHD ratio = 1.22. 10.1029/2023GL105553 5 of 9 Thus, the gasdynamic model over-estimates the density by 14% and the MHD model tends to under-estimate the subsolar magnetopause density by 22%.
Figure 3 shows the difference, in R E , of the observed and MHD predicted subsolar magnetopause locations.Almost all values are less than zero with a median of −0.5 R E .Thus, most of the time, MHD over-estimate the magnetopause distance from the Earth by about 0.5 R E for ∼75% of the events.

Discussion
Nineteen events, including a total of 27 subsolar magnetopause crossings, were surveyed to determine trends in the observed and modeled densities and subsolar magnetopause locations.Figures 2 and 3 and Table S1 in Supporting Information S1 show that the modeled densities and magnetopause locations are not grossly dissimilar from the observations.However, densities derived from global MHD have considerable scatter and are up to about a factor of 2 different from observed densities.MHD magnetopause locations show similar, significant scatter and deviate from observed locations by up to ∼1 R E .Considering the medians in Figures 2 and 3, gasdynamic models tend to over-estimate the subsolar magnetopause density by ∼14% while MHD models tend to under-estimate the density by ∼22%.MHD models tend to predict subsolar magnetopause locations that are ∼0.5 R E further from the Earth (i.e., about 5% of the typical standoff distance of 10 R E ) than the observed locations.About half of this overestimate may be explained by the fact that MHD models are single fluid, while the actual solar wind has higher dynamic pressure because it contains about 4% He 2+ .
The systematic over-estimate of the density using the gasdynamic model results in an under-estimate of the geocoronal density.In Equation 1, an over-estimate of the density of 14% translates into an over-estimate of the ion flux and an under-estimate of the geocoronal density of similar magnitudes.
The systematic under-estimate of the density and over-estimate of the magnetopause location in the MHD model have important consequences for calculating the geocoronal density using Equation 1.Both systematic effects result in an over-estimate of the geocoronal density.The lower density translates directly into a lower ion flux, resulting in an over-estimate of about 22% using an MHD model with the upstream solar wind as the boundary condition.The over-estimate in the calculated geocoronal density that results from the larger magnetosphere in the MHD model is more difficult to quantify.At energies ∼1 keV, the ion flux in Equation 1 primarily comes from the magnetosheath and the contribution from the magnetosphere is much lower (see, e.g., the differences in the magnetosheath and magnetosphere ion fluxes at 1 keV in Figure 1).Therefore, increasing the magnetopause standoff distance by 0.5 R E results in a lower ion flux; however, the magnitude of the reduction depends on the field-of-view and resolution of the ENA camera and the pixel resolution of the MHD simulation.For the simulation in Figure S1 in Supporting Information S1 and the data derived from the IBEX ENA imager with a FOV that covers essentially the entire magnetosheath (see, e.g., Fuselier et al., 2020;Starkey et al., 2022), the reduction in the ion flux caused by the sunward shift of the entire subsolar magnetosheath density profile by 0.5 R E is roughly ∼10%-20%.Thus, combined with the under-estimate of the density, the MHD model over-estimates the geocoronal density by ∼30%-40%.
The scatter in the data in Figures 2 and 3 plays an important role in the estimate of the uncertainty in the geocoronal density using Equation 1.This scatter is about a factor of two larger than the systematic over-and under-estimates of the parameters illustrated in these figures.Therefore, using the gasdynamic and MHD models with the upstream solar wind as a boundary condition results in an uncertainty of about a factor of 2 in the determination of the geocoronal density using Equation 1.
The systematic differences and uncertainties are reduced if the MHD or gasdynamic results are normalized by in situ observations at the subsolar magnetopause.However, ENA observations typically require integration over many minutes and the question then becomes which in situ observations to use.For event #9 in Figure 1, the magnetosheath densities at the inbound and outbound magnetopause crossings differ by more than 50% and for event #19 (including the outbound crossing used in Fuselier et al. (2020)) the densities differ by a factor of 2.
Understanding which observations to use requires reconsidering the terms in Equation 1, in particular J Ion (E,l).With J Ion (E,l) in mind, Figure 4 provides at least a partial answer to the question of which observations to use.Plotted are the phase space densities in the plasma rest frame as a function of velocity for the inbound and outbound crossings for events #6, #9, and #19.Identified in the panels are the center thermal velocities for energy steps 2 through 6 for the IBEX-Hi ENA camera (Funsten et al., 2009).These log-spaced energies span center energies from 0.7 to 4.3 keV.For events #9 and #19, the largest changes in the plasma distributions from inbound to outbound occur below 400 km/s, in the core of the magnetosheath distribution.At higher velocities, including those corresponding to the IBEX energy steps, there is little difference between the plasma phase space densities at the inbound and outbound crossings.Thus, based on these two events, the answer to which observations to use is that, for the IBEX-Hi camera energies, it doesn't matter.Event 6 shows that this is not always the best answer.For this event, it is probably better to use an average of the inbound and outbound phase space densities.
The higher energies are not as affected because this part of the magnetosheath distribution originates at the bow shock.At the shock, a portion of the solar wind distribution reflects off the shock front, gains energy as it propagates back into the upstream, and returns and crosses the shock and enters the downstream (see Gosling & Robson, 1985).The local effects at the subsolar magnetopause are in the core of the distribution and not in the shock-generated higher-velocity parts of the distribution (see also, Starkey et al., 2022).For the IBEX-Hi energies, the break between the core and the reflected/transmitted population occurs between energy steps 2 and 3.  Phase space density versus velocity in the plasma rest frame for three pairs of inbound and outbound magnetopause crossings.For events #9 and #19, the total density is much higher for the outbound crossings (red circles) than for the inbound crossings (blue squares).Despite this density difference, the phase space densities for they two types of crossings are similar at energies above ∼0.7 keV (IBEX-Hi energy step 2).Event #6 shows the opposite at energies >1 keV (IBEX-Hi energy step 3).
In conclusion, in determining the neutral geocoronal density at ∼10 R E from the Earth, a global model of the bow shock, magnetosheath, and magnetosphere is necessary to link line-of-sight ENA measurements with the local in situ ion measurements.The two global models investigated here provide this link up to a point.There are some systematic differences between observations and the models; however, the overwhelmingly important takeaway from the comparison between models and observations in Figure 2 through 4 is that the use of a predictive model introduces an inherent uncertainty of at least a factor of 2 in the final determination of the geocoronal density."The Earth's exosphere and its response to space weather."Support for this study comes from the IBEX mission as a part of NASAʼs Explorer program under Grant 80NSSC20K0719 and from the NASA Goddard Space Flight Center MMS prime contract NNG04EB99C.ChatGPT was used to draft the plain language summary from the abstract.

Figure 1 .
Figure 1.(top to bottom) Ion, electron energy time spectrograms, ion and electron densities, ion velocity components (in GSE), solar wind dynamic pressures convected to the magnetopause, and distances from MMS to the magnetopause when the spacecraft was in the magnetosphere.The ion density at the inbound magnetopause crossing (09:32 UT, left blue arrow) is almost a factor of 2 lower than that at the outbound crossing (09:39:30 UT, right blue arrow).Both magnetopause crossings probably occurred because of the relatively small changes in the solar wind dynamic pressure.The bottom panel shows that, while the spacecraft was in the magnetosphere, the magnetopause moved almost 1 R E sunward by 0938 UT, before moving rapidly earthward and crossing over the spacecraft at 09:39:30 UT.These large changes in the magnetopause location are not reflected in the changes in the solar wind dynamic pressure.

Figure 2 .
Figure 2. Observed/predicted subsolar density for the gasdynamic (top panel) and magnetohydrodynamic (MHD) (bottom panel) models.The gasdynamic model tends to predict higher densities at the subsolar magnetopause than observed, while the MHD model tends to predict lower densities than observed.

Figure 3 .
Figure3.Observed-predicted subsolar magnetopause standoff distance from the Earth using observations and an magnetohydrodynamic (MHD) model.For most events and magnetopause crossings, the MHD model tends to predict a magnetopause that is farther from the Earth than observed.

Figure 4 .
Figure 4. Phase space density versus velocity in the plasma rest frame for three pairs of inbound and outbound magnetopause crossings.For events #9 and #19, the total density is much higher for the outbound crossings (red circles) than for the inbound crossings (blue squares).Despite this density difference, the phase space densities for they two types of crossings are similar at energies above ∼0.7 keV (IBEX-Hi energy step 2).Event #6 shows the opposite at energies >1 keV (IBEX-Hi energy step 3).