THE VECTOR DIRECTION OF THE INTERSTELLAR MAGNETIC FIELD OUTSIDE THE HELIOSPHERE

After crossing the solar wind termination shock in 2004 and 2007, respectively, the Voyager 1 and 2 spacecraft discovered that the downstream flows remain supersonic with respect to the thermal ions. This configuration had been anticipated (Zank et al. 1996a) and occurs because most (≈80%) of the solar wind energy upstream of the shock transfers to a non-thermal population whose presence, although not directly detected, can be inferred from extrapolation of the available data (Decker et al. 2008; Richardson et al. 2008). The shocked solar wind can then simultaneously exhibit subsonic flow with respect to the energetically dominant component of the plasma and supersonic flow with respect to the thermal component.

Pickup ions likely form the bulk of this non-thermal population. Interstellar neutral atoms, unaffected by electromagnetic fields and possessing a velocity characteristic of the solar system's motion with respect to the local interstellar medium (LISM), ≈25 km s⁻¹, can drift inside the termination shock where they encounter the ≈400 km s⁻¹ solar wind. When these atoms ionize, they suddenly come under the influence of the solar wind's magnetic field, are "picked up," and join the wind's outward flow. The high relative velocity of the neutral atoms becomes an effective non-thermal temperature with k_BT = mv²/2 ≈ 1 keV. Although not directly observed at the termination shock, pickup ions have been detected throughout the inner heliosphere.

Models suggest that the Voyager spacecraft will remain within the heliosphere for 10–20 years (see Table 1 of Opher et al. 2006) before finally encountering the heliopause (HP), the boundary between the solar system and local interstellar space. There, the interstellar magnetic field abuts the heliospheric field in much the same way that, within the solar system, the fields of magnetized bodies meet the interplanetary magnetic field at magnetopauses. Magnetic reconnection often occurs at such interfaces and has been observed at, among other locations, Earth, Mars (Eastwood et al. 2008), and Saturn (Huddleston et al. 1998). By analogy, it is also expected to occur at the HP (Fahr et al. 1986).

However, the high energy content of the pickup ions implies a large value for the heliosheath plasma β (where β = 8πnk_BT/B² is the ratio of the thermal and magnetic pressures). Observations of reconnection at the terrestrial magnetopause (Scurry et al. 1994) and in the solar wind (Phan et al. 2009) suggest that, during high-β conditions, reconnection only occurs between anti-parallel magnetic fields, i.e., when the shear angle between the reconnecting fields is ≈180°. Particle-in-cell (PIC) simulations of reconnection have demonstrated why this occurs: diamagnetic drifts in plasmas with a high β but shear angles <180° suppress reconnection (Swisdak et al. 2003).

In this paper, we present PIC simulations showing that diamagnetic stabilization also occurs at the HP and limits the possible sites of reconnection to those where the interstellar magnetic field lies nearly anti-parallel to the heliospheric field. Large-scale MHD simulations of the heliosphere can then determine where HP reconnection should occur.

We find the locus of anti-parallel reconnection sites to be relatively small, so in situ detection of reconnection signatures at the HP by either Voyager spacecraft seems unlikely. However, on several occasions since their launch, the spacecraft have detected bursts of radio emission at 2–3 kHz believed to be associated with the interaction of interplanetary shocks and the HP (Gurnett et al. 2003). We argue that reconnection triggered by this interaction drives the electron beams producing the radio emission. By combining the locations determined for one set of bursts (Kurth & Gurnett 2003) with the requirement that HP reconnection be nearly anti-parallel, we can constrain the magnitude and orientation of the local interstellar magnetic field. In particular, the vector direction of the field (the difference between B and −B) can be easily determined. This is, to our knowledge, the first determination of this parameter.

As part of our investigation we perform large-scale three-dimensional MHD simulations of the heliosphere. The code is based on BATS-R-US, a three-dimensional parallel, adaptive grid code developed by the University of Michigan (Gombosi et al. 1994) and adapted for the outer heliosphere (Opher et al. 2004).

