TURBULENT HEATING OF THE DISTANT SOLAR WIND BY INTERSTELLAR PICKUP PROTONS IN A DECELERATING FLOW

Interstellar pickup ions originate as neutral atoms in the local interstellar medium, which then penetrate into the region of supersonic solar wind as a result of the slow relative motion of the heliosphere through this medium. When one of these atoms is ionized through an encounter with a solar UV photon or a charge-exchange collision with a solar wind ion, it suddenly appears in the reference frame of the solar wind plasma as an energetic charged particle traveling toward the Sun at the solar wind speed, V_sw. These new ions immediately gyrate about the local magnetic field, forming a ring–beam distribution which is unstable to the generation of MHD waves (Wu & Davidson 1972; Lee & Ip 1987; Bogdan et al. 1991). They subsequently scatter to a stable, nearly isotropic shell in velocity space through the cyclotron resonant interaction with the self-generated and ambient waves. At this point, these particles are comoving with the bulk solar wind and are said to be "picked up." As they accumulate in the outflowing solar wind, they form a distinct energetic population of the plasma (see Zank 1999 for an excellent and extensive review).

The primary component of the inflowing neutral gas is hydrogen, with a density at the solar wind termination shock of N_o ∼ 0.1 cm⁻³ (Gloeckler et al. 1997; Bzowski et al. 2009). The resulting pickup protons come to dominate the internal thermal energy of the plasma far from the Sun and are responsible for large-scale observable effects at heliocentric distances r>20 AU. For instance, the momentum loading by the pickup of these particles causes a continual slowing of the solar wind speed (Wallis 1971; Holzer 1972; Isenberg 1986; Lee 1997), which is measured by the plasma instrument on Voyager 2 (Wang et al. 2000; Wang & Richardson 2003).

The Voyager plasma instrument also measures the temperature of the core solar wind distribution. One infers from the observed properties of the termination shock and the insufficient heating in the core distribution across the shock (Richardson et al. 2008; Richardson & Stone 2009) that the core remains substantially distinct from the accumulating pickup protons throughout the supersonic solar wind. Nevertheless, significant ongoing heating of the core is found far from the Sun (Richardson et al. 1995; Richardson & Smith 2003). At distances of 40–80 AU, the only plausible source for this heating is some process that taps the large pool of energy residing in the interstellar pickup protons.

We have constructed a model of this heating that has had some success in reproducing the observed solar wind temperatures. The basic idea, as first presented by Williams et al. (1995), is that the waves generated through the scattering of new pickup protons to near isotropy feed into a turbulent cascade which is ultimately dissipated by heating the core protons. A phenomenological model of turbulent evolution was developed to follow the behavior of large-scale turbulent eddies in the radially expanding solar wind (Zank et al. 1996; Matthaeus et al. 1999b; Smith et al. 2001). The steady-state turbulence was assumed to be driven by shear and compressional motions in the inner heliosphere, and by the pickup proton fluctuations in the outer heliosphere. The turbulent fluctuations were taken to dissipate at the Kolmogorov rate and the dissipative energy increased the temperature of the core solar wind protons. The core heating predicted from this steady-state model followed the general trend of the Voyager observations. A more exact treatment of the pickup proton scattering and wave generation led to quantitative improvements in this agreement (Isenberg et al. 2003; Isenberg 2005). Most recently, the substantial variability of the solar wind plasma parameters, averaged over a solar rotation, was incorporated into the model to test the comparison even further (Smith et al. 2006).

In Smith et al. (2006), the steady-state model was repeatedly applied for each solar-rotation-averaged observation period at Voyager 2. Each Voyager period was matched to a corresponding earlier period at 1 AU, assuming advection of the solar wind plasma at constant speed between the two radial positions. The average solar wind conditions from the National Space Science Data Center (NSSDC) Omnitape compilation at 1 AU were then used as input to the inner boundary of the model, which was propagated back out to the Voyager position. The variable input to the steady-state equations produced an effective response of the solar wind system to net changes on timescales of a solar rotation or longer. The average effect of shorter variations was still represented by the shear and compressional terms driving the turbulence. Core temperature results equivalent to those of Smith et al. (2006) are shown in Figure 1, where our use of three-solar-rotation time averages here is sufficient to illustrate the model.

