Theory Helps Observations: Determination of the Shock Mach Number and Scales From Magnetic Measurements

The Mach number is one of the key parameters of collisionless shocks. Understanding shock physics requires knowledge of the spatial scales in the shock transition layer. The standard methods of determining the Mach number and the spatial scales require simultaneous measurements of the magnetic field and the particle density, velocity, and temperature. While magnetic field measurements are usually of high quality and resolution, particle measurements are often either unavailable or not properly adjusted to the plasma conditions. We show that theoretical arguments can be used to overcome the limitations of observations and determine the Mach number and spatial scales of the low-Mach number shock when only magnetic field data are available.


INTRODUCTION
Collisionless shocks [1] are one of the most ubiquitous phenomena in space plasmas. The unfading interest in collisionless shocks is related to the fact that they are the most efficient accelerators of charged particles in the known universe [2][3][4][5][6][7][8][9][10][11][12][13][14]. A shock is a discontinuity in terms of magnetohydrodynamics (MHD) [15], at which the flow velocity along the shock normal drops while the density increases. In reality, this discontinuity has a finite width, and the electric and magnetic fields vary continuously inside the transition from the upstream region of low density and entropy to the downstream region of higher density and entropy. A collisionless shock efficiently converts the energy of the directed flow into the thermal energy of charged particles, the energy of accelerated particles, and the field energy. The conversion occurs via the interaction between charged particles and the electric and magnetic fields of the shock. Thus, understanding the processes inside the shock requires, first and foremost, knowledge of the fields inside the transition layer together with their dependence on time and space. Observational determination of this is not an easy problem. Direct observational separation of spatial and temporal variations is not possible with singlespacecraft measurements. Multi-spacecraft measurements help (e.g., Russell et al. [16], Dunlop [17]) but to a limited extent. In addition, determination of the Mach number and the spatial scales requires sufficiently good particle measurements. At present, most spacecraft that study shocks are not designed to properly measure parameters of narrow cold beams [18][19][20][21][22]. Therefore, the solar wind is not resolved to the required precision. At the same time, the available magnetic field measurements are typically very good. It would be helpful if the Mach number and the spatial scales could be reasonably estimated using magnetic field data alone. It makes sense to first attempt to develop such a methodology for low-Mach number shocks, which are expected to be nearly stationary and planar and have moderately structured profiles [23][24][25][26][27][28][29][30][31]. There is no strict definition of what low-Mach numbers are. Typically, a shock with the Alfvénic Mach number M ≲ 4 may be expected to suit the above purposes. In contrast, in several theoretical works, dependences of several spatial scales on the Alfvénic Mach number have been derived, and certain tools have been proposed for application to the magnetic field data [32][33][34][35][36][37][38].
The Mercury bow shock is typically a low-Mach number shock and the 20 Hz MESSENGER magnetic field measurements are sufficiently good for the application of the proposed methods [39,40]. This study demonstrates the efficiency and consistency of theoretical predictions by applying the proposed methods to selected shocks. The study is organized as follows: first, we summarize the theoretical estimates of the width of various shock features proposed so far (Section 2). Next, we briefly outline the numerical analysis used as an additional check (Section 3). After that, we analyze in detail two selected shock profiles.

THEORETICAL BASIS FOR ESTIMATES OF THE SHOCK SCALES AND MACH NUMBER
Several theoretical estimates are available for planar stationary shocks when the ion reflection is weak. Note that these estimates allow determining the scales in terms of the ion inertial length but do not allow determining the latter in dimensional units (km). The analysis in the study is done in the normal incidence frame (NIF), that is, the shock frame in which the upstream plasma flow is along the shock normal.

Phase-Standing Whistler Precursor
In the shock with a coherent whistler precursor, the wavelength of the whistler wave train can be estimated if the precursor is assumed to be phase standing in the shock frame [29,36,41,42]: where V u is the NIF upstream plasma velocity, θ Bn is the angle between the shock normal and the upstream magnetic field, M = V u /V A is the Alfvénic Mach number, V A B u / 4πn u m p is the Alfvén speed, ω pi 4πn u e 2 /m p is the proton plasma frequency, n u is the upstream proton number density, and m p is the proton mass. For simplicity, the plasma is assumed to consist of protons and electrons. Possible small addition of α-particles is neglected. Then, conversion of the temporal scale into spatial scale can be done by measuring the time between two successive maxima of the whistler precursor. Equation 2 defines the characteristic whistler length l W together with its relation to the whistler wavelength λ.