The Sun is the center of our coordinate system with the x-axis parallel to the line through, but pointing inward from, the nose of the HP. The z-axis lies parallel to the solar rotation axis and the y-axis completes the right-handed triplet. In order to better describe the different sub-environments of the heliosphere, such simulations model multiple fluids (Zank et al. 1996b; Pogorelov et al. 2006); we include five hydrogen populations, one ionized and four neutral (Opher et al. 2009a). All four neutral populations are described by separate systems of Euler equations with corresponding source terms describing neutral-ion charge exchange. The parameters for the inner boundary (located at 30 AU) were chosen to match those used by Izmodenov et al. (2008): proton density at the inner boundary of n = 8.74 × 10⁻³ cm⁻³, temperature T = 1.087 × 10⁵ K, and a Parker spiral magnetic field with strength B = 2μG at the equator. Except where otherwise noted, the heliospheric field orientation corresponds to that during solar cycle 22 (roughly 1986 to 1996), in which in the x–z plane through the nose of the HP, the field in the northern hemisphere points in the $-\mathbf {\hat{y}}$-direction. We denote this as B_SW,y < 0.

The outer boundary conditions are: n = 0.06 cm⁻³, velocity =26.3 km s⁻¹, and T = 6519 K. The neutral hydrogen in the LISM is assumed to have n = 0.18 cm⁻³ and the same velocity and temperature as the ionized LISM. We use fixed inner boundary conditions for the ionized fluid and soft boundaries for the neutral fluids. We impose outflow outer boundary conditions everywhere except the $-\mathbf {\hat{x}}$ boundary where inflow conditions are used for the ionized and neutral populations from the interstellar medium. The grid's outer boundaries are at ±1500 AU in all three directions and the computational cell sizes range from 0.73 to 93.7 AU. The orientation of the interstellar magnetic field is characterized by the angles α_IS, the angle between the field and the interstellar wind, and β_IS, the angle between the field and the solar equator. Given α_IS and β_IS, the vector direction of the field is still undetermined up to a sign that can be fixed by noting whether the y-component of B_ISM is positive or negative.

In Figure 1, we show an overview of a heliospheric simulation with B_ISM,y>0 and the heliospheric polarity of solar cycle 22. The colors represent |B| and the black lines are streamlines of the flow. The roughly circular surface at a radius of 80 AU from the Sun where the magnetic field strength increases from ≈0.02 nT to ≈0.2 nT is the termination shock. Downstream from the shock, the plasma density and temperature increase and the solar wind velocity decreases, as expected. The neutral hydrogen species do not have a dramatic effect on the overall morphology of the system, although they affect the distances to the termination shock and HP, bringing them closer to the Sun, as can be seen by comparison to previous simulations employing only the ionized fluid (Opher et al. 2006, 2007).

Zoom In Zoom Out Reset image size

Reconnection occurs at the grid scale at a rate determined by numerical details rather than by physical processes in these MHD simulations. Previous modeling has shown that without the addition of some sort of extra-physical diffusive process, e.g., an anomalous resistivity, MHD reconnection is inherently slow (Biskamp 1986). The inclusion of the Hall term in the generalized Ohm's law dramatically enhances the reconnection rate (Birn et al. 2001), but necessitates the resolving of length scales on the order of the ion inertial length d_i = c/ω_pi. For a typical HP density of 0.05 cm⁻³, d_i ≈ 10⁸ cm ≈ 10⁻⁵ AU, several orders of magnitude smaller than the resolution of the global simulation. Since the MHD model does not correctly describe reconnection, we use it as a starting point for kinetic simulations of HP reconnection.

3.1. Numerical Methods

3.2. Simulation Results

The components of any plasma that includes both a magnetic field and a non-parallel pressure gradient undergo a diamagnetic drift given by

$$\begin{equation} \mathbf {v}_{*,j} = -c\frac{\bm {\nabla } p_j \bm {\times } \mathbf {B}}{q_jn_jB^2}{\rm,} \end{equation} \tag{ 1 }$$

