The Structure of Shocks in the Very Local Interstellar Medium

The Voyager 1 spacecraft crossed into the very local interstellar medium (VLISM) after passing the heliopause (HP) on 2012 August 25. Anomalous cosmic rays and termination shock particles, which are dominant in the inner heliosheath (IHS), vanished abruptly and the galactic cosmic-ray flux increased significantly. The Voyager 1 magnetometer detected a change in the magnetic field strength (Krimigis et al. 2013; Stone et al. 2013; and Webber & McDonald 2013). Voyager 1 is now making in situ measurements of the interstellar medium. Although the plasma instrument on Voyager 1 stopped working in 1980, the electron number density can be estimated from the frequency of electron plasma oscillations (Gurnett et al. 2013). A remarkable discovery made by Voyager 1 was the observation of shock waves in interstellar space (Burlaga et al. 2013). The shocks seen in the VLISM are due to interplanetary shocks that propagate through the supersonic solar wind, collide with the heliospheric termination shock (HTS), and then propagate through the inner heliosheath until the HP, where the shocks are partially reflected and transmitted into the VLISM (Whang & Burlaga 1994; Zank & Müller 2003; Washimi et al. 2011; Zank 2015).

The interstellar shocks are surprising in the regard that they possess properties quite different than those of heliosphere shocks (Burlaga et al. 2013). The 2012 VLISM shock was observed (Burlaga et al. 2013) to have a weak jump in the magnetic field, ${\boldsymbol{B}}$ (compression ratio of about B₂/B₁ ≈ 1.4 where 2/1 refers to the downstream/upstream state of the shock), and appeared to correspond to a weak, low beta, and subcritical shock (Mellott & Greenstadt 1984). Burlaga et al. (2013) showed that the shock is quasi-perpendicular (the angle between shock normal and ${\boldsymbol{B}}$ is about 85°). The significant property of the shock is its very smooth nature with very large width, being about 10⁴ times broader than a shock with similar properties at 1 au. Thus far there is no theoretical explanation for such a broad shock structure in the VLISM. Here we present a theory describing the VLISM that accounts for the structure of VLISM shocks, and we model the VLISM shock observed at the end of 2012 by Voyager 1. We show that the electron and proton collisional timescale in the VLISM is much smaller than the characteristic dynamical timescale in the VLISM with corresponding short collisional mean free paths. Consequently, the VLISM is collisional with respect to the thermal plasma. Thermal collisions introduce heat conduction and viscosity (dissipation) into the system, and these processes are responsible for determining the structure of VLISM shocks. Fast and hot neutral atoms created in the supersonic solar wind and IHS are deposited in the VLISM (Zank et al. 1996; Zank et al. 2014, 2016), which, through secondary charge exchange, leads to the creation of suprathermal PUIs in the VLISM. The VLISM PUIs contribute a very small pressure in the VLISM due to their small number density compared to that of the thermal gas. Although the primary dissipation mechanism at the HTS is reflected PUIs (Zank et al. 1996; Richardson et al. 2008), they are insufficient to provide the dissipation for VLISM shocks.

The paper is organized as follows. We first derive the basic equations describing the VLISM environment and the structure of shock waves therein. In Section 3, we model the VLISM shock and compare it with the observations (Burlaga et al. 2013). Finally, we discuss and summarize our results.

The VLISM is comprised of thermal electrons and protons, suprathermal PUIs, and interstellar neutral atoms. In the VLISM, for the scales of interest here, the charge exchange mean free path of protons and hydrogen is sufficiently large that we can neglect their coupling. Zank et al. (2014) discussed collisions between PUIs and thermal VLISM electrons and protons, showing that the collisional mean free path is very large with the implication that PUIs are not equilibrated with the thermal VLISM on scales of less than 75 au. PUIs are scattered by both in situ and self-excited magnetic field fluctuations and therefore experience pitch-angle scattering. A collisionless Chapman–Enskog expansion (Zank et al. 2014, 2016) shows that PUIs can be described by an isotropic pressure, P_I, stress tensor, ${{\boldsymbol{\Pi }}}_{I}$ , and heat flux, q_I, given by