Foot Length
Usually, the expression by Gosling and Thomsen [43] is used to estimate the foot length. This expression is derived for a specularly reflected ion entering the shock with the velocity of the flow. More detailed studies have shown that ion reflection is non-specular and the foot length is substantially smaller and can be estimated as where M is the Alfvénic Mach number and Ω u = eB u /m p c is the upstream gyrofrequency [33,38].

Downstream Magnetic Oscillations
Coherent magnetic oscillations arise due to the gyration of the ion distributions produced at the shock crossing [28,[44][45][46]. If the oscillations are periodic, the effect of the reflected ions is weak. For directly transmitted ions, the distance between two successive maxima can be estimated as [47] where V drift is the component of the downstream flow velocity along the shock normal. In the quasi-perpendicular case, this drift velocity can be approximated as V drift ≈ V u (B u /B d ), which is sufficient for the estimate. The amplitude of the oscillations decreases due to the gyrophase mixing. The decay is faster for higher upstream ion temperatures [45]. If the effect of reflected ions is significant, other peaks/dips may exist between the main maxima/minima [48,49]. In this case, Δ would correspond to the distance between a maximum and the next second maximum, or approximately the twice distance between two adjacent maxima.

Distance From the Overshoot Maximum to the Undershoot Minimum
This distance is more difficult to evaluate. As a rough approximation, it can be estimated as a gyroradius of the ion, which just crossed the ramp, in the downstream magnetic field. Within the narrow shock approximation [47,50] the ion speed upon crossing the shock is where ϕ NIF is the cross-shock potential in the normal incidence frame. For the low upstream temperature, the ion velocity at the entry to the ramp is approximately equal to the flow speed. Therefore, a rough approximation for the distance from the overshoot maximum to the undershoot minimum would be This estimate is less reliable than the others because the gyration occurs in the inhomogeneous magnetic field between the ramp and the undershoot. Note the use of the particle velocity and not the drift velocity because the analysis refers only to a part of the gyration.

Noncoplanar Magnetic Field
In laminar shocks or shocks with weak ion reflection, the noncoplanar magnetic field component inside the ramp is approximately [34,35] B y l W dB z dx This approximation is not valid behind the ramp where the ion distributions begin to gyrate as a whole.

BASICS OF NUMERICAL ANALYSIS
Further check of the estimates will be done using the adjustable test particle analysis [51,52]. We briefly describe the principles of the method. Model magnetic and electric fields are chosen, determined by a small number of parameters. The equations of motion for ions are solved numerically in these fields. A large number of ions from the initial upstream ion distribution are thus numerically traced across the shock front, which provides the ion distribution everywhere. The relevant moments of this distribution are determined numerically and used in the corresponding conservation law, pressure balance (see below). The pressure balance equation is used to derive the magnetic field, which would be consistent with the numerically found ion pressure. The derived magnetic field is compared with the chosen model field. The parameters of the chosen model are varied until a reasonable agreement is achieved. For a model low-Mach number shock, a reasonable model profile is given by the following expressions [45,51]: The parameters R z = B d,z /B u,z and θ Bn are taken from observations. The shock ramp width D is one of the adjustable parameters. The electrostatic field along the shock normal is modeled using E x ∝ dB z dx and whereas the noncoplanar component of the magnetic field is taken from Eq. 7. An incident Maxwellian distribution of ions is traced across the shock. The total ion pressure is numerically determined as a function of the coordinate x along the shock normal. The pressure is inserted in the pressure balance: where p i,xx is the total ion pressure, which is the sum of the dynamic and kinetic pressure. The electron pressure p e is assumed adiabatic, p e ∝ n γ , γ = 5/3. We also assume β e = β i , where β 8πnT u /B 2 u , T u being the upstream temperature. Equation 11 is used to numerically find the magnetic field that would be consistent with the pressure balance. The shock parameters M, s, D, β are varied until the magnetic field magnitude derived from Eq. 11 converges to the initial model magnetic field well beyond the ramp, and the derived and initial ramp slope are in agreement. For low-Mach number shocks,  magnetic oscillations of the derived field around the downstream value are small so that their effect on the ion motion can be neglected. Earlier studies [44,45,51] have shown that, for low-Mach number shocks, the method is not too sensitive to the variation of the parameters within reasonable limits expected by the theory so that the adjustment procedure converges rapidly. Figure 1 shows the shock crossing at 2011/083/12:25:00 at two different scales. The top panel shows the very transition as the magnetic field within ±20 s around the crossing. The bottom panel shows the shock and the magnetic field around the shock, ± 180 s, around the crossing. Figure 2 shows the magnetic field components in MSO coordinates and the magnetic field magnitude. The red vertical lines mark the upstream region used to calculate the normal from the magnetic coplanarity, and the blue vertical lines mark the downstream region. The black vertical line marks the crossing time. The upstream magnetic field vector B u is found by averaging the magnetic field vector over the chosen upstream region. The downstream magnetic field vector B d is determined similarly. In the following, B u = |B u | and B d = |B d |. Figure 3 shows the normalized magnetic field rotated into the shock coordinates: x is along the shock normal and y is the noncoplanarity direction. The fields are normalized on the upstream magnetic field magnitude B u . The rotation is made using the unit vectors:

A LOW-MACH NUMBER SUBCRITICAL SHOCK
The main magnetic field, B z , has a clear monotonic ramp, a whistler precursor, and a barely noticeable overshoot: B d /B u = 1.7 and max |B|/B u = 1.75. There is no visible foot. The region between the blue and the red vertical lines in the Figure 3 is the whistler precursor. The region between the red and the green vertical lines is the ramp.
The black vertical line passes through the maximum of the overshoot. Judging by the above-described features of the profile, the shock should be a low-Mach number subcritical shock with no or negligible ion reflection. The noncoplanar magnetic field exhibits fluctuations at the spatial scale of the ramp and the whistler precursor. The normal component fluctuates inside the ramp. The spatial scales of these fluctuations are substantially smaller than the ramp width. We tend to interpret these deviations from planarity as a small-scale rippling inside the ramp, which propagates along the shock surface. The normal component of the magnetic field does not have fluctuations at the whistler spatial scale outside the ramp, which makes us conclude that the whistler propagates or phase stands along the normal. The angle between the shock normal and the upstream magnetic field is θ Bn = 67°and cos θ Bn = 0.39. Moderate changes of the upstream and downstream intervals did not affect the normal determination noticeably. Figure 4 is a close-up of Figure 3 but plotted using points to show explicitly the resolution of the magnetic field measurements. The horizontal black line marks B = 0. The magnetic field increase in the ramp is nearly linear, and the noncoplanar magnetic field has a rather broad maximum. The behavior is consistent with the relation (7).
The magnetic field is measured as a function of time in the spacecraft frame. Accordingly, in the following, all scales are given as temporal equivalents of the corresponding spatial scales. Proper conversion of temporal durations into spatial lengths should be done by multiplying by the unknown shock speed V sh in the spacecraft frame. This speed cannot be determined in dimensional units (km/s) without density measurements. Inspection of Figure 3 has shown that the time separation of the two successive maxima of the whistler precursor is Δt W = 3.19 s, which gives the spatial-to-temporal correspondence λ = 2πL W = V sh Δt W (see (2)), or l W = V sh Δt W /2π = V sh · 0.5 s. Using the estimate based on the noncoplanar magnetic field from Eq. 7 and d/dt = V sh (d/dx), one gets l W = V sh · 0.6 s. The agreement is quite good. The ramp duration is Δt r = 2.4 s so that the ramp width is D = V sh Δt r = V sh · 2.4 s and FIGURE 3 | The normalized magnetic field rotated into the shock coordinates: x is along the shock normal, and y is the noncoplanarity direction. Colors are as in Figure 2. The region between the blue and the red vertical lines is the whistler precursor. The region between the red and the green vertical lines is the ramp. The black vertical line passes through the maximum of the overshoot.
Frontiers in Physics | www.frontiersin.org March 2022 | Volume 10 | Article 852720 Assuming D ≈ c/ω pi [27], we get an estimate of the Alfvénic Mach number: M A ≈ 1.6. The shock is laminar, indicating low β [23,25,26]. Figure 5 shows the results of the adjustable test particle analysis. The shock angle was taken from the above determination using magnetic coplanarity, and the ramp width was taken as one ion inertial length. In contrast, the Mach number M, the cross-shock potential s, and the ion β i were varied. The best convergence to the downstream magnetic field was found for M = 1.65, β i = 0.1, and s = 0.63. The overshoot in the derived profile is also in agreement with the observations. Figure 6 shows the shock crossing at 2013/047/06:57:00 in a format similar to Figure 1. Figure 7 shows the magnetic field components in MSO coordinates, as well as the magnetic field magnitude, the shock crossing moment, and the upstream and downstream regions used in the analysis. The format is similar to that of Figure 2. Figure 8 shows the normalized magnetic field rotated into the shock coordinates, similar to Figure 3. This shock has a clear overshoot and a noticeable foot, which means that this shock is expected to be a moderately supercritical shock. The downstream-to-upstream main magnetic field ratio is R d = B d /B u = 2.66 and max |B|/B u = 3.54. The angle between the shock normal and the upstream magnetic field is θ Bn ≈ 63°and cos θ Bn ≈ 0.45. Moderate changes of the upstream and downstream intervals did not affect the normal determination noticeably. High-frequency turbulence is present. We remove it using the   wavelet denoising as follows: 1) 2 11 points of data are taken around the crossing time to cover sufficiently the upstream and downstream regions, 2) the Daubechies-10 wavelet transform is applied, 3) five smallest scales are removed, and 4) the inverse wavelet transform is performed. The procedure is done for B x , B y , B z and separately for |B| because, otherwise, the foot region is smeared out. Figure 9 shows the magnetic field magnitude with the high-frequency noise removed. The maximum magnetic field is now R m = |B dn, max |/B u = 3.3. We shall adopt this value as the maximum overshoot magnetic field. Figure 10 shows some meaningful points marked at the magnetic field profile. Among these, the most important for us will be the beginning of the foot at t 1 = − 10 s, the end of the foot and the beginning of the ramp at t 2 = − 1.9 s, with the elevation of the magnetic field of ΔB foot /B u = 0.1, and the overshoot maximum at t 3 = 0.7 s, with the additional elevation of the magnetic field of ΔB ro /B u = 2.2. The minimum of the undershoot occurs at t 5 = 5.8 s, and the magnetic field there drops to B 5 /B u = 2, well below the downstream value. The width of the foot is estimated as L foot ≈ 0.5M(c/ω pi ) = V sh · 8.1 s, where M is the Alfvénic Mach number [33]. Figure 11 shows all three components of the magnetic field (1/V sh ) dB z dt . This behavior implies that the relation (7) may be a good estimate [35]. Using it, we get