where p_j = n_jk_BT_j is the pressure and q_j is the charge of species j. Note that, because of their charges, ions and electrons drift in opposite directions and that separate populations of the same species (in particular, pickup and thermal ions) can drift at different speeds. In our case, the pressure gradient normal to the HP crossed with the guide magnetic field produces drifts parallel to the reconnection outflows.

Using PIC simulations, Swisdak et al. (2003) demonstrated that under such conditions the X-line convects in the ion rest frame with a speed given by |v_*e| + |v_*i|, the sum of the electron and ion diamagnetic drifts. If the drift velocity of the X-line exceeds the speed of the nominally Alfvénic outflows, reconnection is suppressed. Qualitatively, suppression occurs when the X-line propagates fast enough that the plasma does not have sufficient time to establish the necessary flow configuration for reconnection before the X-line passes. Quantitatively, Swisdak et al. (2003) proposed that diamagnetic drifts suppress reconnection when

$$\begin{equation} v_{*,j} > v_{A,j}, \end{equation} \tag{ 2 }$$

where v_A,j is the Alfvén speed for a particular population (i.e., calculated based on n_j) in the reconnecting component of the magnetic field. This condition can be derived from balancing the centrifugal and magnetic tension forces exerted on a fluid streaming with velocity v_*,j toward the X-line (see the Appendix for details). When v_*,j is sufficiently large, the magnetic tension cannot eject the field line from the vicinity of the X-line and reconnection stalls.

Equation (2) can be reformulated as a condition relating the jump in the plasma parameter across the current layer, Δβ, and the shear angle θ between the reconnecting fields:

$$\begin{equation} \Delta \beta _j > \frac{2L_p}{d_i}\tan (\theta /2), \end{equation} \tag{ 3 }$$

where L_p represents a typical pressure scale length near the X-line. Both simulations and observations at Earth's magnetopause suggest that $L_p/d_i \sim \mathcal {O}{(1)}$ (Berchem & Russell 1982; Eastman et al. 1996). According to Equation (3), only anti-parallel reconnection, θ ≈ π, occurs for β ≫ 1. Observations of reconnection at the magnetopause (Scurry et al. 1994) and in the solar wind (Phan et al. 2009) support this conclusion.

The inferred large energy content of pickup ions in the heliosheath implies a substantial value for β, which, when coupled with Equation (3), leads to the conclusion that only anti-parallel reconnection occurs at the HP. To check, we performed PIC simulations of HP reconnection at two shear angles. Cuts through the HP at two places in the simulation of Figure 1 established the initial equilibria for these simulations. Due to the relatively large grid size, the HP in the global simulations has a thickness of approximately 10 AU. This is only an upper limit, however, as the real thickness is certainly closer to the 1d_i ≪ 10 AU observed at the terrestrial magnetopause. Furthermore, in both experiments (Yamada 2007) and simulations (Cassak et al. 2005), reconnection is seen to proceed slowly in thick current sheets before accelerating rapidly when the sheet thickness reaches d_i. In our kinetic simulations, the MHD models provide the asymptotic parameters of the reconnecting plasmas, but we take the width of the current layer to be 1d_i. In Figure 2, we show the initial profiles of the densities and temperatures at the two locations used to provide the asymptotic profiles for the PIC simulations.

Zoom In Zoom Out Reset image size

The initial proton distribution is a superposition of a cold, Maxwellian population representing the solar wind and a much hotter Maxwellian (20% by number) corresponding to the pickup ions. Choosing a thermal distribution for the pickup ion distribution is a somewhat crude approximation. The actual distribution is still a subject of investigation, since neither of the Voyager spacecraft can measure it directly. Close to the Sun, where many of the pickup ions are created, a ring distribution could be appropriate, but as the solar wind expands outward, a number of altering effects have been identified (Isenberg 1987; Chalov et al. 1995) even before the interaction with the termination shock introduces further modifications (Zank et al. 2010). In our kinetic simulations, the particle distribution functions are free to evolve (and do so, particularly in the reconnecting current layer), and so the effects of choosing a different initial distribution are ultimately unclear, albeit of obvious interest.