$$\begin{eqnarray}&&{{\boldsymbol{P}}}_{I}={P}_{I}{\boldsymbol{I}}+{{\rm{\Pi }}}_{I},\end{eqnarray} \tag{ 1 }$$

$$\begin{eqnarray}{{\boldsymbol{\Pi }}}_{I}=\left(\begin{array}{ccc}1 & 0 & 0\\ 0 & 1 & 0\\ 0 & 0 & -2\end{array}\right)\displaystyle \frac{{\eta }_{k{\ell }}}{2}\left(\displaystyle \frac{\partial {U}_{k}}{\partial {x}_{{\ell }}}+\displaystyle \frac{\partial {U}_{{\ell }}}{\partial {x}_{k}}-\displaystyle \frac{2}{3}{\delta }_{k{\ell }}\displaystyle \frac{\partial {U}_{m}}{\partial {x}_{m}}\right),\end{eqnarray} \tag{ 2 }$$

$$\begin{eqnarray}&&{q}_{I}(x,t)=-\displaystyle \frac{1}{2}{{\boldsymbol{K}}}_{I}.\bigtriangledown{P}_{I},\end{eqnarray} \tag{ 3 }$$

where ${{\boldsymbol{P}}}_{I}$ is the PUI pressure tensor, η is the PUI viscosity, ${{\boldsymbol{K}}}_{I}$ is a spatial PUI diffusion tensor, and ${\boldsymbol{U}}$ is the bulk flow velocity.

Consider now the thermal VLISM electrons and protons, and assume both have a nearly Maxwellian distribution. The collisional streaming or frictional timescale for a charged particle a colliding with a stationary background population of charged particles b is given by (e.g., Zank 2014)

$$\begin{eqnarray}&&{({\tau }_{s}^{{ab}})}^{-1}=\displaystyle \frac{{n}_{b}{q}_{a}^{2}{q}_{b}^{2}\mathrm{ln}{\rm{\Lambda }}}{2\pi {\epsilon }_{0}^{2}{m}_{a}^{2}{V}_{{Ta}}^{3}}\displaystyle \frac{{m}_{b}}{{m}_{a}}\left(1+\displaystyle \frac{{m}_{a}}{{m}_{b}}\right)\displaystyle \frac{{T}_{a}}{{T}_{b}}\displaystyle \frac{G({x}_{b})}{{x}_{a}},\end{eqnarray} \tag{ 4 }$$

where G(x) is the Chandrasekhar function,

$$\begin{eqnarray}&&G(x)\equiv \displaystyle \frac{f(x)-{xf}^{\prime} (x)}{2{x}^{2}},\end{eqnarray} \tag{ 5 }$$

x_a/b ≡ v/V_Ta/b, and f(x) is the error function. The n, q, m, and T for particles a and b denote number density, charge, mass, and temperature, respectively. ${V}_{{Ta}/b}=\sqrt{2{K}_{{\rm{B}}}{T}_{a/b}{m}_{a/b}^{-1}}$ is the thermal speed, K_B is Boltzmann's constant, and $\mathrm{ln}{\rm{\Lambda }}$ is the Coulomb logarithm.

The VLISM is collisionally equilibrated, with electrons and protons having the same temperature of about 7500 K (McComas et al. 2015). The electron and proton thermal speeds are 477 km s⁻¹ and 11 km s⁻¹, respectively. The Voyager 1 magnetometer measured the magnetic field magnitude upstream of the VLISM shock to be about 0.38 nT (Burlaga et al. 2013). Based on electron plasma oscillations observed in the vicinity of the 2012 VLISM shock, the electron density upstream of the shock is about 0.06 cm⁻³ (Gurnett et al. 2013, 2014). For these parameters, we can evaluate the proton–proton (pp), electron–electron (ee), electron–proton (ep), and proton–electron (pe) collisional frictional timescales and the mean free paths (Zank 2014).

These are illustrated in Figure 1, which shows that the ep and pe collisional timescales are constant and independent of the velocity,

$$\begin{eqnarray}{\tau }^{\mathrm{ep}} & = & 1.5\times {10}^{5}\,{\rm{s}},{\lambda }^{\mathrm{ep}}=0.47\,\mathrm{au};\\ {\tau }^{\mathrm{pe}} & = & 2.7\times {10}^{8}\,{\rm{s}},{\lambda }^{\mathrm{pe}}=20.4\,\mathrm{au},\end{eqnarray} \tag{ 6 }$$