A LOW-MACH NUMBER SUPERCRITICAL SHOCK
Note that this relation is not valid behind the ramp because of the strong non-gyrotropy of the ion distribution [32]. Together with the estimate of L foot , we get the estimate of the Mach number Another estimate of the Mach number can be obtained from [50] where s 2eϕ NIF /m p V 2 u is the normalized NIF cross-shock potential. Taking s = 0.5, one gets M ≈ 4.1. The two estimates  agree very well. There should be no illusions though because all these are approximations. Using M ≈ 4, we get The distance between the maximum of the overshoot and the minimum of the undershoot is 5 s. A rough estimate of this distance is which is not bad at all. Consistence of all these estimates encourages to conclude that the chosen s = 0.5 is not far from reality. The absence of downstream magnetic oscillations, that is, only one overshoot and undershoot, indicates high β [45]. Figure 12 shows the results of the adjustable test particle analysis. The best convergence to the downstream magnetic field was found for M = 4, β i = 0.75, and s = 0.5. The overshoot in the derived profile is also in agreement with the observations. The position of the undershoot is close to the predicted. However, the magnetic field in the undershoot of the derived profile is somewhat higher than the observed one probably because the overshoot modifies the ion motion. The initial profile used for adjustment was a simple tanhlike profile. A better agreement could potentially be achieved by using a more sophisticated model with an overshoot and undershoot. We leave this issue for further studies.