Zoom In Zoom Out Reset image size

The average temperatures in Figure 1 show that the model (thin blue line) and the Voyager data (circles) agree reasonably well inside 40 AU, but there are clear discrepancies beyond that point. The disagreement seen between r = 40 AU and 60 AU has been explained by Smith et al. (2006) as due to the latitudinal structure of the solar wind flow during solar minimum. Voyager 2 travels to increasingly high latitude as it moves beyond 30 AU, and the well-organized structure of the wind during solar minimum means that the average solar wind parameters measured at Earth are not a good representation of the wind that reaches the spacecraft. We will justify this interpretation further in Section 3. (The thicker red line in Figure 1 will be discussed in Section 4.)

However, when Voyager passed 60 AU, solar maximum conditions had returned and it is more reasonable to assume that the Omnitape data provide valid input parameters for the model. During that time, the model predicted a core solar wind temperature which is about twice that measured at the spacecraft. It is important to ask what physically justifiable changes can be made to this model to lower the predicted temperature in this region.

Smith et al. (2006) considered several modifications in an effort to address this discrepancy. Various internal parameters that model the turbulent properties could be adjusted, but plausible new values either result in small differences in the distant temperature or in unacceptable disagreements with observations closer to the Sun. A smaller value of N_o can bring the model core temperature in line with the Voyager measurements, but only if this inflowing density is approximately halved. Such a low hydrogen density at the edge of the supersonic solar wind is not consistent with observations (Bzowski et al. 2009). Isenberg (2005) extended the earlier steady-state model to include the additional effect of second-order Fermi acceleration of the pickup protons, which would reduce the wave energy available to drive the turbulence. The consequent reduction in the model core temperature was also shown to be small. Finally, Ng et al. (2010) recently considered the consequences of replacing the Kolmogorov dissipation in the model with an Iroshnikov–Kraichnan formulation. Their initial results were encouraging, but more study is planned. We refer the reader to the cited papers for more detailed discussion of these efforts.

In this paper, we investigate the effect of including the deceleration of the solar wind which must accompany the pickup proton heating. The previous models took V_sw to be constant in r, but momentum loading from the pickup process requires the solar wind to slow down. The speed of a steady-state solar wind will fall linearly with heliocentric radius at large r (Holzer 1972; Isenberg 1986; Lee 1997), and this trend is seen in the Voyager data during solar maximum time periods. Incorporation of the solar wind slowdown is especially appealing in the context of the excess temperatures obtained in our model, since a lower energy per particle (∼V²_sw) in the initial ring beam would result. A lower initial energy could reduce the pickup proton driving of the turbulence, and thus lower the dissipative heating of the solar wind core.

However, in this paper we show that a decelerating solar wind leads to a hotter solar wind proton core, contrary to the above expectation. In the next section, we derive the modifications to the phenomenological turbulence model when the solar wind speed is not constant. Section 3 describes the specific deceleration resulting from the pickup of interstellar hydrogen. Section 4 gives the results of our model calculations when the solar wind deceleration is included, and explains this effect. In Section 5, we discuss these results and present our conclusions.

The evolution of solar wind turbulence can be described by a simplified phenomenological model which follows the intensity of the energy-containing eddies in a radially expanding steady solar wind (Zank et al. 1996; Matthaeus et al. 1999b; Smith et al. 2001, 2006; Isenberg et al. 2003; Isenberg 2005; Breech et al. 2008). These large-scale eddies are driven by sources due to compressions and shears in the inner heliosphere, and by pickup proton isotropization in the outer heliosphere. The turbulent energy is assumed to proceed to dissipative scales through a well-developed cascade, without specifying the details of the nonlinear couplings that actually transport the energy in the inertial range. These models are conceptually similar to the "engineering models" of hydrodynamical turbulence (Bradshaw et al. 1981), which originated in the pioneering work of von Kármán & Howarth (1938).

In these models, the turbulence is characterized by two quantities: Z² = 〈δv²〉 + 〈δb²/4πρ〉, the turbulent intensity at the energy-containing scales in Elsässer units; and λ, the correlation length of the turbulent fluctuations. The cross helicity of the fluctuations is assumed here to be zero, appropriate for the equatorial region of the solar wind beyond ∼10 AU. Extensions of this model to include non-zero cross helicity have been presented (Matthaeus et al. 2004; Breech et al. 2005, 2008), but will not be considered here.