where τ is the collisional timescale and λ is the mean free path. The timescales for ee and pp collisions are not constant and depend on the velocity. The collisional heat conduction and viscous coefficients depend on the collisional scattering timescales, τ_s, (see, for example, Section 4.8 of Zank 2014). Therefore, we need to estimate τ_s for the pp and ee collisions and we approximate an appropriate collisional time, $\langle {\tau }_{s}\rangle$, as follows. The coefficient of viscosity, η, can be estimated as (Zank et al. 2014)

$$\begin{eqnarray}&&\eta \equiv \displaystyle \frac{4\pi }{3}\int {c}^{2}{\tau }_{s}{f}_{0}{{mc}}^{2}{dc},\end{eqnarray} \tag{ 7 }$$

where $c=v-U$ is the random velocity. ${\tau }_{s}^{\mathrm{ee}}$ and ${\tau }_{s}^{\mathrm{pp}}$ can be approximated as fifth degree polynomial curves from Figure 1 and f₀ is a Maxwellian distribution. Integrating from 0 to 2.5V_T, we obtain the viscosity coefficient, η, for ee and pp collisions. Using the pressure upstream of the VLISM shock and the relation $\eta \simeq P\langle {\tau }_{s}\rangle$, we derive the most probable estimates

$$\begin{eqnarray}{\tau }^{\mathrm{pp}} & = & 4.2\times {10}^{6}\,{\rm{s}},{\lambda }^{\mathrm{pp}}=0.31\,\mathrm{au};\\ {\tau }^{\mathrm{ee}} & = & 7.5\times {10}^{4}\,{\rm{s}},{\lambda }^{\mathrm{ee}}=0.24\,\mathrm{au}.\end{eqnarray} \tag{ 8 }$$

The VLISM is collisional with respect to the thermal plasma since the electron and proton collisional mean free paths are small compared to the almost featureless VLISM. Therefore, dissipation terms such as thermal heat conduction, K_g, and viscosity, η_g, associated with Coulomb collisions (see Zank 2014) determine the structure of VLISM shocks. In the Appendix, we show that the structure of the observed VLISM shock in the absence of the collisional viscous term is not smooth, meaning that heat conduction alone is insufficient to determine the shock structure and an extra dissipation term (viscosity) is necessary to insure a smooth solution (see also Mostafavi et al. 2017a, 2017b).

Zoom In Zoom Out Reset image size

To describe shock waves in the VLISM, we ideally need to consider a multi-fluid system of equations for thermal electrons, protons, and non-thermal PUIs. Solving such a system of equations is complicated, but we can reduce the model to a single-fluid-like system of equations while retaining the effects of an equilibrated thermal plasma and non-equilibration with the PUI components. Since PUIs are not equilibrated with the background plasma in the VLISM (Zank et al. 2014, 2016), they should be treated as a separate pressure component, P_I, in the system. The reduced single-fluid model comprises the continuity, momentum, and energy equations, together with Maxwell's equations, (see Zank et al. 2014 and Mostafavi et al. 2017a, 2017b)

$$\begin{eqnarray}&&\displaystyle \frac{\partial \rho }{\partial t}+{\rm{\nabla }}\cdot (\rho {\boldsymbol{U}})=0;\end{eqnarray} \tag{ 9 }$$

$$\begin{eqnarray}&&\rho \left(\displaystyle \frac{\partial {\boldsymbol{U}}}{\partial t}+{\boldsymbol{U}}\cdot {\rm{\nabla }}{\boldsymbol{U}}\right)\\ &&\quad =\,-{\rm{\nabla }}({P}_{g}+{P}_{I})-{\rm{\nabla }}\cdot ({{\boldsymbol{\Pi }}}_{I}+{{\boldsymbol{\Pi }}}_{g})+{\boldsymbol{J}}\times {\boldsymbol{B}};\end{eqnarray} \tag{ 10 }$$