VERIFICATION WITH MAGNETOSPHERIC MULTISCALE
There were no particle data for the analyzed MESSENGER shocks. The above analysis has been done using magnetic field   measurements alone. It is desirable to verify this analysis with shocks for which the Mach number and the scales could be determined by the above methods and, independently, by the conventional methods involving additional measurements. For this task, an MMS1 shock was selected with sufficient magnetic features to apply the above approach. Two independent analyses have been performed: one applied the theoretical estimates to the magnetic profile without utilizing any other information and the other was done in the standard way. We present the comparison below. The shock crossing occurred at 2020-11-12 14:36:04. Figure 13 shows a part of the shock, |B|/B u , in GSE coordinates, with the time set to zero at the shock crossing. The black line shows the normalized magnetic field magnitude. The sampling rate is 16 measurements per second. The upstream magnetic field is determined by averaging over about 20 first seconds of the Figure 13. The downstream magnetic field is determined by averaging over about 50 last seconds of the figure.
The shock normal is determined using magnetic coplanarity. The found shock normal isn (0.627 5, 0.778 6, −0.000 5). The red line shows the denoised magnetic field. The denoising (i.e., removing high-frequency fluctuations to only retain what is assumed to be the stationary profile) is done by applying discrete wavelet transform. The Daubechies 10 wavelet transform was applied to the 4,096 data points, and five finest scales were removed. The vertical blue lines show the chosen beginning and end of the foot used for further analysis. The duration of the foot is 10.16 s. The denoised magnetic field is used to determine R d ≈ 2.5 and R m ≈ 3.6. Figure 14 shows the noncoplanar and main magnetic field components in the vicinity of the ramp. Using Eqs 3, 7, the Mach number is estimated as M ≈ 3.35, and the ion inertial length corresponds to the duration ≈ 6.0 s. Using the distance between the two adjacent downstream maxima gives a result that is inconsistent with other estimates. Using twice this distance [48,49], together with the foot length, also gives M ≈ 3.35. Note that the precision of the determination of the durations used for the estimates and the application of the wavelet transform is not sufficient to ensure the above excellent agreement of the Mach number estimates. The shock was independently analyzed using the MMS magnetic field and particle (density and velocity) measurements. The upstream and downstream magnetic field vectors in GSE coordinates are B u = (1.48, − 6.53, 4.52) nT and B d = (7.93, − 12.94, 12.4) nT, respectively. The upstream and downstream velocity Note that the model shock normal and the normal found earlier from coplanarity are quite similar. In order to avoid confusion, velocities V 1 and V 2 are the bulk plasma velocities relative to the spacecraft in the upstream and downstream regions, respectively. Value V u is the NIF upstream plasma speed. The Mach number is The derived Mach number is M = 3.68 with the shock speed of V sh = 20.3 km/s. Because the spacecraft instruments are not quite appropriate for catching the cold solar wind, the OMNI density [54] is often used to replace the upstream ion density measured by the spacecraft. In this case, the OMNI density is N u,OMNI = 7.4 cm −3 , slightly lower than the MMS measured value, and the corresponding Mach number is M = 3.36, with the shock speed of V sh = 10 km/s. The Mach numbers obtained in the two approaches differ by less than 10%, which is within the precision of the determination of the shock parameters.