The first site has nearly anti-parallel reconnecting fields (shear angle θ = 165°), while in location 2 the fields have a shear angle of θ ≈ 130°. In order to specify the plane of reconnection for these two-dimensional simulations, we use the prescription derived in Swisdak & Drake (2007) for finding the X-line orientation based on the asymptotic plasma parameters. At locations 1 and 2, the values of Δβ are ≈2 and 4, respectively. Equation (3) then predicts that reconnection should be diamagnetically stabilized at location 2, but not at location 1.

In Figure 3(a), we show a snapshot of the out-of-plane current density from the anti-parallel case. The bulge of the separatrices into the low field strength heliosheath is typical of asymmetric reconnection and has been observed in simulations of the terrestrial magnetopause (Krauss-Varban et al. 1999; Nakamura & Scholer 2000). Since the solar wind has the lower field strength in both the terrestrial and heliospheric systems, the separatrices bulge toward the Sun in both cases.

Zoom In Zoom Out Reset image size

We show a similar image from the non-anti-parallel reconnection case in Figure 3(b). Because of the diamagnetic drift, the X-line has drifted to the right from its initial position at approximately the velocity given by Equation (1). A comparison of the size of the islands downstream from the X-lines in the two simulations indicates that essentially no flux has reconnected in Figure 3(b). The normalized reconnection rate as a function of time for the two simulations can be seen in Figure 3(c). The anti-parallel case (solid line) clearly exhibits reconnection at a rate of ≈0.06, typical of "fast" Hall reconnection (Birn et al. 2001). In the second case (dashed line), where reconnection is not anti-parallel, very little flux reconnects.

The plasma wave instruments on the Voyager spacecraft have detected several bursts of radio emission at frequencies between 2 and 3 kHz during their more than 25 years of operation. A strong interplanetary shock passed by the spacecraft before each burst, which led to the suggestion that the interaction of global merged interaction regions (GMIRs) and the HP boundary produced the emission (Gurnett et al. 2006). An application of several direction-finding techniques to one set of bursts demonstrated that the sources lie at the radial distance of the HP and stretch along a line parallel to the galactic plane and passing near the nose of the heliosphere.

Cairns & Zank (2002) have suggested a detailed mechanism for producing the radio emission. Just beyond the HP, lower hybrid waves can resonantly accelerate ambient electrons that then serve as a seed population for shock acceleration. The enhanced electron beams that form then excite Langmuir oscillations, which subsequently decay, via nonlinear effects, into the observed electromagnetic radiation. This scenario successfully explains, first, why the Voyager spacecraft do not observe radio emissions (or the precursor Langmuir waves and electron beams) when passed by shocks in the inner heliosphere, and second, why the 2–3 kHz signal seems to be radially confined to a thin layer near the HP.

Reconnection at the HP can also explain these observations. Within the heliosphere, the current layers composing the heliospheric current sheet are broad enough (∼100d_i) that, despite occasional observations of X-lines (Gosling 2007), reconnection proceeds slowly on average. Instead, the most significant reconnection should occur when the plasma encounters a boundary, as happens at planetary magnetopauses and, presumably, the HP. The superposition of solar wind structures in a GMIR means that the magnetic field strength can be larger, by a factor of ≈4, than the ambient value (Whang et al. 1995). At the HP, this enhanced field will produce stronger reconnection and therefore naturally accounts for the association between the arrival of the GMIRs at the HP and the bursts of emission. Second, observations and simulations have demonstrated that reconnection can produce both beams of electrons traveling at the electron Alfvén speed (Cattell et al. 2005) and more isotropic populations at relativistic energies (Øieroset et al. 2002; Hoshino et al. 2001; Drake et al. 2006). The electromagnetic decay mechanisms described in Cairns & Zank (2002) can then generate the 2–3 kHz radiation from these electrons. Reconnection, then, naturally localizes the emission at the HP, produces the necessary energetic electrons, and correlates the emission with the arrival of GMIRs.