To incorporate a spatially dependent solar wind speed into the phenomenological turbulent evolution model, we return to the earlier analysis of Zhou & Matthaeus (1990) and Matthaeus et al. (1996). The steady-state behavior of the turbulent intensity and the correlation length in a solar wind flow V_sw is described by

$$\begin{equation} {\bm V}_{\rm sw} \cdot \bnabla Z^2 + \left({\frac{{\bnabla \cdot {\bm V}_{\rm sw} }}{2} + M\sigma _D } \right)Z^2 = F_{{\rm shear}} + F_{{\rm pickup}} - \frac{{\alpha Z^3 }}{\lambda } \end{equation} \tag{ 1 }$$

$$\begin{equation} {\bm V}_{\rm sw} \cdot \bnabla \lambda - M\sigma _D \lambda = G_{{\rm shear}} + G_{{\rm pickup}} + \beta Z. \end{equation} \tag{ 2 }$$

The product Mσ_D represents the effects of mixing between the large-scale expansion of the solar wind flow and the small-scale fluctuations. Here, the mixing coefficient M is a geometry-dependent term coupling the spatial inhomogeneity and the turbulent motions (Breech et al. 2008), σ_D is the normalized energy difference between the kinetic and magnetic fluctuation components, σ_D = (r_A − 1)/(r_A + 1), and r_A is the Alfvén ratio of these energies, r_A = E_v/E_b. The F and G terms represent the forcing effects on Z² and λ due to the turbulent sources of shear and pickup. The parameters α and β are set by considerations of local turbulence theory, where we choose α = 2β, corresponding to self-similar decay of the fluctuations (Hossain et al. 1995; Matthaeus et al. 1996).

In the supersonic region beyond 1 AU, we take the solar wind flow to expand spherically in the radial direction. We assume (Zank et al. 1996; Matthaeus et al. 1996, 1999b; Smith et al. 2001) that the turbulence is composed of transverse fluctuations which are uniformly distributed in the plane perpendicular to the average magnetic field direction. We then take that direction to be the Parker spiral angle, Ψ, where tan Ψ = B_ϕ/B_r = rΩ/V_sw, and Ω is the angular rotation frequency of the Sun. Under these assumptions, the mixing coefficient M becomes

$$\begin{equation} M = \frac{{V_{\rm sw} \cos ^2 \Psi }}{r} - \frac{{\cos 2\Psi }}{2}\frac{{dV_{\rm sw} }}{{dr}}. \end{equation} \tag{ 3 }$$

This result corrects the expression given in Zank et al. (1996), and is consistent with the formalism of Breech et al. (2008).

The forcing terms, F and G, are parameterized in the same manner as the earlier versions of this model, taking the compression and shear effects to fall as r⁻¹. We then obtain the equations

$$\begin{eqnarray} &&\frac{{dZ^2 }}{{dr}} + \left({1 + \sigma _D \cos ^2 \Psi } \right)\frac{{Z^2 }}{r} + \frac{{1 - \sigma _D \cos 2\Psi }}{2}\frac{{V^{\prime}}}{{V_{\rm sw} }}Z^2\nonumber\\ &&\quad = C_{{\rm sh}} \frac{{Z^2 }}{r} + \frac{Q}{{V_{{\rm sw}} }} - \frac{{\alpha Z^3 }}{{\lambda V_{{\rm sw}} }} \end{eqnarray} \tag{ 4 }$$

$$\begin{eqnarray} &&\frac{{d\lambda }}{{dr}} - \sigma _D \left({\frac{{\cos ^2 \Psi }}{r} - \frac{{\cos 2\Psi }}{2}\frac{{V^{\prime}}}{{V_{{\rm sw}} }}} \right)\lambda\nonumber\\ &&\quad = - \hat C_{{\rm sh}} \frac{\lambda }{r} - \frac{\beta }{\alpha }\frac{{\lambda Q}}{{V_{{\rm sw}} Z^2 }} + \frac{{\beta Z}}{{V_{{\rm sw}} }}, \end{eqnarray} \tag{ 5 }$$

