Observations of the Outer Heliosphere, Heliosheath, and Interstellar Medium

The Voyager spacecraft have left the heliosphere and entered the interstellar medium, making the first observations of the termination shock, heliosheath, and heliopause. New Horizons is observing the solar wind in the outer heliosphere and making the first direct observations of solar wind pickup ions. This paper reviews the observations of the solar wind plasma and magnetic fields throughout the heliosphere and in the interstellar medium.

A schematic diagram of the interaction of the heliosphere and LISM looking down at the solar equatorial plane. b) The trajectories of V1, V2, and New Horizons and the locations of the V1 and V2 TS and HP crossings Fig. 2 Plasma observations from 1 AU into the LISM plasma speed, density, and temperature from 1 AU into the LISM. Data inside the HP are 25-day running averages from the V2 plasma science (PLS) instrument and the densities outside the HP are from the V1 and V2 plasma wave subsystem (PWS). The radial speed varies between 300 and 850 km/s; variations smooth out at larger distances as the streams interact Elliott et al. 2019).
At 20 AU, V2 observed solar minimum at low-latitudes with slow, dense cool solar wind. After solar minimum a roughly 1.3-year oscillation was observed at V2 and also at spacecraft throughout the heliosphere (Richardson et al. 1994;Gazis et al. 1995). A similar period was observed in cosmic rays, auroras, and convection speeds on the Sun. At about 50 AU, V2 observed the 2006 solar minimum when it was at 20°S heliolatitude. At these higher latitudes a significant fraction of fast wind was mixed in with the slower equatorial wind, giving a rise in speed to over 500 km/s and an increase in temperature.
Outside of 50 AU the solar wind is slowed down by pickup ions (PUIs); speeds are consistently below 400 km/s at both V2 and NH even though the solar cycle activity was quite different 2008a;Elliott et al. 2019). PUIs originate from LISM neutrals that are ionized and accelerated in the solar wind. The PUIs initially have thermal energies comparable to solar wind flow energies and by 15 AU dominate the plasma thermal pressure. The PUI fraction was ∼17% at the TS giving a 18% slowdown of the solar wind, so roughly 1/3 of the solar wind flow energy is transferred to the PUIs before the TS.
The solar wind and heliosphere vary with the solar cycle. At solar minimum the Sun's dipole tilt is small and the coronal holes accelerate solar wind at high-latitudes. Corotating interaction regions (CIRs) form when high speed solar wind streams overtake slower streams, with a forward-reverse shock pair forming at the interface; these are most prevalent in the declining phase of the solar cycle. The largest solar transients are coronal mass ejections (CMEs), which propagate through the heliosphere as interplanetary CMEs (ICMEs). ICMEs expand outward to about 15 AU until they reach equilibrium with the surrounding solar wind (Richardson et al. 2012). ICMEs propagate through the heliosphere to the HP, where they can drive shocks in the VLISM (Kim et al. 2017;Liu et al. 2014;Gurnett et al. 2013). More ICMEs occur at solar maximum. These ICMEs merge, with slow ICMEs overtaking fast ones, to form merged interaction regions (MIRs). Shocks driven by ICMEs are observed in the solar wind out to the TS, with the last large shock, observed at 78 AU (Richardson et al. 2008a), showing jumps in B, V R , density N, a large increase in the temperature T, and enhanced energetic particle fluxes.
The TS is described in detail below; at the TS the speed drops by a factor of 2, the density goes up by a similar factor, and T increases by an order of magnitude (Richardson et al. 2008b). V R drops slowly across the heliosheath. N decreased after the TS at solar minimum, then increased at ∼105 AU when solar maximum reaches the heliosheath, then was fairly constant until the HP. T decreased from 120,000°K after the TS to 50,000°K outside ∼90 AU. The MIRs driven by ICMEs are still observed in the heliosheath (see section below) with two examples at ∼100 and ∼109 AU. Outside the HP, the Voyager PWS instruments determine the plasma density when plasma oscillations are observed (Gurnett and Kurth 2019;. These densities increase away from the HP. PLS data suggest the plasma T is high, 30,000°-50,000°K, after the HP ); T must eventually decline to the LISM value of 6400°K ).

The Energetic Particle Heliosphere
In addition to the plasma heliosphere, Voyager has also observed the energetic particle heliosphere. Figure 3 shows the distance and heliolatitude of V1 and V2 in panel (a), 140-220 keV ions in panel (b), >70 MeV GCRs in panel (c), and the sunspot number in panel (d). The keV particles peak near the Sun, where they are accelerated by shocks driven by ICMEs or corotating interaction regions (CIRs). Their intensities increase with distance, with occasional increases at large events such as the GMIR in Sept 1991. At 65-75 AU the intensities reach a minimum. Then, almost 10 AU before the TS crossing, the intensities increased due to particles accelerated at the TS leaking inwards into the supersonic wind. These particles are observed when field lines connect the TS and the spacecraft; the intensities increase up to the TS.
At the TS, where these particles are accelerated, the intensity jumps, then remains fairly steady in the heliosheath with a slow decrease as the spacecraft approach the HP. Note that even though V1 and V2 are separated by over 100 AU, the intensities are very similar. At the HP the intensities rapidly decrease to background levels (see HP section). The GCRs come from outside the heliosphere and their intensities generally increase with distance. They are modulated by the solar cycle, with less at solar maxima when magnetic fields are higher. At the TS the slope increases, then the intensity takes a final jump at the HP (for more details see Rankin et al. 2022).