DISCUSSION AND CONCLUSION
In this study, we applied existing theoretical estimates to the magnetic field measurements to determine the Aflvénic Mach number and the scale parameters of two low-Mach number shocks. One of these shocks has a very low overshoot and a clear whistler precursor. Another one possesses a substantial overshoot and a foot. In both cases, we estimated the Mach number using at least two independent theoretical approaches and found good agreement between the various methods. In addition, this allowed us to determine the correspondence of the duration of the measurement of a particular feature to its physical spatial scale in terms of the upstream convective gyroradius and/or ion inertial length. As always, the determination of the shock parameters requires making some assumptions, such as stationarity and planarity. Although the methods were applied to rather clean shocks with classical profiles at this stage, they will possibly allow extension to less favorable cases, in part by comparison with the success of the present study. The method has been tested with an MMS observed shock, for which sufficiently good particle measurements are also available. The Mach number obtained with the magnetic measurements and theory and the Mach number obtained with both magnetic field and particle measurements differ by less than 10%, which is encouraging.

DATA AVAILABILITY STATEMENT
Publicly available datasets were analyzed in this study. These data can be found at: https://pds-ppi.igpp.ucla.edu/mission/MESSENGER/ MESS/MAG https://lasp.colorado.edu/mms/sdc/public/ AUTHOR CONTRIBUTIONS MG, EG, CR, and AD participated in the analysis and paper writing.