where V' = dV_sw/dr, and Q is the fluctuation source due to the pickup proton isotropization. The pickup source is modeled as

$$\begin{equation} Q = \zeta \frac{{V_{{\rm sw}} ^2 }}{n}\frac{{dN}}{{dt}}, \end{equation} \tag{ 6 }$$

where n is the solar wind density, dN/dt is the rate at which new pickup protons are created, and ζ is the fraction of initial pickup proton energy which is given up in the isotropization process, generating new fluctuations which drive the turbulence.

These turbulent evolution equations are supplemented by an energy equation to yield heating of the core solar wind protons by the turbulent dissipation

$$\begin{equation} \frac{{dT}}{{dr}} = - \frac{4}{3}\frac{T}{r} - \frac{2}{3}T\frac{{V^{\prime}}}{{V_{{\rm sw}} }} + \frac{m}{{3k_B }}\frac{{\alpha Z^3 }}{{\lambda V_{{\rm sw}} }}, \end{equation} \tag{ 7 }$$

where m and T are the proton mass and core temperature, respectively, while k_B is the Boltzmann constant. This equation is identical to that used in the previous works, with the addition of the middle term on the right-hand side representing the adiabatic heating (or cooling) when the spherically symmetric solar wind flow speed changes with r.

The pickup proton-driven turbulent heating, given by the last term in Equation (7), increases with distance from the Sun and dominates the other effects beyond ∼25 AU. The magnitude of the pickup proton source Q is determined by the details of the pitch-angle scattering of new pickup protons from their initial ring distribution in the azimuthal average spiral field to a closed shell distribution in velocity space. These details are parameterized by the coefficient ζ. In the initial models (Williams et al. 1995; Matthaeus et al. 1999a; Smith et al. 2001), the new pickup protons were assumed to scatter to a bispherical distribution (Galeev & Sagdeev 1988; Williams & Zank 1994), represented by ζ = V_A/V_sw where V_A is the local Alfvén speed. The bispherical distribution, introduced in the context of the intense fluxes found in cometary pickup, corresponds to the velocity-space shell which results when pickup ions efficiently damp the resonant stable modes and only scatter through the interaction with the self-generated unstable waves. In the case of interstellar pickup protons, this assumption was shown to inject too much energy to the turbulent fluctuations, leading to solar wind temperatures far in excess of the observed values (Smith et al. 2001).

The revised analysis by Isenberg et al. (2003) and Isenberg (2005) recognized that the ionization that creates new interstellar pickup protons in the outer heliosphere is much slower than the expected rate for redistribution of wave energy by nonlinear turbulent processes. Thus, a power-law spectrum of dispersive ambient waves was taken to govern the scattering of the new pickup protons. The cyclotron–resonant interaction was modeled using the four parallel-propagating transverse wave modes, right or left circularly polarized, propagating in either direction along the magnetic field. The waves were described by the standard dispersion relation in a cold electron–proton plasma, with intensities proportional to k^−5/3 and equally distributed in all four modes. This "dominant turbulence" analysis resulted in particular shapes for the fully scattered pickup proton shells, which retained more energy in the particles than those of the bispherical distribution, reducing the energy available to drive the turbulence. Additional reductions were obtained by including the effect of non-resonant spreading of the initial pickup ring distribution due to the fluctuating field direction in the ambient turbulence (Isenberg 1996; Németh et al. 2000). The resulting energy partition coefficient, ζ, is a function of V_A/V_sw and of the non-resonant pitch-angle spread, denoted by Δ. This function yielded solar wind heating much closer to the observed values than the bispherical model. This partition function was then used by Smith et al. (2006) in the variable-input model shown in Figure 1. The details of this dispersive, dominant turbulence function ζ are given in Isenberg (2005).

The decrease in solar wind speed by the momentum loading of the interstellar pickup protons can be clearly seen in Figure 2. The circles are the solar wind speeds measured at Voyager 2, averaged over three solar rotations, divided by the averaged, time-shifted speed measured at 1 AU from the NSSDC Omnitape data. In this single figure, the time shift between Voyager 2 and 1 AU is calculated by taking the observed average wind speed at Voyager to be a constant, as in Smith et al. (2006). This procedure neglects the slight effect of the radial deceleration on the choice of the appropriate Omnitape time period. In all the analysis to follow, this deceleration will be included in the determination of the time shift.