$$\begin{eqnarray}\displaystyle \frac{\partial {P}_{g}}{\partial t}+{\boldsymbol{U}}\cdot {\rm{\nabla }}{P}_{g}+{\gamma }_{g}{P}_{g}{\rm{\nabla }}\cdot {\boldsymbol{U}} & = & \displaystyle \frac{1}{3}{\rm{\nabla }}\cdot ({{\boldsymbol{K}}}_{g}\cdot {\rm{\nabla }}{P}_{g})\\ & & -({\gamma }_{g}-1){{\boldsymbol{\Pi }}}_{g}\,:({\rm{\nabla }}{\boldsymbol{U}});\end{eqnarray} \tag{ 11 }$$

$$\begin{eqnarray}\displaystyle \frac{\partial {P}_{I}}{\partial t}+{\boldsymbol{U}}\cdot {\rm{\nabla }}{P}_{I}+{\gamma }_{I}{P}_{I}{\rm{\nabla }}\cdot {\boldsymbol{U}} & = & \displaystyle \frac{1}{3}{\rm{\nabla }}\cdot ({{\boldsymbol{K}}}_{I}\cdot {\rm{\nabla }}{P}_{I})\\ & & -({\gamma }_{I}-1){{\boldsymbol{\Pi }}}_{I}\,:({\rm{\nabla }}{\boldsymbol{U}});\end{eqnarray} \tag{ 12 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial }{\partial t}\left(\displaystyle \frac{1}{2}\rho {{\boldsymbol{U}}}^{2}+\displaystyle \frac{{P}_{I}}{{\gamma }_{I}-1}+\displaystyle \frac{{P}_{g}}{{\gamma }_{g}-1}+\displaystyle \frac{1}{2{\mu }_{0}}{B}^{2}\right)\\ &&+\,{\rm{\nabla }}\cdot \left(\displaystyle \frac{1}{2}\right.\rho {\boldsymbol{U}}{U}^{2}+\displaystyle \frac{{\gamma }_{I}}{{\gamma }_{I}-1}{P}_{I}{\boldsymbol{U}}+\displaystyle \frac{{\gamma }_{g}}{{\gamma }_{g}-1}{P}_{g}{\boldsymbol{U}}\\ &&+\,\displaystyle \frac{1}{{\mu }_{0}}{B}^{2}{\boldsymbol{U}}-\displaystyle \frac{1}{{\mu }_{0}}{\boldsymbol{U}}\cdot {\boldsymbol{BB}}+{{\boldsymbol{\Pi }}}_{I}\cdot {\boldsymbol{U}}-\displaystyle \frac{1}{{\gamma }_{I}-1}\displaystyle \frac{{{\boldsymbol{K}}}_{I}}{3}\cdot {\rm{\nabla }}{P}_{I}\\ &&+\,\left.{{\boldsymbol{\Pi }}}_{g}\cdot {\boldsymbol{U}}-\displaystyle \frac{1}{{\gamma }_{g}-1}\displaystyle \frac{{{\boldsymbol{K}}}_{g}}{3}\cdot {\rm{\nabla }}{P}_{g}\right)=0;\end{eqnarray} \tag{ 13 }$$

$$\begin{eqnarray}&&{\boldsymbol{E}}=-{\boldsymbol{U}}\times {\boldsymbol{B}};\ \displaystyle \frac{\partial {\boldsymbol{B}}}{\partial t}=-{\rm{\nabla }}\times {\boldsymbol{E}};\ {\mu }_{0}{\boldsymbol{J}}={\rm{\nabla }}\times {\boldsymbol{B}};\\ &&{\rm{\nabla }}\cdot {\boldsymbol{B}}=0.\end{eqnarray} \tag{ 14 }$$

${\boldsymbol{B}}$ and ${\boldsymbol{E}}$ denote the magnetic and electric fields, P_g is the thermal gas (protons and electrons) pressure, ρ is the total density, and ${\boldsymbol{U}}$ is the bulk flow velocity. We denote all quantities pertaining to PUIs and thermal gas with the subscripts I and g. γ_g/I are the adiabatic indices of the thermal plasma/PUIs. The collisional thermal plasma stress tensor and PUI stress tensor in a 1D form may be expressed as ${{\rm{\Pi }}}_{g/I}=-\tfrac{1}{3}{\eta }_{g/I}\,\tfrac{{{dU}}_{x}}{{dx}}$, where η_g is the thermal collisional viscosity coefficient and η_I is the collisionless PUI viscosity coefficient associated with magnetic fluctuations (Jokipii & Williams 1992; Zank et al. 2014).