A significant test of this model is whether the locations of the 2–3 kHz radiation, as determined by Kurth & Gurnett (2003) for a series of events between 1992 and 1994 (during solar cycle 22), correspond to the regions where reconnection should occur. The diamagnetic stabilization of reconnection facilitated by pickup ions described in Section 3 implies that reconnection will only occur where the reconnecting fields are anti-parallel. The locations where such reconnection occurs depend on the orientation of the interstellar magnetic field, which will not be directly measured until one of the Voyager spacecrafts passes the HP. However, inferences from indirect measurements, including backscattered Lyα emission and the difference in the flow directions of interstellar hydrogen and helium (the hydrogen deflection plane), do present some constraints (see Opher et al. 2007, and references therein). Other constraints have been derived based on the heliospheric asymmetries and heliosheath flows measured by Voyager 2 (Opher et al. 2009b). If reconnection is the source of the 2–3 kHz radiation, the requirement that the sources lie near the locus of anti-parallel reconnection can provide a strong further constraint on the vector direction of the interstellar magnetic field.

We used the results of the global MHD model described in Section 2 to compare the locations of the radio sources determined by Kurth & Gurnett (2003) and the locations of anti-parallel reconnection at the HP. We first empirically define the HP in the simulation as the surface where log(T) = 10.9; other values shift the HP somewhat but do not significantly affect our conclusions. Assuming the field in the solar wind, B_SW is purely azimuthal upstream of the HP, and hence neglecting any draping of the heliospheric field, we calculate where $\mathbf {B}_{{\rm SW}}\boldsymbol{\cdot }(\mathbf {B}_{{\rm ISM}} \boldsymbol{\times }\mathbf { \hat{n}})/| \mathbf {B}_{{\rm ISM}}||\mathbf {B}_{{\rm SW}}|$ vanishes, where B_ISM is the field in LISM and $\mathbf {\hat{n}}$ is the normal to the HP. (We use this metric rather than the simpler ($\mathbf {B}_{{\rm ISM}}\boldsymbol{\cdot }\mathbf {B}_{\rm SW})/| \mathbf {B}_{{\rm ISM}}||\mathbf {B}_{\rm SW}| = -1$ to eliminate the effects of spurious normal magnetic fields at the HP.) This quantity vanishes for both anti-parallel and parallel field configurations so we must combine it with another measure (we use the field shear angle) to identify the anti-parallel locations.

In Figure 4(a), we show the map produced by the interstellar field parameters–α_IS = 30°, β_IS = 60°, |B| = 4.4 μG, B_ISM,y>0—that provide the best fit. These parameters are consistent with the values determined by other methods. Our models suggest that the locus of anti-parallel reconnection is more sensitive to α_IS and |B| than it is to β_IS. Changing α_IS by 10° or |B| by 20% moves the anti-parallel locus significantly away from the sources, while shifts in β_IS of up to 20° keep the concordance acceptable. Panels 4(b) and 4(c) display the map for α_IS = 30°, β_IS = 80° and α_IS = 20°, β_IS = 60°, parameters at the edge of what we deem the acceptable range.

Zoom In Zoom Out Reset image size

Perhaps of most interest is the dependence of our results on the vector sign of the local interstellar magnetic field. Figure 4 corresponds to an interstellar field with B_ISM,y>0. Reversing the sign of B_ISM leads to a HP with no anti-parallel reconnection sites near the nose. If anti-parallel reconnection³ causes the 2–3 kHz radio bursts at the HP, then the existence of these bursts fixes the vector direction of the interstellar magnetic field.

The polarity of the heliospheric field reverses due to the 11 year solar cycle. A location that exhibits anti-parallel reconnection during one cycle will, in the next 11 years, find that the heliospheric and interstellar fields are parallel, and hence unlikely to reconnect.⁴ A similar change happens at Earth's magnetosphere when the interplanetary magnetic field rotates (on a much faster timescale) from, for instance, southward to northward. However, just as in the magnetospheric case, even when anti-parallel reconnection is not possible near the nose of the HP it can occur, albeit more weakly, at much higher latitudes near the cusps of Earth's field. This effect cannot be directly tested with the data from Figure 4 because all of the source locations were determined during a two-year period (1992–1994) of solar cycle 22. The Voyager spacecraft did detect radio bursts during solar cycles 21 and 23, but could not determine the locations of their sources, specifically if they occurred at the cusps as our model suggests. Intriguingly, these events were substantially weaker than those of solar cycle 22, which agrees with the observation that in the terrestrial case cusp reconnection is weaker than its equatorial counterpart.