Zoom In Zoom Out Reset image size

The large relative speeds found between 40 and 55 AU occurred during solar minimum, when the interplanetary magnetic field was highly organized in heliolatitude. During this time, the solar wind seen at the high-latitude position of Voyager 2 is not well related to the equatorial wind at Earth. Conversely, we see that three-rotation averages of solar wind speed during the remaining time periods are reasonably well correlated over the latitude difference between Voyager 2 and Earth, implying that our purely radial model is justified in taking the Omnitape data for its inner boundary condition.

The dashed line in Figure 2 is the speed ratio obtained from the three-fluid model of Isenberg (1986), taking V_sw = 450 km s⁻¹ at 1 AU, the inflowing hydrogen density N_o = 0.1 cm⁻³, and a proton ionization rate at 1 AU ν_o = 7.5 × 10⁻⁷ s⁻¹leading to an ionization cavity radius of L = 5.6 AU in the upwind direction. These values and the resulting speeds are consistent with similar findings by Wang & Richardson (2003).

We represent this deceleration in the turbulent heating model of this paper by taking the functional form

$$\begin{equation} \begin{array}{ll} V_{\rm sw} (r) = V_{\rm 1AU} - 1.4(r - 18) & r \ge 18\;{\rm AU} \\\\[-6pt] \phantom{V_{\rm sw} (r)}= {\rm }V_ {\rm 1AU} & r < 18\;{\rm AU}, \end{array} \end{equation} \tag{ 8 }$$

where V_sw and V_1AU are in units of km s⁻¹ and r is in units of AU. The dependence (8) is shown in Figure 2 by the solid line. Figure 3 shows the three-solar-rotation-averaged solar wind speed measured at Voyager 2, compared to the predicted speed using Equation (8) to obtain the time-shifted values of V_1AU. Apart from the period of solar minimum, the agreement is very good.

Zoom In Zoom Out Reset image size

We repeatedly integrate the model equations (4)–(7) in r from 1 AU out to the three-solar-rotation-averaged position of Voyager 2, using the time-shifted three-solar-rotation-averaged solar wind conditions at 1 AU from the NSSDC Omnitape as input boundary values for the model. Here, and in all the following analysis, the deceleration (8) is included in the calculation of the time shift. Since we are primarily interested in the heating from the pickup protons, taking place beyond 10 AU, we assume in our model that the average magnetic field is azimuthal, Ψ = π/2, and that the Alfvén ratio is a constant at r_A = 1/2 (Roberts et al. 1987). Finally, we take C_sh = 1.4, $\hat C_{{\rm sh}}$ = 0.1, α = 0.8, and β = 0.4, roughly consistent with the parameter choices of Breech et al. (2008).

In Figure 4, we show the model results for the fluctuating intensity of one transverse component of the magnetic field,

$$\begin{equation} \big\langle \delta B_N^2 \big\rangle = \frac{{2\pi \rho }}{{1 + r_A }}Z^2, \end{equation} \tag{ 9 }$$

and the correlation length, λ, along with the corresponding values obtained from analyzing the Voyager observations, shown as circles.

Zoom In Zoom Out Reset image size

The Voyager magnetic field data are obtained from the NSSDC and the variance of the N-component (in the standard R, T, N heliospherical coordinate system) is computed over 10 hr intervals. These variances are then averaged over 0.5 AU, which provides greater smoothing than the values presented in Smith et al. (2006). Voyager/MAG data are available out to and beyond the termination shock, but we limit this analysis to observations prior to the shock crossing at r = 84 AU. Voyager/PLS observations are analyzed in a similar manner using different codes with secondary averages over three solar rotations rather than 0.5 AU. This accounts for the greater fluctuation level in the temperature plots. Values for the correlation length of the fluctuations are computed using either the integration or the e-folding definition of this length, as described in Smith et al. (2001), and both are shown in Figure 4(b). The correlation length analysis uses 120 hr maximum lag, and these values are then averaged over 5 AU. This analysis was limited to r < 56 AU due to frequent data gaps beyond that distance.