A 1D steady-state model in which all physical quantities depend on the x Cartesian coordinate is used and U = (U_x, 0, 0) and B = (0, 0, B_z). The 1D form of Equations (10)–(12) in normalized form is

$$\begin{eqnarray}{\gamma }_{g}{M}_{s1}^{2}\displaystyle \frac{\partial y}{\partial \bar{x}} & = & -\displaystyle \frac{\partial }{\partial \bar{x}}({P}_{g}^{{\prime} }+{P}_{1}{P}_{I}^{{\prime} })-{y}_{{\rm{B}}}{\gamma }_{g}{M}_{s1}^{2}{B}_{z}^{{\prime} }\displaystyle \frac{\partial {B}_{z}^{{\prime} }}{\partial \bar{x}}\\ & & +\displaystyle \frac{1}{3}{\gamma }_{g}{M}_{s1}^{2}\left({\mathrm{Sch}}_{g}+{\mathrm{Sch}}_{I}\displaystyle \frac{{K}_{I}}{{K}_{g}}\displaystyle \frac{{\rho }_{I}}{{\rho }_{g}}\right)\displaystyle \frac{{\partial }^{2}y}{\partial \bar{x}};\end{eqnarray} \tag{ 15 }$$

$$\begin{eqnarray}{{yP}}_{1}\displaystyle \frac{\partial {P}_{I}^{{\prime} }}{\partial \bar{x}} & + & {\gamma }_{I}{P}_{I}^{{\prime} }{P}_{1}\displaystyle \frac{\partial y}{\partial \bar{x}}=\displaystyle \frac{1}{3}{P}_{1}\displaystyle \frac{{K}_{I}}{{K}_{g}}\displaystyle \frac{{\partial }^{2}{P}_{I}^{{\prime} }}{\partial {\bar{x}}^{2}}\\ & & +\displaystyle \frac{1}{3}({\gamma }_{I}-1){\gamma }_{g}{M}_{s1}^{2}{\mathrm{Sch}}_{I}\displaystyle \frac{{K}_{I}}{{K}_{g}}\displaystyle \frac{{\rho }_{I}}{{\rho }_{g}}{\left(\displaystyle \frac{\partial y}{\partial \bar{x}}\right)}^{2};\end{eqnarray} \tag{ 16 }$$

$$\begin{eqnarray}&&y\displaystyle \frac{\partial {P}_{g}^{{\prime} }}{\partial \bar{x}}+{\gamma }_{g}{P}_{g}^{{\prime} }\displaystyle \frac{\partial y}{\partial \bar{x}}=\displaystyle \frac{1}{3}\displaystyle \frac{{\partial }^{2}{P}_{g}^{{\prime} }}{\partial {\bar{x}}^{2}}+\displaystyle \frac{1}{3}({\gamma }_{g}-1){\gamma }_{g}{M}_{s1}^{2}{\mathrm{Sch}}_{g}{\left(\displaystyle \frac{\partial y}{\partial \bar{x}}\right)}^{2},\end{eqnarray} \tag{ 17 }$$

where ${M}_{s1}^{2}={\rho }_{1}{U}_{x1}^{2}/{\gamma }_{g}{P}_{g1}$ is the thermal gas Mach number of the flow, $\bar{x}\equiv x/L$, L ≡ K_g/U₁, $\overline{x}$ is distance normalized to the diffusion length scale, and ${\mathrm{Sch}}_{g/I}\equiv {\eta }_{g/I}/({\rho }_{1g/I}{K}_{g/I})$ is the Schmidt number of the thermal plasma/PUI. The inverse compression ratio is y = U_x/U_x1 and ${P}_{g}^{{\prime} }={P}_{g}/{P}_{g1}$, ${P}_{I}^{{\prime} }={P}_{I}/{P}_{I1}$, ${B}_{z}^{{\prime} }={B}_{z}/{B}_{z1}$ and ${P}_{1}={P}_{I1}/{P}_{g1}$. The inverse Alfvén Mach number squared of the flow far upstream is ${y}_{{\rm{B}}}\,={B}_{z1}^{2}/({\rho }_{1}{U}_{x1}^{2}{\mu }_{0})={M}_{{Az}0}^{-2}$.