The Solar Wind in the Outer Heliosphere
Many features of the evolution of the solar wind are discussed in previous reviews (Burlaga et al. 1996;Richardson and Stone 2009;von Steiger and Richardson 2006;Burlaga et al. 2021b); therefore, we provide only a high-level overview. New Horizons (NH) provides the most recent solar wind data from the outer heliosphere. NH does not have a magnetometer; however, the plasma instrument provides excellent H + and He ++ data and the first pickup ion measurements from the outer heliosphere. Figure 4 shows the V2 and NH solar rotation averaged plasma speeds, densities, temperatures, dynamic pressures, and thermal pressures versus distance. The bottom panel shows NH remains at low heliolatitudes whereas V2 moved south after the Neptune encounter at 30 AU. NH speeds are close to, but slightly below the V2 speeds. The difference grows larger outside ∼43 AU; this increase is probably because V2 observed solar minimum conditions at these distances and, since it was at 20°S, it observed a mixture of fast and slow solar wind. Both spacecraft observed an R −2 density fall-off as expected. The NH densities are slightly higher than those at V2 inside 43 AU and significantly higher outside 43 AU, giving similar fluxes. The NH interstellar pickup ion (PUI) densities shown in purple are fairly constant, so the percentage of PUIs increases with distance; at 50 AU about 8% of the ions are PUIs (McComas et al. 2017. The solar wind thermal ion temperatures decrease to 20-25 AU, then increase due to energy transfer from the PUIs (Smith et al. 2006). The NH PUI temperatures are fairly constant from 10-46 AU at near 4 × 10 6 K; the PUIs are also heated as they move outward (see Zirnstein et al. 2022). The NH thermal pressure plot shows that the PUI pressure dominates the thermal pressure outside ∼12 AU. The physics of the plasma changes when this occurs, with wave speeds determined by the PUIs and shocks preferentially heating the PUIs.

Interstellar Pickup Ion Effects
The first observational evidence for interstellar PUIs were observations of pressure balance structures in the solar wind; outside of 20 AU a hot PUI component was required if the structures were in pressure equilibrium (Burlaga et al. 1994). The increase in temperature outside 20-25 AU discussed above is another effect of the PUIs; the initial PUI ring distributions are unstable and drive magnetic waves which heat the thermal plasma (Smith et al. 2006;Isenberg et al. 2005). The energy to heat the PUIs comes from the solar wind flow energy; slowing of the solar wind is first observed at about 30 AU Wang and Richardson 2003;Elliott et al. 2019) and the slowdown is about 17% ahead of the TS (Richardson et al. 2008a). The solar wind slowdown is related to the percentage of PUIs by n pui /n sw = (3γ − 1)/(2γ − 1)(δv/v 0 ), where n pui is the PUI density, n sw is the thermal solar wind proton density, δv is the solar wind speed decrease, v 0 is the unperturbed solar wind velocity, and γ is the ratio of specific heats (Lee 1995;. For γ = 5/3, n pui /n sw = 7/6(δv/v 0 ) and a 17% slowdown implies that PUIs make up 20% of the solar wind (although we note below that γ may not be constant). Elliott et al. (2019) tested the  equation by using pickup ion and solar wind densities measured with solar wind to calculate n pui /n sw and using the speed at NH and speed measurements at 1 AU from ACE and STEREO to estimate the amount of slowing (δv/v 0 ). They found that when they used a constant polytropic index the two sides of the Richardson et al. 1995 equation did equate. With investigation it became clear that the relationship between the solar wind density and solar wind temperature slowly and systematically evolved with distance. Elliott et al. (2019) fit plots of solar wind density and temperature measurements for given distance ranges using 3 methods. Method 1 fit a simple scatter plot with a line for a given distance range. Method 2 binned the temperature into density bins for a given distance range and fit the binned results. Method 3 fit sets of 8 adjacent points and then binned the polytropic index results by speed for a given distance and determined an index since the index did not vary much by speed. All 3 methods produced a solar wind polytropic index that decreased with increasing distance (Fig. 5). Theoretically, the polytropic index usually refers to the entire plasma and not only the solar wind portion. However, the physical effect encapsulated by the radially decreasing polytropic index determined for the solar wind portion of the data is that the solar wind temperature and density relationship is slowly evolving with increasing distance as the addition of interstellar pickup ions heats and slows the solar wind.
When Elliott et al. (2019) substituted the formula for the polytropic index as a function of distance into the formula derived by , both sides of the equation J.D. Richardson et al. Fig. 6 The average measured solar wind slowing with standard deviations (solid grey line and grey bar) and, in purple, the estimated amount of slowing determined from the NH measured n pui /n sw in the solar wind using the  formula, and the method 3 solar wind polytropic index radial variation with distance determined from solar wind observations ( Fig. 15 from Elliott et al. 2019) agreed and the measured slowing was consistent with the measured PUI density (Fig. 6). Additional measurements are required to quantify how the solar wind is heated and slowed by picking up interstellar material. One needs to measure solar wind parcels as they move away from the Sun, the interstellar neutrals, and the pickup ions simultaneously. Although such measurements do not exist, the interaction between the solar wind and interstellar material can be simulated in a more self-consistent way than is done in the work of  and Elliott et al. (2019) in order to better understand how the solar wind is heated and slowed as interstellar material is picked up by the solar wind on its journey through the heliosphere.

Solar Wind Magnetic Fields Between 7 and 87 AU
This section describes observations of B in the solar wind by Voyager 1 between 7 AU and 80 AU on scales from one to 128 days (Burlaga and Vinas 2005a). Relatively large clusters of strong magnetic fields, MIRs, move past a spacecraft during a fraction of a solar rotation. MIRs form near one AU and can persist through the heliosheath to the heliopause. During the declining phase of the solar cycle, MIRs often occur quasi-periodically, as shown in data from 1980 in Fig. 7(a). More commonly, MIRs occur sporadically, as illustrated in Fig. 7(d) (Burlaga and Vinas 2005a). MIRs can merge to form larger clusters of strong magnetic fields, called global merged interaction regions, GMIRs, that encircle the sun and extend to high latitudes. GMIRs typically form near ∼20-30 AU, and persist throughout the heliosphere and heliosheath as shown in Fig. 1(b).
All of these interaction regions have "filamentary structures" associated with jumps in B of various sizes over a wide range of scales. Jumps in B and filamentary structures, seen in all of the panels in Fig. 7, are fundamental features of the heliosphere and heliosheath, related to the fact that the solar wind is a driven, non-linear, non-equilibrium system. These features are characteristic of the multifractal structure of the solar wind and heliosheath magnetic field (Burlaga 1995(Burlaga , 2001(Burlaga , 2004. Boltzmann-Gibbs statistical mechanics, with an additive entropy function, cannot describe non-equilibrium physical systems with large variability and multifractal structure. Tsallis (1988Tsallis ( , 1994Tsallis ( , 2004a introduced a generalization of Boltzmann-Gibbs statistical mechanics in which the entropy is non-additive (non-extensive) and he introduced a probability distribution function, the Tsallis distribution function, that can describe observations of the magnetic field in the heliosphere and heliosheath. The Tsallis distribution is proportional to the q-exponential function exp q (−B q E) = (1 − (q − 1)B q E) −1/(q−1) . When q = 1 the qexponential function is a Gaussian distribution. The q-exponential function approaches the form of an exponential function for small E and a power law function for large E. The 2nd moment of the q-exponential function is finite for q < 5/3 and diverges for 5/3 ≤ q ≤ 3. The kappa function, which has been used in plasma physics to model speed distribution functions (Olbert 1968;Vasyliunas 1968), is the q-exponential distribution with q = 1 + 1/κ.
One can quantitatively describe the jumps and filamentary structure by studying the normalized increments of B, dBn(t i ) = (B(t i ) + t n ) -B(t i ))/ B(t i ) , on scales t n =2 n days, where n = 0, 1, 2, 3, 4, 5, 6, and 7 and t i is time. We consider scales from 1 to 128 days. The distribution of increments of dBn(t i ) is described by the Tsallis distribution R q (x) = A(1 + (q − 1) B q (dBn) 2 (−1/(q − 1)) , which is the fundamental distribution function in non-equilibrium statistical mechanics (Tsallis 1988). The non-extensivity parameter q for a Gaussian distribution is q = 1, and in general q > 1. Figure 8 shows the distribution of the increments dBn(t i ) (shown by the plus signs) for a given lag in units of days, ranging from one day to 128 days, as shown on the right side of each panel. Each distribution function was fit with the Tsallis distribution shown as solid curves for each lag t n = 2 n days, where n = 0, 1, 2, 3, 4, 5, 6, and 7.
During 1991, Voyager 1 was in the supersonic solar wind at 45 AU. The observed distributions of increments of B (shown by the plus signs) were highly non-Gaussian. (A Gaussian curve would appear as a parabola on the log scale in Fig. 8). This result clearly demonstrates that the large-scale magnetic fluctuations of B in the supersonic wind cannot be described by Boltzmann statistics. But the distributions can be described by the Tsallis distribution of non-extensive statistical mechanics, which are shown by the solid curves. The column on the J.D. Richardson et al.  Fig. 8 shows observations made by Voyager 1 during 2002 near 87 AU, 7 AU before the termination shock crossing. Non-Gaussian (q > 1) Tsallis distributions are observed on scales from 1 to at least 4 days, but at larger scales the distributions are Gaussian. Figure 8 shows data from the CIRs in panel (a) and from a GMIR in panel (b) of Fig. 7. The very large jumps in B/ B are the cause of the non-Gaussian behavior. This behavior is common near 20 -25 AU. However, at larger distances (during 2002 for example) jumps in B were present only at relatively small scales, as expressed by the non-Gaussian fits with a Tsallis distribution on scales from 1 to 4 days. Figure 9 shows the kurtosis for the observations (solid circles) and the Tsallis fits (the exes, x) closer to the sun (panels (a) and (b)) and farther from the Sun (panels (c) and (d)). The solid curves are fits to the observations with an exponential growth curve. The values of the kurtosis were high during 1980 and 1991 when Voyager was within 45 AU of the Sun, decreasing from K = 15 to a plateau for scales greater than 10 days with a relatively large value, K = 7. During 2011 and 2012, when Voyager 1 was at 85 -87 AU, the kurtosis was much smaller, decreasing from ∼K = 6 at a lag of one day to K ∼ 3 (a Gaussian distribution) for lags of 8 to 128 days.
This striking change in the behavior of the kurtosis at larger distances was also observed for the non-extensivity parameter, q. The Tsallis distribution had a finite variance when q < 5/3 and a divergent moment when q > 5/3. For the 1980 and 1991 data (near 8 AU and 45 AU, respectively), q > 5/3 at all scales from one day to 128 days, owing to large tails of the probability distribution functions (PDFs) caused by large fluctuations and jumps in B(t). In Tsallis statistical mechanics, the PDF describes the metastable or quasi-stationary state with the parameter q stat . By contrast, for the 2001 and 2002 data q < 5/3, and q approached 1 (the PDF approached a Gaussian) at scales greater than the solar rotation period. The values of q < 5/3 are a consequence of the smaller tails of the PDFs (smaller jumps in B) during 2001 Fig. 9 The kurtosis of the observations (solid circles) and of the Tsallis fits (x's) closer to the sun (panels (a) and (b)) and farther from the sun (panels (c) and (d)) and 2002. During 2001 and 2002, q tends to approach unity indicating that the distributions are becoming Gaussian. The transition from q > 5/3 at <45 AU to q < 5/3 at >80 AU suggests the possibility of a "phase transition" between 45 AU and 80 AU.
As shown in Fig. 8, the Voyager observations of B showed that the increments of the magnetic field strength by dBn(t i ) = (B(t i ) + t n ) − B(t i ))/ B(t i ) on scales t n = 2 n days, where n = 0, 1, 2, 3, 4, 5, 6, and 7 and t i is time can be described by the Tsallis distribution from 1 AU to 87 AU.
A MHD model which includes pickup ions was used to calculate the radial evolution from 1 AU to 80 AU of the distribution of increments of the field strength B on scales from 1 to 64 days (Burlaga et al. 2007). Plasma and magnetic field data from 1 AU were input into a 1-D MHD model. Figure 10 (left panel) shows the predictions of the model at 80 AU. The distributions of increments dBn of B are shown by the dots and the fit to the observations with a Tsallis distribution is given by the dashed curves. The model distributions are welldescribed by Tsallis distributions on scales t n = 2 n days, where n = 0, 1, 2, 3, 4, 5, and 6, i.e. on scales from 1 to 64 days. The model correctly predicts that there is a significant deviation of the observations from the Tsallis distribution at the largest scale, 2 7 days = 128 days. Thus, the 1-D MHD model predicts that the increments of B at 80 AU are described by the Tsallis distribution on scales from 1-64 days.
The MHD model also predicts that the entropic index q decreases from q > 5/3 at R ∼ 45 AU to q < 5/3 at R ∼80 AU (Fig. 11a), indicating a change from a divergent to a J.D. Richardson et al.

Fig. 10
Model predictions and observations of the magnetic field distributions convergent 2nd moment of the Tsallis distribution, as observed by Voyager 1. This result suggests the possibility of a "phase transition" or a relaxation effect from q > 5/3 at < 45 AU to q < 5/3 at > 80 AU. The model also predicts that the width w n = 1/(dBn) 1/2 of the Tsallis distribution increases linearly with increasing lag (Fig. 11b).

The Termination Shock
The solar wind is super-magnetosonic in the inner heliosphere, so as it approaches the heliopause a termination shock forms where the plasma slows down, is compressed, heated, and becomes subsonic. Voyager 1 crossed the TS in 2004 at 94 AU and V2 crossed in 2007 at 84 AU; these crossings revealed the size of our heliosphere. TS effects were observed upstream of the crossings. Voyager 1 discovered the termination foreshock region starting at 85 AU, in which the energetic particle intensities increase when magnetic field lines are connected to the TS. The particles energized at the TS flow inward along these field lines and populate the foreshock. Voyager 2 entered the foreshock region at 75 AU, about 9 AU before the TS crossing. Electron plasma oscillations were observed weeks prior to the shock, as would be expected in an electron foreshock associated with the TS (see Fig. 14). The foreshock is a dynamic region since changes in the plasma pressure and magnetic field direction change the pathlength between the spacecraft and TS and thus the particle intensities.
V1 crossed the TS during a tracking gap; therefore, it did not observe the actual shock, but several observations indicated that V1 had entered the heliosheath. The magnetic field strength increased by a factor of ∼3, consistent with expectations for the TS, and the variability of B increased dramatically (Burlaga et al. 2005). The low-energy particle intensities jumped at the TS, but the higher energy anomalous cosmic rays (ACRs) did not (Decker   ). The first and last V2 TS crossing occurred in data gaps, but the TS moved across V2 three times while V2 had DSN tracking (Fig. 12), providing the first and only in situ TS data. The first V2 termination shock precursor in the plasma data were the three step-like speed decreases, shown in Fig. 12, that started 90 days before the TS. V R decreased from 400 to 300 km/s, removing ∼40% solar wind flow energy ahead of the shock. The last decrease, from 350 to 310 km/s, was coincident with a low-energy particle intensity increase; the inward pressure gradient from the particle increase is sufficient to produce the reduction in the solar wind speed (Florinski et al. 2009). Thus, the V2 TS is the first example of a particle-mediated shock, where the particles accelerated at the shock move upstream and change the shock structure.
The TS is dynamic and may reform on the scale of several hours. The TS looked like a classic super-critical, quasi-perpendicular shock in 2 of the 3 crossings (Burlaga et al. 2008); one is shown in Fig. 12. It has a foot region, formed by reflected ions, with larger B and lower V. The most rapid changes in the plasma and magnetic field were in the ramp region, which was followed by an overshoot of B and then fluctuations downstream. The high-resolution magnetic field data resolved the ramp structure, which had quasi-periodic fluctuations with length scales of about 1000 km (Burlaga et al. 2008). The other TS crossing was different ( Fig. 12), with two ramp-like features that suggest the TS was reforming (Burlaga et al. 2008). The shock strength was 2-3 and the TS speed was 60-100 km/s, similar to the speeds of planetary bow shocks (Richardson et al. 2008a) A surprise from the V2 TS was that the thermal plasma downstream of the shock was only heated to 10 5 K, compared with expectations of over 10 6 K ( Fig. 13). Only 20% of the solar wind flow energy heated the thermal solar wind plasma; most of the energy heated the PUIs (Richardson et al. 2008a). The solar wind thermal ions gain so little energy they remain supersonic in the heliosheath; the PUIs determine the sound speed, which is subsonic as required. At the TS, almost all the thermal solar wind ions pass directly through the shock potential (Zank et al. 1996), unlike at planetary bow shocks where up to 50% of the thermal ions are reflected (Richardson 1987). The PUI distribution is much broader in energy than the thermal ions and some of the PUIs are reflected. The flow energy of the solar wind goes mainly to the pickup ions. The PUI gain about 65% of the energy, thermal solar wind ions gain 20%, and more energetic ions gain 15% (Decker et al. 2008). However, energetic ions were accelerated only to lower energies. ACRs were not created at either the V1 or V2 TS crossing.
The TS is not symmetric in radial distance. The first evidences of asymmetry were the two crossings of the inner foreshock boundary. V1 entered the foreshock at 85 AU. When V1 crossed the TS at 94 AU, V2 was entering the foreshock region at 75 AU. The V2 TS crossing was at 84 Au, also 10 AU closer than the V1 crossing (Burlaga et al. 2008;Decker et al. 2008;Richardson et al. 2008a;Stone et al. 2008). Changes in the solar wind dynamic pressure account for only 2-3 AU of this difference. Models can reproduce the asymmetry if the VLISM magnetic field were of order 3-5 nT, much larger than the ∼1 nT pre-Voyager estimates but comparable to the observed VLISM B. The observed VLISM B near the HP was larger at V2 than V1, consistent with B being more compressed in the south due to the draping of the magnetic field. The field must also be at the right orientation to the LIC flow (Izmodenov et al. 2005;Opher et al. 2007;Pogorelov et al. 2009) to produce the asymmetry. As of 2021 (35 AU beyond the TS), the V1 magnetic field remains solar-like and the pristine field direction is not known. IBEX observations show the pressure is highest south of the HP nose, consistent with a closer TS and HP in the south (McComas and Schwadron 2014).

Plasma Waves Associated with the Termination Shock
The Voyager PWS instruments observed plasma waves associated with the TS Kurth 2005, 2008). Voyager 1 and 2 both observed electron plasma oscillations prior (upstream) to the TS crossing. At V1, these were observed up to ten months prior to the TS crossing whereas Voyager 2 observed them only about one month prior to the TS. The reason for this variation is not definitively known, but could be due to variations in the relative motion of the TS, the energy of the electron beams underlying the generation of the waves, and the specific geometry of the magnetic field connecting the TS to the spacecraft at the time the waves were observed. For both spacecraft, the plasma oscillations were detected in the 311 Hz spectrum analyzer channel that has a bandwidth of ±15%. The plasma oscillations occur at the electron plasma frequency f pe = 8980(n e ) 1/2 , where frequencies in Hz and electron densities n e in cm −3 . Hence, for 311 Hz, n e is 1.2 × 10 −3 cm −3 and is typical of plasma densities in the solar wind at 90 AU. An example of the Voyager 2 observations is given in Fig. 14 from Gurnett and Kurth (2008).
While Voyager 1 crossed the TS during a data gap, Voyager 2 observed multiple crossings, allowing the detection of plasma waves associated with the shock. An example of the V2 TS plasma waves is shown in Fig. 15 (Gurnett and Kurth 2008). As in planetary bow shocks, the wave signature of the TS is a broad spectrum below the local electron plasma frequency. Voyager does not have a wave magnetic field sensor; however, by analogy with planetary bow shocks, it is likely this feature primarily comprises electrostatic waves. For two of the Voyager 2 TS crossings the broadband signatures corresponded well with the ramp of the shock, suggesting the currents in the ramp drove these waves.

The Heliosheath
Investigating the global nature of the heliosheath, including its thickness and overall pressure balance, is fundamental to understanding the heliosphere as a whole. V1 and V2 found that the heliosheath pressure is dominated by suprathermal particles, which was verified by the Voyager Low Energy Charged Particle (LECP) and plasma subsystem (PLS) data , a combination of Voyager LECP particle data and Cassini/INCA energetic J.D. Richardson et al. The broadband signature starting at 00:11 corresponds to the ramp of the shock and is similar to the spectra observed at planetary bow shocks (from Gurnett and Kurth 2008) neutral atoms (ENAs) data , and recent MHD models (Opher et al. 2020). Rankin et al. (2019) calculate the total effective heliosheath pressure and compare it with IBEX observations to provide additional evidence that the heliosheath dynamics are driven by suprathermal energetic processes.
Voyager 1 moved through the heliosheath at northern heliolatitudes and Voyager 2 was in the south. Although they were both launched in 1977, they arrived at the termination shock at different times, 2004.95 and 2007.49, respectively. The magnetic field magnitudes and directions are shown in Fig. 16 (Burlaga et al. 2021a). The azimuthal angle, λ, of the magnetic field is scattered about 90°and 270°and the elevation angle, δ, is scattered about Fig. 16 Daily averages of the V1 (left) and V2 (right) magnetic field magnitude, azimuthal angle, and elevation angle in the heliosheath 0°, as predicted by the Parker magnetic field model. A sector structure, with alternating inward and outward magnetic fields, was observed by both spacecraft, but Voyager 2 observed a long period of unipolar magnetic fields when it was below the heliosphere current sheet (HCS). The average B did not change until very close to the HP crossing at either V1 or V2. The magnetic field intensity was strongly influenced by solar activity. Voyager 1 observed more MIRs and GMIRs than Voyager 2. Relatively few weak MIRs were observed by Voyager 2 when it was at lower latitudes than the HCS. Relatively strong magnetic fields were observed at Voyager 1 and Voyager 2 for several months as the spacecraft approached the heliopause, when they moved through the magnetic barrier. Figure 17 shows the V2 plasma velocity, density, and temperature across the heliosheath. The average speed |V| is roughly constant at 145 km/s until just before the HP, with excursions between 120 and 180 km/s. V R decreases slowly from 130 km/s near the TS to 80 km/s before the HP boundary region, where it decreases more rapidly and then goes to zero at the HP. V T increases away from the TS; it becomes hard to measure with PLS when the flow angle V surpasses 50°(beyond the field of view of the instrument) but the LECP data show it reaches 100 km/s near the HP. The RT angle increases to a plateau at about 55°in 2012 that extends to the HP. The flow angle out of the solar equatorial plane ρ changes from −10°n ear the TS to −25°near the HP with a maximum flow angle of −30°. The density initially falls by a factor of 2 after the TS in 2009, then increases in 2011. This may be a solar cycle affect since the lower density region corresponds to solar minimum in the outer heliosphere. The average density then stays fairly constant from 2011 until it increases before the HP. The T drops from 150,000°K near the TS to 50,000°K by 90 AU and then stays, on average, constant to the HP. MIRs are observed at 99, 109, and 116 AU with increases in B, V R , N, and T. These are driven by solar transients and propagate through the heliosheath. After the 116 AU MIR, a series of pressure pulses are observed with increases in N, T and the energetic particle intensity. Figure 18 compares the changes in N, T, keV ions, MeV ions, and GCRs. The correlations are very good for the all the parameters, suggesting the heliosheath is traversed by pressure waves with alternating compression and rarefaction regions. Before the HP there is a rise in B, N, T, and a decrease in V R , as discussed in the HP section. These changes suggest that pressure waves move through the heliosheath compressing magnetic flux tubes and the plasma and energetic ions within by similar amounts.

Flow Speeds in the Heliosheath
Since V1 does not have a working plasma instrument, the plasma density and temperature are not known and speeds are derived from particle anisotropies using the Compton-Getting J.D. Richardson et al. Fig. 17 The V2 plasma speed, radial speed, flow angles in the RT and RN planes, density and temperature in the heliosheath (CG) effect (Compton and Getting 1935;Gleeson and Axford 1968;Forman 1970;Ipavich 1974). The CG method uses energetic particle fluxes from different look directions to derive plasma flow speeds. Krimigis et al. (2011) report that V1 entered a stagnation region in 2010, about 8 AU before the HP, where V R went to 0 and was sometimes negative (inward) and V T and V N were small. The V1 CRS instrument observed > 0.5 MeV anisotropies during spacecraft rolls and derived V R consistent with the LECP V R , also showing a stagnation region near the HP ). This stagnation region was unexpected; one possible explanation is that HP instabilities distort the flow field (Borovikov and Pogorelov 2014). An additional puzzle was that, despite the speed decrease, B did not increase and the magnetic flux was not conserved .
The V2 PLS V R profile was also surprising. V R decreased only slowly until the HP crossing. The LECP V2 CG V R profile derived from 28-53 keV ion anisotropies differs from the PLS V R (see Fig. 19); PLS shows a slow steady decrease while LECP shows large, >100 km/s, excursions in V R . Starting in 2016, about 8 AU before the HP, the LECP CG V R averaged about half that observed by PLS.  suggest this slowdown region is analogous to the V1 stagnation region. The V2 CRS V R using >0.5 MeV ions are again consistent with the LECP V R . These authors suggest that the CG calculations are flawed since they disagree with the direct plasma observations; however, another possibility is that PLS speeds are not correct.  show that if the PLS V R profiles are used to calculate the magnetic flux, a conserved quantity (Parker 1963), then the magnetic flux at 1 AU and at V1 and V2 are similar. However, if the CG V R profiles are used then the magnetic flux in the heliosheath differs greatly from that at 1 AU. They also show that if the V1 V R profile were the same as that at V2, then magnetic flux conservation would hold at V1 as well. It is not understood why the different energy particle Fig. 18 A comparison of the changes of plasma density, temperature, keV and MeV heliosheath ions, and GCRs. For each quantity x, (x − x )/ x is plotted where x is the running one-year average. GCR changes are multiplied by 10 but the other plots have no normalization anisotropies from LECP and CRS give the same incorrect V R values, but the V2 speeds and magnetic flux conservation argue against the presence of a stagnation region.

Reconnection in the Heliosheath
The compression of the heliosheath as it approaches the HP brings oppositely directed magnetic field lines together and could drive reconnection in the heliosheath (Lazarian and Opher 2009;Drake et al. 2010;Pogorelov et al. 2009Pogorelov et al. , 2013, but the importance of this process is not clear. Opher et al. (2011) suggest that near the HP the heliosheath is comprised of magnetic bubbles formed by reconnection. Observations of the particle intensities in Fig. 20 show that in the unipolar (no reconnection) zone V2 LECP electron fluxes drop out and ion fluxes decrease (Hill et al. 2014). These fluxes are high in the sector zone where oppositely directed magnetic fields may reconnect and create magnetic bubbles that trap particles and keep the intensities high. Drake et al. (2017) predict B and N are correlated in the sector zone, but not in the unipolar zone, and data substantiate this prediction. However, the agreement of the magnetic flux values at 1 AU and in the heliosheath and the paucity of observations of magnetic D-sheets (Burlaga and Ness 2010;Burlaga et al. 2021b) argue against significant reconnection. The role of reconnection in the heliosheath is still not well understood.

Transients Reaching Beyond the Heliosphere
The Sun is very active and drives transient structures with effects that persist through the heliosphere and into the LISM. The largest of these events are coronal mass ejections (CMEs) Fig. 19 Left: the top panel shows the V2 radial speed measured by PLS (black) and derived from the CG effect using LECP (blue) and CRS (red) data. The middle and bottom panels compare the magnetic flux at 1 AU (Omni data set) with that calculated at V2 using the PLS and LECP CG speeds, respectively. Right: the top panel shows the V1 radial speed derived from the CG effect using LECP (blue) and cosmic ray subsystem (CRS) (red) data. The middle and bottom panels compare the magnetic flux at 1 AU (Omni data set) with that calculated at V1 using the PLS speeds and using the V2 LECP CG profile, respectively The white regions are where V2 was in a unipolar magnetic field region (below the heliospheric current sheet) and the gray regions are where V2 was in the sector region and saw both toward and away magnetic field directions. V1 is always in the sector region (Hill et al. 2014) which propagate through the heliosphere as interplanetary CMEs (ICMEs). ICMEs expand outward to about 15 AU until they reach equilibrium with the surrounding solar wind (von Steiger and Richardson 2006). ICMEs have been tracked through the heliosphere to the TS. These ICMEs can drive shocks which persist out to the HP (Kim et al. 2017;Liu et al. 2014;Gurnett et al. 2013). Corotating interaction regions (CIRs) form when high speed solar wind streams overtake slower streams, with a forward-reverse shock pair forming at the interface.
Both ICMEs and CIRs accelerate particles; ICMEs are more likely to drive large events, MIRs and GMIRS, that influence the outer heliosphere. More ICMEs occur at solar maximum. These ICMEs merge, with faster ICMEs overtaking slower ones, to form MIRs. These MIRs are quasi-periodic in the outer heliosphere near solar maximum, occurring roughly twice a year, with pressure enhancements of factors of 5-10. These MIRs run into the TS and move it outward, sending pressure pulses through the heliosheath. Figure 17 shows that these pressure pulses are characterized by highly correlated increases in the plasma density, pressure, and keV to MeV particle intensities as the whole region is compressed. The exception is B, whose magnitude is only weakly correlated with the pressure changes. These heliosheath pressure pulses eventually reach the HP, push the HP outward, and drive shocks into the LISM. These shocks generate electron beams which drive the plasma oscillations observed by PWS.

Thickness of the Heliosheath
The thickness of the heliosheath from the TS to the HP was 28 AU at V1 and 36 AU at V2. Steady state MHD simulations give widths of 55-65 AU (see Kleimann et al. 2022, this journal). This difference implies that time dependence is important or that pressure is missing from the heliosheath; mechanisms suggested to resolve this problem are solar cycle variations (Izmodenov et al. 2005(Izmodenov et al. , 2008, removal of hot heliosheath ions by charge exchange (Malama et al. 2006;Opher et al. 2020), inclusion of thermal conductivity (Izmodenov et al. 2014), and escape of ACRs across the HP (Guo et al. 2018). All produce insufficient heliosheath thinning, as discussed in Kleimann et al. (2022, this journal)

Heliosheath Magnetic Fields and Plasma: Solar Maximum
To contrast the magnetic observations in the heliosheath with those in the supersonic solar wind we discuss V2 observations from 2015, when V2 moved from 106.81 AU to 109.93 AU, latitude 30.9°S to 31.3°S and longitude 217.7°to 217.9° . The average B during 2015 was 0.126 nT, which is higher than the heliosheath average, because in 2015 Voyager 2 was observing solar maximum. At least 3 MIRs were observed during the first 140 days of the year and a GMIR was observed from day 260 to at least day 305. It is noteworthy that such strong interaction regions, with B exceeding 0.2 nT, survive out to distances approaching that of the heliopause.
In the GMIR in Fig. 21, the density and temperature were enhanced significantly and the speed increased, suggesting that the GMIR was still growing in strength. The GMIR was related to a large increase in the dynamic pressure, NV 2 , that might have produced a shock in the VLISM. The GMIR produced a large decrease in >70 MeV/nuc GCRs, as typically observed through the heliosheath and heliosphere beyond 10 AU (Burlaga et al. 1985). Figure 22 shows that the distribution of the magnetic field strength during 2015 was lognormal, which is characteristic of the distribution of B that is observed when the sun is active and produces ejecta and magnetic clouds that interact with one another to produce the strong magnetic fields in MIRs and GMIRs. The distribution of the azimuthal angle λ has two nearly equal height peaks at λ = 90°and 270°, corresponding to the Parker spiral magnetic field sector directions. The two peaks of λ indicate that the heliospheric current sheet extended to large latitudes in the southern hemisphere, relative to the latitude of Voyager 2, and was warped, producing four sectors in the magnetic field. The well-defined sector boundary that V2 crossed near day 70, 2015 was thick; the change in the direction of B lasted ∼7 days. A minimum variance analysis showed that B rotated smoothly through a current sheet in a plane that was tilted 47°with respect to the solar equatorial plane. The velocity did not change when Voyager 2 crossed the current sheet, consistent with a tangential discontinuity. The distribution of the elevation angle, δ, peaked at δ ∼ 0°, as expected for a spiral magnetic field. There was a change in δ to near −65°during an extended interval between days 150 and 210. These unusual angles might be attributed to a ripple in the The small-scale structure of the V2 2015 magnetic field can be quantified by plotting the increments of B(t) on a scale of one hour, dB1h (t) = B (t + 1 hr) − B(t). Figure 23(b) shows that dB1h(t) for the one-hour of B was very bursty. The dots in Fig. 22(a) show the distribution of the increments of dB1h(t). The Tsallis distribution provides an excellent fit to the distribution of increments (the coefficient of determination R 2 = 0.999), with a non-extensivity parameter q = 1.66 ± 0.03. Similar results for the distribution of daily averages give a fit to the Tsallis distribution with R 2 = 0.979 and q = 1.60 ± 0.17. Thus, on these small scales, the non-extensively parameter during the undisturbed conditions shown in Fig. 23 is similar to that observed during the quiet solar wind conditions in 2000 discussed earlier.

Heliosheath Magnetic Fields and Plasma: Solar Minimum
The non-extensive statistical mechanics of (Tsallis 1988(Tsallis , 2004a(Tsallis , 2004b(Tsallis , 2009) provide a description of the driven, open, non-equilibrium systems in the supersonic solar wind and heliosheath. This section discusses V2 observations of the 2010 "quiet" heliosheath magnetic field, when the spacecraft moved between 91.02 and 94.5 AU at latitudes 28.8°S to 29.3°S . The magnetic fields were carried by solar wind that left the southern coronal hole on the sun during 2008 and 2009, when solar activity was historically low (Ahuwalia and Ygbuhay 2011) and B at 1 AU was only 4.0 nT. Voyager 2 was south of the heliospheric current sheet (which was relatively close to the solar equatorial plane, since it was solar minimum) and it measured a relatively weak average magnetic field of 0.08 ± 0.04 nT. Thus, Voyager 2 observed a "minimum energy state" of the heliosheath during 2010. J.D. Richardson et al. The full width at half maximum was 0.060 nT for daily averages and 0.079 nT for the hourly averages. Both B distributions in Fig. 24 are accurately described by a Gaussian distribution (with coefficients of determination R 2 = 0.94 and 0.97, respectively), which is characteristic of B in the solar wind and the heliosheath when solar activity is very low. Tsallis (2004a,b) introduced the concept of the "q-triplet" = (q stat , q sen , q rel ) associated with the "non-extensively index" q stat obtained from the Tsallis distribution of increments of B, the parameter q rel related to the correlation coefficient, and the parameter q stat derived from the multifractal spectrum. Burlaga and Vinas (2005a) were the first to compute the q-triplet for a physical system, using solar wind data. They obtained the values q stat = 1.75 ± 0.06, q sen = −0.6 ± 0.02, and q rel = 3.8 ± 0.034 for observations made in the solar wind near 30 AU in 1989 and near 85 AU in 2002. They showed that the q-triplet had the properties associated with non-extensive statistical mechanics. For example, the deviation of the 3 q's from unity is a measure of the departure from thermodynamic dynamic equilibrium. Nevertheless, the fluctuations in the solar wind and heliosheath are in a metastable, quasistationary state.
The one-day increments of B, dB1d = B(t + 1 day) − B(t), plotted in Fig. 25, show that the increments are "spiky" or intermittent. An excellent fit (R 2 = 0.98) to the observations (the squares in Fig. 25b) was obtained with the Tsallis (q-Gaussian) distribution, shown by the solid curve, and the value q = 1.6 ± 0.01 was obtained from the fit. This value of q is the same as that observed during solar maximum conditions during 2015 discussed in Sect. 2, even though in 2010 the sun was inactive, with few sunspots! As a measure of the significance of q = 1.6, we note that q = 1 for a Gaussian distribution and q = 1.7 is associated with the onset of chaos in the z = 2 logistic map (Tsallis 2009(Tsallis , p. 1964. The non-extensively parameter q is seldom much larger than 1.6 in the heliosheath and solar wind. The parameter q derived from the distribution of increments of B is the same as the parameter q stat . In Boltzmann Gibbs statistical mechanics, the PDF describes the thermal equilibrium state characterized by the temperature T. In Boltzmann Gibbs statistical mechanics there is an exponential relaxation time to thermal equilibrium (exponential decay, with a relaxation time τ ). Tsallis statistics are associated with the q-exponential relaxation of microscopic quantities toward thermal equilibrium (q-exponential decay), with a relaxation parameter q rel that can be computed from the "correlation function" C(τ ) = q rel = (B(ti + τ ) − B(ti) ) × ((B(ti + τ ) − B(ti) ) / (B(ti) − B(ti) )2 that is plotted as a function of τ on a log-log scale in Fig. 25. The result is accurately described (R 2 = −0.99) by a linear fit with a slope s = −0.35 ± 0.02 on scales from 1 to 16 days. Thus, we find that q rel = 1 − 1/s = 3.9 (Tsallis 2009, p. 156). This result is in good agreement with the results that q rel = 3.98 in 1989 and q rel = 3.54 in 2002 obtained by Burlaga and Vinas (2005b). In Boltzmann-Gibbs statistical mechanics, the correlation function is an exponential, rather than a power law, related to q cor in Fig. 26 as predicted by non-extensive statistical mechanics which we find in the Voyager 2 data.
In Boltzmann Gibbs statistical mechanics there is exponential sensitivity to the initial conditions (strong chaos, described by an exponential growth characterized by zero Lyapunov exponents and a growth parameter q sen ). In Tsallis statistics the system is related to q-exponential sensitivity to the initial conditions described by q-exponential growth. This parameter is derived from the "multifractal spectrum" (Burlaga et al. 1993(Burlaga et al. , 2006Halsey et al. 1986;Hentschel and Procaccia 1983;Macek 2007;Macek and Szczepaniak 2008;Nauenberg 2003;Pirraglia 1993;Stanley and Meakin 1988;Tel 1988;Tsallis 2004b, 2009, and Tsallis and Brigatti 2004. The multifractal spectrum of daily averages of B observed by Voyager 2 during 2010 was computed using the methods described by Burlaga (1995) and Burlaga et al. (2006), and the results are plotted as points in Fig. 27 for 1 < τ < 8 days, where f is the multifractal distribution function, which is a function of the parameter α. The relevant quantities are the values of α at f = 0, but the observational points do not extend to f = 0. The observations were extrapolated to f = 0 by means of the quadratic and a cubic fit. Using the cubic fit one obtains the minimum αm = 0.73 ± 0.03 and the maximum αM = 1.35 ± 0.02, similar to the values obtained by Burlaga and Ness (2012) in the heliosheath. From these numbers one obtains q sen from the relation 1/(1 − q sen ) = 1/α min −1/α max (Lyra and Tsallis 1998), hence q sen = −0.5 ± 0.3, compared with the average value for the distant heliosphere q sen = −0.6 ± 0.2 obtained by Burlaga and Vinas (2005a).
In summary, during 2010 Voyager 2 observed a 1) q-Gaussian distribution of daily increments of B with q = 1.6, 2) a power law correlation (rather than an exponential correlation) on scales from 1 to 16 days with q rel = 3.9, and 3) a multifractal structure of B on scales from 1 to 8 days with q sen = −0.5 ± 0.3.) Figure 28   shows that the length of the curve (B(t), τ ) (which is treated as a fractal) as a function of lag τ in a log-log plot is a straight line (R = ±0.93) on scales from τ = 1 to 100 days. Thus, the length of the curve is a power law L ∼ τ −s , which corresponds to a frequency spectrum ∼ f −α where α = 3 − 2s = 1.1 ± 0.10 consistent with a "1/f" spectrum on scales from 1 to 100 days. This spectrum corresponds to an equal distribution of energy at all scales. This 2010 observation is the first demonstration of an f −1 spectrum associated with B in the heliosheath .

The Heliopause Region
The HP is the boundary between the solar wind plasma and the interstellar plasma. The location and nature of the HP were unknown before the Voyager crossings; it was thought to be a sharp tangential discontinuity in analogy with planetary bow shocks. This section describes the HP and its precursors upstream in the heliosheath. Figure 29 shows overviews of the HP regions at V1 and V2. V1 crossed the HP on day 238 of 2012 at 121.7 AU. V2 crossed it on day 318 of 2018 at 119.0 AU. The first clear evidences of the V1 HP crossing were an abrupt increase in magnetic field strength B, a decrease of heliosheath energetic particles, and an increase in the GCR counting rate at 2012.56 Krimigis et al. 2013;). These changes were followed by a recovery, a larger abrupt change, another recovery, then a final crossing after which B remained high and steady, the energetic particles disappeared, and the GCR intensities plateaued (see region T near 2012.6 in Fig. 29). Figure 31 shows that the direction of B did not change across the HP, contrary to expectations that the VLISM field would be tilted with respect to the solar field. The steady V1 B direction caused controversy about whether this boundary was the HP or a transition from closed field lines to field lines connected to the LISM, allowing heliosphere ions to escape into the LISM and GCRs to enter the heliosphere. The strong uniform magnetic fields, absence of energetic particles, the plateau of the GCRs, and the subsequent observation that the density was 0.08 cm −3 (Gurnett et al. 2013) were strong evidence that V1 had crossed the HP, but other hypotheses persist (Gloeckler and Fisk 2015). The V2 HP crossing had similar features, an increase in GCR intensity, an increase in B, and a decrease in heliosphere ion intensity . The magnetic field direction again did not change. V2 has a working plasma instrument which observed a sharp change in the plasma flux at the HP . The higher densities in the VLISM were again confirmed by PWS .
V1 and V2 observed HP boundary regions and complex structure before the HP crossings. Figure 30 compares the V1 and V2 HP crossings to illustrate which features are common and which may depend on time and/or position. The first precursor of the HP may be the decrease in CG V R speed at both V1 and V2 shown in Fig. 18 about 8 AU before the HP (Krimigis et al. 2013. As discussed above, neither B nor the PLS V R at V2 change at this boundary; the cause of the CG speed decrease is not known but may be a HP boundary effect. J.D. Richardson et al.

The Heliopause Boundary Region and Magnetic Barrier
The V1 HP crossing is preceded by the HP boundary region, a ∼1.3 AU wide region that started at day 125 of 2012 with a step increase in the GCR rates and a decrease in heliosheath electrons (Fig. 29) and ended at the HP with a second step increase of GCRs, an increase in B, and a dropout of low-energy ions to near background levels Krimigis et al. 2013;. The magnetic field lines in this region may be reconnected to the VLISM magnetic field enabling GCRs to enter the heliosheath and heliosheath electrons to exit.
V2 observed a similar boundary region with an increase in the GCR intensity and decrease in heliosheath electron flux on day 229 of 2018, so open field line boundary regions may be a common feature at the HP. At V2, the GCR increase was coincident with the beginning of the magnetic barrier, a region of increased B not observed at V1. The average B in the V2 magnetic barrier was ∼ 0.40 ± 0.06 nT, significantly greater than the average B of ∼0.13 nT from days 1-229 of 2018. The magnetic barrier had the strongest magnetic fields observed in the heliosheath and persisted for 80 days, from day 229 to day 309; these fields were comparable in strength to the magnetic field in the VLISM observed by Voyager 1. Neither the plasma parameters nor the LECP ion intensity changed at the start of the magnetic barrier. The V2 magnetic field increased to 0.7 nT at the HP, significantly larger than the 0.4 nT field observed at V1. The V2 magnetic field was predicted to be stronger than that at V1 to account for the difference in the TS locations.

The Plasma Boundary Region
V2 provided the first plasma data from the HP region. Figure 31 shows a plasma boundary layer began on day 150 of 2018, 1.5 AU ahead of the HP, where the plasma speed decreases by 30%, the density increases by a factor of 2, and T increases by 30%. The LECP 53-85 keV ion intensity increase correlates well with the plasma density suggesting the whole plasma flux tube is compressed, similar to the other pressure pulses observed in the heliosheath. The GCR slope increased slightly at the beginning of this region (Fig. 31), but the magnetic field did not change at the start of the plasma boundary layer but did decrease slowly from days 150-215 as N and T increased. The origin of the plasma boundary region is not known; compression of plasma as it approaches the HP could cause an increase in the density, temperature, and LECP ion flux as observed, but would also produce an increase in B which was not observed.

The Heliopause Transition Layer
The HP transition layer started on day 302, 2018, about 0.06 AU inside of the HP. In this boundary layer the density increases by a factor of two, the magnitudes of V R and V N decrease, the flow angle in the RT plane increases, and B increases by 30% . This region contains solar wind plasma that is modified by the HP boundary. The V2 HP was a sharp boundary where the outward plasma flux dropped to background levels, B increased, the GCR intensity increased, and the heliospheric particle intensity decreased in the 8-hour data gap on day 309. Figure 32 shows the particle intensity profiles across the HP . At V1 several HP precursors were observed which looked like partial HP crossings or encounters J.D. Richardson et al. Fig. 31 The plasma boundary layer. The panels show daily averages of the plasma speed, density, temperature, B, and GCR intensity. The dashed lines show the beginning of the plasma boundary layer, the start of the magnetic barrier, and the HP (from  with flux tubes from the VLISM moving into the heliosheath. After the V1 HP, ACRs were observed in the VLISM for about 25 days. At V2 the leakage from the heliosheath into the VLISM was much more pronounced, especially at higher energies. A drop off in ACRs occurred 65 days after the HP but these ions did not drop to background levels until about 120 days after the HP. The GCRs at V1 very quickly approached their asymptotic value, but at V2 the GCR increase took about 30 days. The region where magnetic field lines are connected outside the HP was much larger at V2 than V1, perhaps because V2 was further from the nose of the heliosphere.

The Heliopause
The V1 HP crossing was marked by an increase in the magnetic field strength, a decrease of heliosheath energetic particle intensities, and an increase in the GCR counting rate. After the HP, B remained high and steady, the heliosheath particles disappeared, and the GCR intensities plateaued. However, the direction of the magnetic field did not change. PWS data confirmed that this boundary was the HP, observing the higher densities, ∼0.06 cm −3 , expected in the VLISM (Gurnett et al. 2013). Several HP precursors were observed at V1, with smaller decreases in B and the energetic particles and increases in the GCRs centered on days 212 and 230 of 2012 Krimigis et al. 2013;. The precursors may be flux tubes moving from the VLISM into the heliosheath.
The V2 crossing in Fig. 30 did not have precursors like those at V1. On day 309 of 2018, B sharply increased, the heliosheath energetic particle intensity decreased, the GCR counting rate increased, and the radially outward plasma currents dropped to background levels. PWS observations from days 35-55 of 2019 confirmed the density had increased to The biggest surprise of the V1 HP crossing was that the direction of B did not change Burlaga and Ness 2014). At V1 the magnetic field direction near the HP was nearly constant but different from the Parker field direction of 270°by about 20°in azimuth and 18°in elevation angle. Models disagreed on whether this lack of B rotation was a coincidence of the geometry where the HP was crossed or if the rotation of the VLISM B toward the Parker spiral direction were an intrinsic HP feature. This issue was resolved (Fig. 29) when B did not change direction at the V2 HP crossing and thus this lack of a rotation in B seems a usual, but not understood, HP characteristic .
At V2, the magnetic field azimuthal angle was very close to 270°and the elevation angle was about 20°. The strong uniform magnetic fields, the absence of energetic particles, the plateau of the GCRs and the cessatiopn of radial flow were strong evidence that V2 had crossed the HP. The lack of field rotation is still a puzzle and discussed in the modeling chapter (Kleimann et al. 2022). 6 The LISM

The LISM Temperature
The PLS instrument normally observes LISM currents in only one energy channel in one detector. This gives one observation but three unknowns, the plasma density, temperature, and speed into the cup. However, there are two times PLS can estimate the plasma temperature; when PWS measures the density and when the spacecraft rolls, giving data from different look directions. When the density is known from PWS, the measured currents in the 10-30 eV ion detector limit T and V to a 2D box in TV space. Using a reasonable range J.D. Richardson et al.

Fig. 33
Radio emissions were observed by both Voyager PWS instrument some time after the Saturn encounters. They consist of a low-frequency component (< 2.4 kHz) that shows little frequency drift and a series of higher frequency emissions that drift in frequency with time. Gurnett et al. (1993) developed the basic explanation of these emissions as radio waves mode converted from electron plasma oscillations at the electron plasma frequency in the source. Further, they argued for the source plasma density to be high enough to explain the frequencies of the emissions, the source had to be beyond the heliopause where a strong increase of density would be expected to balance the pressure from the much lower density and hot plasma in the heliosheath for V, T is estimated to be 30,000-50,000 K. When the spacecraft rolls PLS (sometimes) observes the variation of the currents with look angle. Fitting the currents vs. angle gives T, V, and N; the fit T is again on the order of 35,000 K . The temperatures are on the upper end of those expected from MHD models (Zank et al. 1996) and may suggest heating by reconnection in the PDL (Fuselier and Cairns 2017).

The Electron Density in the VLISM
The Voyager PWS instruments began making remote measurements of the density of the VLISM as early as 1983 through the detection of radio emissions as shown in Fig. 33 in the frequency range of 1.8 to 3.5 kHz from the vantage point of heliocentric radial distances of greater than about 13 AU (Kurth et al. 1984(Kurth et al. , 1987. At the time of these measurements, the nature of the source was not understood, although sources at or beyond the heliopause were hypothesized. The basic explanation for these emissions was formulated by Gurnett et al. (1993) as a second group of radio emissions was detected beginning in 1992. Both the events starting in 1983 and in 1992 showed similar structures, a low-frequency (<2.4 kHz) component with little drift in frequency and a higher frequency component displaying frequency drifts of order 1 kHz/yr. Another, significantly weaker event commenced in 2002. The temporal spacing of the events was similar to the 11-year solar cycle, which implied a solar influence in the occurrence of the emissions. The 1983 (1992) event was first detected within 412 (408) days of a very strong Forbush decrease. A Forbush decrease is a decrease in the flux of galactic cosmic rays at Earth caused by disturbances in the interplanetary magnetic field that impede the access of the cosmic rays to the inner heliosphere. This  Gurnett et al. (1993) to explain the frequency drift of the radio emissions detected by Voyager inside the heliosphere. The propagation of a shock through a density gradient beyond the heliopause causes electron plasma oscillations at the plasma frequency that mode couple into radio emissions and propagate into the heliosphere, where they are detected by the Voyagers correlation led to the hypothesis that transients from the Sun propagate through the heliosphere, form GMIRs which hit the HP and drive shocks that propagation into the VLISM. Using observed shock speeds from in situ observations of transients at a range of heliocentric distances by Ulysses, Pioneers 10 and 11, and the Voyagers, and some simple assumptions about the thickness of the heliosheath and propagation speeds therein, a time-of-flight calculation gave an estimate for the distance to the heliopause of 116 to 177 AU (Gurnett et al. 1993). This estimate encompasses the distance of the observed heliopause crossings in the range of 120 AU by both Voyagers. Gurnett et al. (1993) explain the frequency drift of the radio emissions as owing to the propagation of the shocks through a density gradient beyond the heliopause (Fig. 34). The shocks accelerate electron beams which excite electron plasma oscillations at the local electron plasma frequency. The plasma oscillations then mode-convert into radio emissions that propagate freely away from the source. These radio waves were detected by Voyager outside 10 AU. Assuming the primary radio emission is at the electron plasma frequency f pe , the frequency of the radio emission gives the plasma frequency at the source. Hence, the frequency drift is a result of shocks propagating through the VLISM. With n e related to the electron plasma frequency (Hz) by n e = (f pe /8980) 2 , the VLISM electron densities at the radio sources range from about 0.04 cm −3 to 0.15 cm −3 . Voyager PWS measurements in the VLISM substantiated the earlier interpretations of the kHz radio observations. The Voyagers were at the right time and location to observe electron plasma oscillations driven by shocks propagating into the VLISM from solar transients. As of 2022, eight intense electron plasma oscillations events have been observed by Voyager 1 as shown in Fig. 35 (Gurnett and Kurth 2019). The frequencies of these emissions show that the VLISM electron density increases from 0.04 cm −3 near the heliopause to 0.12 cm −3 20 AU beyond the heliopause. Electron plasma oscillations are associated with an electron foreshock not unlike those ahead of planetary bow shocks in the solar wind; however, these shocks are associated with solar transients that are transmitted through the heliopause and propagate through the VLISM (Gurnett et al. 2015). These have been observed at a rate loosely averaging once per year.
The V2 PWS has detected two plasma oscillation events since it crossed the heliopause in 2018. The V2 instrument no longer produces high resolution wideband data, but the PWS spectrum analyzer channels allow determination of the VLISM plasma density (Kurth and Gurnett 2020). The two V2 measurements are plotted with PLS and PWS measurements J.D. Richardson et al.

Fig. 35
Detections of electron plasma oscillations by Voyager 1 after it crossed the heliopause. The increasing frequency of the events with increasing time and radial distance shows a radial gradient in the VLISM (Gurnett and Kurth 2019). These plasma waves were predicted to be present in the VLISM in the early 1980's and are responsible for the observed radio emissions from upstream of the termination shock up to beyond 146 AU in Fig. 36. A clear radial density increase is observed by both Voyagers despite the broad separation in their locations. Figure 37 shows that starting in 2015 a very weak line has been observed in the V1 PWS wideband observations at the electron plasma frequency that is likely due to thermal plasma oscillations (Burlaga et al. 2021b;Ocker et al. 2021). The densities obtained from the plasma oscillation events and the thermal plasma line are given in Fig. 36. As noted by Burlaga et al. (2021b), the plasma densities from the thermal plasma line show density increases at two pressure fronts recognized in the Voyager 1 magnetometer data. The ratios of the plasma densities and the magnetic fields before and after the shock in 2014 and at both pressure fronts are remarkably similar.

The Propagation of Solar Transients into the VLISM
As discussed above, solar transients are the primary drivers of the plasma oscillations observed in the VLISM. CMEs and CIRs are the two primary types of solar transients. Note that the propagation time of a solar transient from the Sun to the heliopause is several hundred days, during which many solar transients can occur. Therefore, interactions between the solar transients during the transit time and formation of an MIR occur in the outer heliosphere. A radio wave event observed by Voyager in the VLISM is most likely caused by merging of multiple transients.
MIRs in the outer heliosphere can form from the merging of a series of ICMEs (e.g., Burlaga 1995;Wang and Richardson 2002;. The Voyager measurements in the VLISM enable timing analysis of the propagation of solar transients into the VLISM. The radio and shock events observed at Voyager 1, as shown in Fig. 34, also allow investigations of the ultimate destiny of solar transients as they approach the outer heliosphere and how they affect the VLISM. The April-May 2013 radio wave event in Fig. 35 helped determine that Voyager 1 had crossed the heliopause, providing the first evidence that the density had increased to the ex- Fig. 36 A summary of plasma densities measured by the V2 PLS and the V1 and V2 PWS instruments beginning just inside the termination shock. The two plasma oscillation events detected by V2 after it crossed the heliopause show a radial density gradient remarkably similar to that observed by V1. Given that the two Voyagers are separated by 70°as viewed from the sun, the radial gradient is a large-scale feature of the VLISM, at least toward the nose of the heliosphere (Kurth and Gurnett 2020) pected VLISM levels. This event is hypothesized to have been produced by the 2012 March CMEs hitting the heliopause (Gurnett et al. 2013). Liu et al. (2014) made the first attempt to establish the timing of the propagation of solar transients into the VLISM, combining J.D. Richardson et al.

Fig. 38
MHD propagation of solar wind streams from the Earth to 120 AU. The upper left panel shows the observed solar speed at Wind, and the other panels show the predicted speeds at various distances. The shaded regions in the first panel represent the ICME intervals at the Earth, and the shaded region in the last panel indicates the period of the radio emission observed by Voyager 1. The predicted magnetic field is plotted at 80 AU to show the size of the MIR. Reproduced from Liu et al. (2014) wide-angle imaging observations from STEREO, in situ measurements, and MHD propagation of the measured solar wind disturbances. In 2012 March, the active region NOAA AR 11429, one of the biggest active regions in solar cycle 24, exhibited extraordinary activity (Liu et al. 2014). The active region emitted a series of large CMEs as it rotated with the Sun from the east to west. A cluster of shocks and ICMEs were observed near the Earth, and their propagation outward was predicted using an MHD model (see Fig. 38). The transient stream interaction results in the formation of a large MIR preceded by a shock in the outer heliosphere. The predicted arrival time of the shock and MIR at 120 AU is April 22 2013, which agrees with the April-May 2013 radio emission period and the time of a transient disturbance in galactic cosmic rays (Gurnett et al. 2013;Krimigis et al. 2013).
This MHD model does not include the transition across the termination shock and the heliosheath. Liu et al. (2014) use passage of a shock through the Earth's bow shock and magnetosheath as an analogy, to determine the speed of the shock in the heliosheath. Subsequent MHD simulations are fully three-dimensional (3D) and multi-fluid (e.g., Fermo et al. 2015;Kim et al. 2017;Guo et al. 2021). These recent models have the advantage of including the termination shock, heliosheath and heliopause, but face difficulties such as how to determine the 3D solar wind input, and reproduce correct locations of the heliospheric structures and precise timing with the Voyager 1 measurements. Nevertheless, Kim et al. (2017) reproduce some of the shocks in the magnetic field observed at Voyager 1, assuming an ad hoc 3D solar wind input with a prescribed polar coronal hole geometry. Their results indicate that merging of CIRs may also have played a role in the formation of some of the MIR shocks. Observationally, the timing of events at V1 and V2 can also be compared . V1 crossed the heliopause into the VLISM in 2012 when Voyager 2 was still in the heliosheath. From 2012.5-2016.5, solar maximum conditions persisted in the heliosheath with five MIRs observed at Voyager 2. These MIRs occur at a similar frequency to the transients observed in the VLISM by Voyager 1. Figure 39 shows that the timing between Voyager 1 and 2 measurements indicates that the transients observed in the VLISM by Voyager 1 may have been driven by the MIRs observed at Voyager 2.

Summary
The Voyager spacecraft provided the first survey of the heliosphere from the Earth into the LISM, making the first measurements of the outer heliosphere where PUIs dominate the thermal pressure, of the termination shock, the heliosheath, heliopause, and the VLISM. These observations provided many surprises and have left us with many puzzles for future missions to answer. The Voyager spacecraft have sufficient power to operate all instruments until the mid-2020s; after this time, the instruments will be turned off serially which will extend the useful life of the spacecraft to at least 2030. New Horizons is providing another view of the outer heliosphere and may cross the TS and enter the heliosheath, making the first direct measurements of the PUIs. IBEX (and soon IMAP) provide a global view of the heliosphere through ENA measurements. But more missions are needed! J.D. Richardson et al. Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.