The model core proton temperature results using the decelerating solar wind (8) are given by the thick red line in Figure 1, while the thin blue line shows the results using a constant solar wind speed, V_1AU, for each Voyager point. These results exhibit an increased proton temperature when the solar wind deceleration is included, contrary to the expectation of lower input energy from pickup. The temperature increase is more clearly illustrated using steady input conditions at 1 AU. Figure 5 shows the model temperatures with (thick line) and without (thin line) deceleration for constant 1 AU input values of V_sw = 450 km s⁻¹, V_A = 65 km s⁻¹, $Z^2 {\rm } = {\rm }700$ (km s⁻¹)², λ = 0.027 AU, and T = 7 × 10⁴ K. The larger temperature for decelerating wind is evident.

Zoom In Zoom Out Reset image size

We find that this behavior is due to the sensitive dependence of the partition function ζ on the value of V_A/V_sw, which overcomes the effect of smaller pickup energy per particle in slower wind. This effect is amplified by the fact that, in the azimuthal field of the outer heliosphere, a decelerating solar wind leads to an increasing Alfvén speed. In steady state we have BV_swr = constant and nV_swr² = constant, so V_A ∼ V^−1/2_sw and V_A/V_sw ∼ V^−3/2_sw. Figure 6 shows the value of V_A/V_sw as a function of r for the steady input solutions in Figure 5. The dotted line corresponds to the constant speed model, and the solid line shows the decelerating model. We see that the value of this parameter rises sharply for r>18 AU in the decelerating wind, leading to increased energy transferred to the fluctuations from isotropizing the pickup protons.

Zoom In Zoom Out Reset image size

Of course, a decelerating flow will also be heated by adiabatic compression, given by the second term on the right-hand side of Equation (7). However, we find that the dominant effect here is the increased turbulent driving by the pickup proton scattering. The dashed curve in Figure 6 shows the ratio of the second and third terms of the right-hand side of Equation (7) for the decelerating steady input solution of Figure 5, and it is seen that the compressional heating is less than 10% of the turbulent dissipative heating.

In this paper, we have extended our model of turbulent heating driven by pickup of interstellar hydrogen in the outer heliosphere to include the effect of the solar wind deceleration which must accompany the pickup. One motivation for this step was to see if pickup in the slower wind would reduce the model solar wind temperatures to agree with the Voyager 2 data beyond 60 AU. However, we found that the larger values of V_A/V_sw produced by the deceleration cause enough of an increase in the energy partition function ζ to overcome the reduction in the energy/particle due to the factor of V²_sw. Thus, the energy added to the model turbulence actually increases in the decelerating wind, leading to hotter temperatures.

Consequently, we are still seeking possible changes in this model to improve the predictions of core proton temperatures in the outer heliosphere. Several obvious mechanisms remain to be considered, even within the context of this turbulent evolution picture. For instance, our model currently places all the dissipative energy into heating the core protons. However, if a substantial fraction of this energy were to heat electrons instead (Breech et al. 2009; Cranmer et al. 2010), there could be better agreement in the proton core temperatures.

Another possible improvement follows from recognizing that the pickup proton driving appears in quasi-parallel propagating waves, while the turbulence formalism used here is really appropriate to two-dimensional quasi-perpendicular fluctuations. Such an improved model is being developed along these lines by Oughton and colleagues (Oughton et al. 2006, 2010). This "two-component" model describes the evolution of distinct, but interacting, fluctuating fields: a quasi-two-dimensional turbulent field and a separate "wave-like" field with different phenomenological couplings. We are presently investigating the predictions of this model for describing the observed solar wind properties in the outer heliosphere.

We thank the Voyager/MAG team for data used in this study. We are also grateful for valuable conversations with E. Möbius and S. Oughton. P.A.I. and C.W.S. were supported in part by NSF grant ATM 0635863 and NASA Guest Investigator grants NNX07AH75G and NNX08AJ19G. C.W.S. was further supported by Caltech subcontract 44A1085631 for the ACE/MAG instrument that contributes to the Omni data set. W.H.M. was supported by NASA grant NNX08A147G, and J.D.R. was supported under NASA contract 959203 from JPL to MIT.