Based on the thermal collisional term scales derived above from the expressions in Zank (2014), the dominant collisional terms in a collisional magnetized plasma can be estimated. The small mass of electrons compared to protons means they do not contribute significantly to either the thermal heat conduction or thermal viscosity. Proton–electron collisions are negligible since the mean free path is very large compared to that of proton–proton collisions. The dominant collisional terms derive from proton–proton collisions. The diffusion coefficient associated with proton–proton collisions, K_pp, upstream of the observed VLISM shock is about 10¹⁵ m² s⁻¹, and the thermal viscosity is 2.5 × 10⁻⁸ kg m⁻¹ s⁻¹, giving a thermal Schmidt number of about 0.25. Energetic PUIs outside the heliosphere are generated by secondary charge exchange, and their number density is very small in the VLISM (n_p ≈ 5 × 10⁻⁵ cm⁻³; Zirnstein et al. 2014), therefore they contribute a very small pressure compared to the thermal plasma (P_I1/P_g1 ≈ 0.05). A characteristic value for the PUI heat conduction, K_I, in the VLISM is ∼10 ¹⁴ m² s⁻¹. The pitch-angle scattering timescale in the VLISM is taken to be $3\times 2\pi \times {{\rm{\Omega }}}_{g}^{-1}$, giving 5.17 × 10² s, which gives a value for the PUI viscosity, η_I, as ∼3.57 × 10⁻¹³ kg m⁻¹ s⁻¹, making the PUI Schmidt number about 10⁻³. The PUI Schmidt number is much smaller than that of the thermal gas (${\mathrm{Sch}}_{g}\gg {\mathrm{Sch}}_{I}$), thus PUIs and the dissipation associated with them cannot mediate the VLISM shock, as shown by Mostafavi et al. (2017a, 2017b). Consequently, we can neglect PUI dissipation terms in favor of the corresponding thermal gas terms. Unlike the HTS, PUIs will therefore behave adiabatically, ${P}_{I}{\rho }^{-{\gamma }_{I}}$ is constant through a broad VLISM shock, and the thermal gas will be heated preferentially. Recall that PUIs were predicted (Zank et al. 1996) and observed (Richardson et al. 2008) to be heated preferentially at the HTS, with the thermal plasma component remaining relatively unheated. The physical processes determining the structure of the HTS and the VLISM shocks are therefore different.

We use the 1D form of the model equations to determine the structure of the VLISM shock observed by Burlaga et al. (2013). Equations (15)–(17) can be rewritten as a second-order equation,

$$\begin{eqnarray}&&\displaystyle \frac{{d}^{2}y}{d{\overline{x}}^{2}}+\displaystyle \frac{3}{{\mathrm{Sch}}_{g}}\left[\left(\displaystyle \frac{{\gamma }_{I}}{{\gamma }_{g}}\displaystyle \frac{1}{{M}_{s1}^{2}}{\left(\displaystyle \frac{1}{y}\right)}^{{\gamma }_{I}+1}-1\right)-{\mathrm{Sch}}_{g}y+\displaystyle \frac{{y}_{{\rm{B}}}}{{y}^{3}}\right]\displaystyle \frac{{dy}}{d\overline{x}}\\ &&=\,\displaystyle \frac{9}{{\mathrm{Sch}}_{g}}\displaystyle \frac{1}{{M}_{s1}^{2}}(1-y)Z(y),\end{eqnarray} \tag{ 18 }$$

where

$$\begin{eqnarray}Z(y) & \equiv & {M}_{s1}^{2}\displaystyle \frac{{\gamma }_{g}+1}{2}\left(y-\displaystyle \frac{{\gamma }_{g}-1}{{\gamma }_{g}+1}\right)\\ & & -\left(1+{P}_{1}-\displaystyle \frac{({\gamma }_{I}-{\gamma }_{g})}{{\gamma }_{g}({\gamma }_{I}-1)}\displaystyle \frac{(1-{y}^{1-{\gamma }_{I}})}{(1-y)}\right.\\ & & \left.-{\gamma }_{g}{M}_{s1}^{2}\displaystyle \frac{{y}_{{\rm{B}}}}{2y}\left(\displaystyle \frac{{\gamma }_{g}-2}{{\gamma }_{g}}-y\right)\right).\end{eqnarray} \tag{ 19 }$$