Figure 5 displays results for the solar cycle 23 polarity of the solar dipole (B_SW,y>0) for a case otherwise identical to that shown in Figure 4(c)—β_IS = 60°, α_IS = 20°, |B| = 4.4 μG, B_ISM,y>0. If plotted in the same manner as Figure 4, the entire figure would be white because the interstellar and heliospheric fields have shear angles <90°. Instead, the figure shows that B_SW and B_ISM are nearly parallel over the entire face of the HP. In this case, we expect anti-parallel reconnection to only happen at high heliospheric latitudes; however, the inadequate resolution of the global models there does not allow us to accurately predict the favorable locations.

Zoom In Zoom Out Reset image size

Future data sources may allow us to test this dependence. This model predicts a transition to energetic reconnection near the nose in solar cycle 24 as the polarity of the heliospheric current sheet reverses again. However, current observations suggest that the Sun may be entering a deep solar minimum with very little magnetic activity. If true, it may be accompanied by a concomitant decrease in HP reconnection. In any case, in situ samplings of reconnection by the Voyager spacecraft, aside from being unlikely due to the small odds that a given outward trajectory through the HP will pass near an X-line, will not be possible for at least another decade. Indirect measurements of reconnection at the HP are more likely in the near future. Observations (Lin et al. 2003) and simulations (Drake et al. 2009) have demonstrated that reconnection can energize ions (this occurs in addition to the electric energization necessary for the 2–3 kHz radio emission). When these energetic ions interact with the surrounding plasma, they can undergo charge exchange and create energetic neutrals that can perhaps be sensed remotely by such missions as the Interstellar Boundary Explorer (IBEX). Such a signal would appear in addition to the band of energetic neutrals recently reported by IBEX (McComas et al. 2009) which is, in this view, unrelated to HP reconnection.

We thank the reviewer for helpful comments that led to significant improvements in the paper. Computations were carried out at the NASA AMES Research Center and National Energy Research Supercomputing Center. M.O. and F.A.B. acknowledge the support of NASA-Voyager Guest Investigator grant NNX07AH20G and the National Science Foundation CAREER grant ATM-0747654.

To more formally derive the condition for diamagnetic suppression of reconnection given in Equation (2), begin by considering the outflow from an X-line. In the simplest case, with no diamagnetic drifts, reconnection produces bent magnetic field lines that accelerate away from the X-line due to the $\mathbf {J}\boldsymbol{\times }\mathbf {B}$ force. Specifically, in the coordinate system described in Section 3.1 where outflow is parallel to $\pm \mathbf {\hat{x}}$, the force is proportional to J_zB_y.

In the case with diamagnetic drifts along the $\mathbf {\hat{x}}$-axis, one of the outflows from the X-line will be in the direction opposite to the drifting plasma. The change in momentum that a bent field line can cause in a time Δt over a box of dimensions Δ_x by Δ_y is

$$\begin{equation} J_zB_y\Delta _x\Delta _y\Delta t/c. \end{equation} \tag{ A1 }$$

The unbending of the field line, and hence reconnection, can only occur if this quantity is large enough to overcome the momentum of the plasma in the box traveling with the diamagnetic velocity

$$\begin{equation} -\rho _j v_{*,j}\Delta _y (v_{*,j}\Delta t), \end{equation} \tag{ A2 }$$

where ρ_j and v_*,j are the density and diamagnetic velocity of species j. Combining these expression leads to a condition for the suppression of reconnection. By using Ampere's law to substitute for J_z = (∂B_x/∂y)/(4π/c) and Gauss's law to equate ∂B_x/∂x and −∂B_y/∂y, we arrive at the suppression condition

$$\begin{equation} \rho _j v_{*,j}^2 > \frac{B_x^2}{4\pi }. \end{equation} \tag{ A3 }$$

Equation (A3) is equivalent to Equation (2) of the main text.