Equation (18) is the shock structure equation in the presence of both the thermal gas heat flux and viscosity.

To model the VLISM shock observed by Burlaga et al. (2013), we assume the shock is perpendicular. Upstream of the VLISM shock, the Alfvén speed is V_A1 = 34 km s⁻¹ and the sound speed is C_s1 = 14 km s⁻¹ (V_A > C_s). The magnetoacoustic speed of the ambient medium is about 36.7 km s⁻¹, similar to that estimated by Gurnett et al. (2013), and yields a subcritical shock. The upstream plasma beta β₁ =P_g1/P_B1 ≃ 0.2 < 1. Such a quasi-perpendicular, low beta, low magnetosonic Mach number shock has traditionally been designated a laminar shock (e.g., see Figure 1 of Mellott 1985). Burlaga et al. (2013) showed that Voyager 1 (which is moving relative to the Sun at 17 km s⁻¹) moved past the broad 2012 VLISM shock in about 8.7 days. We need to estimate the shock propagation speed to simulate the VLISM shock. Since a subcritical shock has a magnetosonic Mach number less than about two, the upstream flow speed in the shock frame should be less than 72 km s⁻¹ and therefore, the shock speed in the stationary Sun frame should be less than 52 km s⁻¹ (assuming the post-bow shock/wave VLISM flows toward the Sun with a speed of 20 km s⁻¹; Zank et al. 2013). Kim et al. (2017) used the Multi-Scale Fluid-Kinetic Simulation Suite (MS-FLUKSS) to model shock propagation beyond the HP using near-Earth solar wind data, finding that the 2012 VLISM shock had a compression ratio of about 1.65 and shock speed of about 50 km s⁻¹ in a frame in which the Sun is stationary (T. Kim 2018, private communication). Here we adopt a shock speed of 40 km s⁻¹ with respect to the stationary Sun and, therefore, Voyager 1's traversal of the shock in 8.7 days corresponds to a thickness of about 0.12 au. The thickness of the shock is much greater than that expected for a shock with similar properties in the solar wind at 1 au (Burlaga et al. 2013). The thermal heat conduction length scale, K/U₁, and the viscous length scale, $\eta /{\rho }_{1}{U}_{1}$, associated with proton–proton collisions is about 0.115 au and 0.03 au, respectively.

To model the VLISM shock, we use the upstream parameters M_s1 = 4.17, y_B = 0.31, P_I1/P_g1 = 0.05, γ_I = γ_g = 5/3 and solve Equation (18) in the stationary frame of the shock. Figure 2(a) shows the inverse compression ratio as a function of unnormalized distance (au), being a smooth transition connecting the upstream to the downstream state. The VLISM shock compression ratio is ∼1.67. The shock thickness is about 0.12 au and is determined by proton–proton collisions. Figure 2(b) illustrates the change in magnetic field through the shock. The normalized pressure is shown in the Figure 3. The shock is mediated by the thermal gas and magnetic field pressure, which are the dominant components downstream of the shock. PUIs do not contribute a large downstream pressure. The fast magnetosonic Mach number shows that the flow changes from supersonic to subsonic (Figure 4). Since PUIs do not contribute a large pressure through the shock, their inclusion in the sound speed (i.e., the long wavelength sound speed $=\sqrt{{a}_{g}^{2}+{a}_{I}^{2}};$ Zank et al. 2014) scarcely change the magnetosonic Mach number.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The VLISM is collisional with respect to the thermal plasma since the electron and proton collisional mean free paths are small compared to the almost featureless VLISM. The dominant proton–proton collisions yields a thermal collisional heat flux and viscous dissipation terms. VLISM PUIs are created by secondary charge exchange, experience pitch angle scattering by magnetic fluctuations, and contribute a collisionless heat flux and collisionless viscosity. However, the thermal gas Schmidt number is much greater than the PUI Schmidt number which ensures that we can neglect the PUI dissipation terms. We have shown that the structure of the interstellar shock observed by Burlaga et al. (2013) is determined by interstellar proton–proton collisions, making the shock the first in situ observed example of a classical collisional subcritical shock. Thus, the shock structure is not controlled by wave–particle interactions but by particle collisions. The weak shock is dominated by the magnetic field and thermal gas pressure and PUIs do not contribute a large pressure through the shock transition. The overall thickness of the shock transition is about 0.12 au, which corresponds to the VLISM heat conduction length scale of K/U₁ ∼ 0.115 au.

P.M. acknowledges the support of a NASA Earth and Space Science Fellowship Program grant 16-HELIO16F-0022. G.P.Z. acknowledges the support of NASA grants NNX14AC08G, NNX15AI65G, NNX17AB04G, NNX14AJ53G, and the NSF EPSCoR RII-Track-1 Cooperative Agreement OIA-1655280. We thank Dr. Tae Kim for useful discussions.

If we neglect thermal viscosity, that PUIs behave adiabatically, and that the only dissipation term in the VLISM is the thermal heat flux, then the 1D shock structure equation is now

$$\begin{eqnarray}&&\displaystyle \frac{{dy}}{{dx}}=({\gamma }_{g}-1)\displaystyle \frac{{U}_{x1}}{3{K}_{g}}(y-1)(y-{y}_{A})\displaystyle \frac{M(y)}{N(y)},\end{eqnarray} \tag{ 20 }$$

where

$$\begin{eqnarray}M(y) & = & \displaystyle \frac{1}{2}{y}_{A}{y}_{{\rm{B}}}(1-y)-\displaystyle \frac{1}{2}(y+1){(y-{y}_{A})}^{2}\\ & & +\displaystyle \frac{1}{{M}_{I1}^{2}({\gamma }_{I}-1)}\displaystyle \frac{({y}^{1-{\gamma }_{I}}-1)}{(1-y)}{(y-{y}_{A})}^{2}\\ & & -\displaystyle \frac{1}{{\gamma }_{I}{M}_{I1}^{2}}\displaystyle \frac{{\gamma }_{g}}{{\gamma }_{g}-1}y\displaystyle \frac{({y}^{-{\gamma }_{I}}-1)}{(1-y)}{(y-{y}_{A})}^{2}\\ & & -\displaystyle \frac{1}{({\gamma }_{g}-1){M}_{s1}^{2}}{(y-{y}_{A})}^{2}+\displaystyle \frac{{\gamma }_{g}}{{\gamma }_{g}-1}y{(y-{y}_{A})}^{2}\\ & & +{y}_{{\rm{B}}}(y-{y}_{A})-\displaystyle \frac{1}{2}\left(\displaystyle \frac{{\gamma }_{g}}{{\gamma }_{g}-1}\right){y}_{{\rm{B}}}y(y+1-2{y}_{A});\end{eqnarray} \tag{ 21 }$$

$$\begin{eqnarray}&&N(y)={(y-{y}_{A})}^{3}-{y}_{{\rm{B}}}{(1-{y}_{A})}^{2}-\displaystyle \frac{1}{{M}_{I1}^{2}}{(y-{y}_{A})}^{3}{y}^{-{\gamma }_{I}-1}.\end{eqnarray} \tag{ 22 }$$

Equation (20) admits two solutions. One regime yields smooth transitions from an upstream state to a downstream state, i.e., solutions mediated completely by the thermal gas heat flux. The second regime corresponds to double-valued solutions that are not physical, and a smooth transition from upstream to downstream is impossible. Thermal heat conduction is insufficient to smooth the shock, and an additional thermal dissipation (thermal viscosity) is required for a smooth shock transition. Figure 5 plots ${M}_{s1}^{-2}$ as a function of y_B for quasi-perpendicular shocks (${y}_{A}={B}_{x}^{2}/({\rho }_{1}{U}_{x1}^{2}{\mu }_{0})\to 0$), identifying the smoothed and sub-shock regimes. The observed VLISM shock parameters, ${M}_{s1}^{-2}=0.06$ and y_B = 0.31, show that the 2012 VLISM shock is located in the sub-shock region. Thermal viscosity is therefore necessary to ensure a smooth shock transition.

Zoom In Zoom Out Reset image size