SOLAR WIND PROTON, ELECTRON, AND ION HEATING

# Proton, Electron, and Ion Heating in the Fast Solar Wind from Nonlinear Coupling Between Alfvénic and Fast-Mode Turbulence

Steven R. Cranmer and Adriaan A. van Ballegooijen Harvard-Smithsonian Center for Astrophysics, 60 Garden Street, Cambridge, MA 02138, USA
###### Abstract

In the parts of the solar corona and solar wind that experience the fewest Coulomb collisions, the component proton, electron, and heavy ion populations are not in thermal equilibrium with one another. Observed differences in temperatures, outflow speeds, and velocity distribution anisotropies are useful constraints on proposed explanations for how the plasma is heated and accelerated. This paper presents new predictions of the rates of collisionless heating for each particle species, in which the energy input is assumed to come from magnetohydrodynamic (MHD) turbulence. We first created an empirical description of the radial evolution of Alfvén, fast-mode, and slow-mode MHD waves. This model provides the total wave power in each mode as a function of distance along an expanding flux tube in the high-speed solar wind. Next we solved a set of cascade advection-diffusion equations that give the time-steady wavenumber spectra at each distance. An approximate term for nonlinear coupling between the Alfvén and fast-mode fluctuations is included. For reasonable choices of the parameters, our model contains enough energy transfer from the fast mode to the Alfvén mode to excite the high-frequency ion cyclotron resonance. This resonance is efficient at heating protons and other ions in the direction perpendicular to the background magnetic field, and our model predicts heating rates for these species that agree well with both spectroscopic and in situ measurements. Nonetheless, the high-frequency waves comprise only a small part of the total Alfvénic fluctuation spectrum, which remains highly two-dimensional as is observed in interplanetary space.

magnetohydrodynamics (MHD) — plasmas — solar wind — Sun: corona — turbulence — waves
journal: The Astrophysical Journal, 2012, in press

## 1. Introduction

The energy that heats the solar corona and accelerates the solar wind originates in convective motions beneath the Sun’s surface. However, even after many years of investigation, the physical processes that transport a fraction of this energy to the corona and convert it into thermal, magnetic, and kinetic energy are still not understood. In order to construct and test theoretical models, a wide range of measurements of relevant plasma parameters must be available. In the low-density, open-field regions that reach into interplanetary space, the number of plasma parameters that need to be measured is larger because the plasma becomes collisionless and individual particle species (e.g., protons, electrons, and heavy ions) can exhibit divergent properties. Such differences in particle velocity distributions are valuable probes of kinetic processes of heating and acceleration.

The spectroscopic instruments aboard the Solar and Heliospheric Observatory (SOHO)—e.g., the Ultraviolet Coronagraph Spectrometer (UVCS) and Solar Ultraviolet Measurements of Emitted Radiation (SUMER)—have measured several key collisionless plasma properties for a variety of solar wind source regions (Kohl et al., 1995, 1997, 2006; Wilhelm et al., 1995, 1997). These observations augment decades of in situ plasma and field measurements that show similar departures from thermal equilibrium in the collisionless solar wind (e.g., Neugebauer, 1982; Marsch, 1999, 2006; Kasper et al., 2008). In the high-speed solar wind, both coronal and heliospheric measurements point to the existence of preferential ion heating and acceleration, as well as protons being hotter than electrons. There are also marked departures from Maxwellian velocity distributions for protons and other ions, with the temperature measured in directions perpendicular to the background magnetic field often exceeding the temperature parallel to the field (i.e., ).

A large number of different processes have been suggested to explain the measured proton and ion properties. Many of these processes are related to the dissipation of magnetohydrodynamic (MHD) waves, and many involve multiple steps of energy conversion between waves, reconnection structures, and other nonlinear plasma features. It was noticed several decades ago that the damping of ion cyclotron resonant Alfvén waves could naturally give rise to many of the observed plasma properties (see reviews by Hollweg & Isenberg, 2002; Hollweg, 2008). The problem in the solar corona, though, is how these extremely high-frequency ( Hz) waves could be generated from pre-existing MHD fluctuations that appear to have much lower frequencies ( Hz).

One likely source of high-frequency waves and kinetic dissipation is an MHD turbulent cascade. There is ample evidence that turbulence provides substantial heat input to the plasma in interplanetary space (see Coleman, 1968; Goldstein et al., 1995; Tu & Marsch, 1995; Matthaeus et al., 2003). Furthermore, self-consistent models of turbulence-driven coronal heating and solar wind acceleration have begun to succeed in reproducing a wide range of observations without the need for ad hoc free parameters (e.g., Suzuki & Inutsuka, 2006; Cranmer et al., 2007; Rappazzo et al., 2008; Breech et al., 2008; Verdini et al., 2010; Bingert & Peter, 2011; van Ballegooijen et al., 2011; Chandran et al., 2011). The general scenario is that convection jostles open magnetic flux tubes that are rooted in the photosphere and produces Alfvén waves that propagate into the corona. These waves undergo partial reflection, and the resulting “colliding wave packets” drive a turbulent cascade which heats the plasma when the eddies reach small enough spatial scales.

It has been known for many years that Alfvénic turbulence in a strong magnetic field produces a cascade to small scales mainly in the two-dimensional plane perpendicular to the field (Montgomery & Turner, 1981; Shebalin et al., 1983), and thus is not likely to produce high-frequency ion cyclotron waves. In other words, MHD turbulence leads to eddies with large perpendicular wavenumbers and not large parallel wavenumbers . Under typical plasma conditions in the corona and inner heliosphere, the linear dissipation of high- Alfvén waves would lead to the preferential parallel heating of electrons (Leamon et al., 1999; Cranmer & van Ballegooijen, 2003; Gary & Borovsky, 2008). This apparently disagrees with the observational evidence for perpendicular heating of positive ions.

There have been several proposed solutions to the apparent incompatibility between the predictions of MHD turbulence and existing measurements (see also Cranmer, 2009a). For example, turbulent fluctuations may be susceptible to various instabilities that cause ion cyclotron waves to grow (Markovskii et al., 2006; Vranjes & Poedts, 2008) or they may induce stochastic perpendicular motions in ions if they reach nonlinear magnitudes (Voitenko & Goossens, 2004; Wu & Yang, 2007; Chandran, 2010). Nonetheless, heliospheric measurements have provided several pieces of evidence for the existence of ion cyclotron resonance that gives rise to perpendicular ion heating in the solar wind (e.g., Marsch & Tu, 2001b; Bourouaine et al., 2010; He et al., 2011; Smith et al., 2012). The most direct solution to the problem still appears to be for turbulence to transport some fraction of the fluctuation energy to high- cyclotron resonant waves.

The goal of this paper is to investigate the idea proposed by Chandran (2005) for the turbulent generation of ion cyclotron waves. In this scenario, nonlinear couplings between Alfvén waves and other modes such as fast magnetosonic waves produce an enhancement in the high- power-law tail of the Alfvénic fluctuation spectrum. This is made possible by the ability of fast-mode waves to cascade nearly isotropically in wavenumber space. Thus, the gradual nonlinear generation of ion cyclotron waves may provide enough heat to protons and other ions in the corona and inner solar wind (see also Luo & Melrose, 2006; Chandran, 2008a; Yoon & Fang, 2008).

We note that it is not currently possible to produce a rigorous model that contains a fully self-consistent description of MHD wave transport (from the corona to 1 AU), turbulent cascade, mode coupling, and dissipation. In order to make some progress in trying to understand this complex system, we have created models that include a range of simplifying assumptions. One key approximation is that we divide the modeling into two separate components: (1) a large-scale model of the radial dependence of fluctuation energy densities, and (2) a small-scale description of how the “local” fluctuations at each radius evolve in wavenumber space and heat the plasma. Feedbacks from the latter to the former are not included, and we discuss their potential importance in Section 7.

We model the plasma conditions in a representative magnetic flux tube that is rooted in a polar coronal hole and that exhibits a steady-state fast solar wind outflow. In Section 2 we describe a model of background plasma conditions and large-scale wave transport in this flux tube. We take an empirical approach to the solar generation of Alfvén, fast, and slow mode MHD waves by specifying their amplitudes as free parameters at a lower coronal boundary height of 0.01 solar radii () above the photosphere. Section 3 gives a summary of how we model the small-scale transport of cascading wave energy in wavenumber space, and Section 4 describes our treatment of the nonlinear coupling between high-frequency Alfvén and fast-mode waves. In Section 5 we apply quasilinear kinetic theory to predict the net rates of particle heating from the cascading waves. Section 6 presents a selection of results for the collisionless rates of proton, electron, and heavy ion heating. Finally, Section 7 concludes this paper with a brief summary of our major results, a discussion of some of the wider implications of this work, and suggestions for future improvements.

## 2. Large-Scale Model of Coronal Hole Conditions

We wish to better understand the global energy budget of MHD waves and turbulence from the lower solar corona out to the interplanetary medium. The work of this section builds on many earlier models of the radial evolution of Alfvén waves in the fast solar wind (e.g., Hollweg, 1986; Tu & Marsch, 1995; Cranmer & van Ballegooijen, 2005; Chandran & Hollweg, 2009) and extends it to describe the likely behavior of fast and slow magnetosonic waves as well. Below, we describe an empirical model of how the time-steady plasma properties vary with heliocentric distance (Section 2.1) as well as a large-scale view of the dispersion, propagation, and dissipation of linear waves in such a system (Sections 2.22.4).

We model the plasma properties along an open magnetic flux tube rooted in a polar coronal hole. At solar minimum, large unipolar coronal holes are associated with superradially expanding magnetic fields and the acceleration of the high-speed solar wind. Because we only consider a field line along the polar axis of symmetry, we do not need to include the rotational generation of azimuthal magnetic fields (e.g., the Parker spiral effect; see Weber & Davis, 1967; Priest & Pneuman, 1974) or other geometrical effects of streamer-like flux tube curvature (Li et al., 2011). We do not distinguish between dense polar plumes and the more tenuous interplume regions between them. The radial dependence of plasma parameters is described as a function of either the heliocentric distance () or the height above the solar photosphere ().

To specify the radial variation of the time-steady magnetic field strength , mass density , and solar wind outflow speed , we used the empirical description of Cranmer & van Ballegooijen (2005). This model combined a broad range of observational constraints with a two-dimensional magnetostatic model of the expansion of thin photospheric flux tubes into a supergranular network canopy. At AU in this model, the solar wind outflow speed is 781 km s and the proton density is 2.56 cm. This model also specifies the Alfvén speed , which decreases from a maximum value of 2890 km s at down to 31 km s at 1 AU. There is a local minimum in at that is the result of the assumed shape of network “funnels” that expand superradially into the corona.

We also need to know the plasma temperature in order to determine the relative importance of gas pressure versus magnetic pressure. Despite observational evidence for different particle species having different temperatures (and departures from Maxwellian velocity distributions), we generally assume that the majority proton-electron magnetofluid is close enough to thermal equilibrium that strong plasma microinstabilities are not excited (e.g., Gary, 1991; Marsch, 2006). Thus, we specify a one-fluid temperature that is assumed to be equal to both the proton temperature and the electron temperature , and we assume temperature isotropy () for both species.

We used the polar coronal hole model of Cranmer et al. (2007) as a starting point to describe , but this model was modified in two ways. First, we moved the sharp transition region (TR) down from a height of 0.01 to 0.003 solar radii () to better match the conditions of semi-empirical models (e.g., Fontenla et al., 1990; Cranmer & van Ballegooijen, 2005; Avrett & Loeser, 2008). Thus, in the adopted model, at the temperature has risen to 0.48 MK, and it continues to rise to a maximum value of 1.36 MK at . We also increased the temperature slightly at distances greater than 0.2 AU in order to better agree with the mean of the in situ and measurements of Cranmer et al. (2009). At AU, MK and it declines as . The one-fluid sound speed is defined as , where is the monatomic ratio of specific heats, is Boltzmann’s constant, and is the hydrogen atomic mass.

Figure 1 shows the radial dependence of a selection of the background plasma properties defined above. It also shows the dimensionless plasma beta parameter, which is usually defined as the ratio of gas pressure to magnetic pressure, with

 β0=PgasPmag=2γ(csVA)2. (1)

However, we will often use a simpler dimensionless parameter given by

 β=(csVA)2=γβ02, (2)

where and differ only by a factor of 1.2 when . The range of heights shown in Figure 1 extends down into the solar chromosphere, but the wave models discussed below start at a lower boundary condition in the low corona; i.e., they specify the wave and turbulence properties only for .

### 2.2. Linear Properties of MHD Waves

In this section we briefly summarize the dispersion properties of linear MHD waves (i.e., phase and group speeds for the Alfvén mode and the fast and slow magnetosonic modes) and the partitioning between fluctuations in kinetic, magnetic, and thermal energy. In Sections 2.32.4 we assume that all three types of MHD waves are present, and we vary their relative strengths arbitrarily in order to match the observations.

The phase speed is defined in terms of the frequency and the magnitude of the wavenumber . In general, is a function of the Alfvén speed, the sound speed, and the angle between the background field direction and the wavevector . We follow the standard convention of defining a Cartesian coordinate system with the background magnetic field along the axis and the vector having components only in the - plane. Also, for now we express and in the frame comoving with the solar wind. For Alfvén waves,

 V2ph=V2Acos2θ (3)

and for the magnetosonic modes,

 V2ph=V2A+c2s2(1±Σ) (4)

applies with the upper sign corresponding to the fast mode and the lower sign corresponding to the slow mode, and with

 Σ=√1−σcos2θ,σ=4β(1+β)2 (5)

(see, e.g., Whang, 1997; Goedbloed & Poedts, 2004). In Section 2.3 we also need to know the component of an MHD wave’s group velocity in the direction parallel to the background magnetic field. We call this quantity , and for the Alfvén mode it is identically equal to no matter the value of . For the fast and slow modes,

 Vgz=Vphcosθ(1∓σsin2θ2Σ(1±Σ)). (6)

MHD waves excite oscillations in the plasma parameters. We denote the root-mean-square (rms) fluctuation amplitudes in velocity as , in magnetic field as , and in density as . We ignore fluctuations in the electric field because their contribution to the total energy density tends to be negligible when . The kinetic, magnetic, and thermal energy densities associated with each type of fluctuation are given as

 Ki=ρ0v2i2,Mi=B2i8π,Θ=βB208π(δρρ0)2 (7)

respectively, with . For linear Alfvén waves, the total energy density is divided equally between transverse kinetic and magnetic fluctuations along the axis, with and

 KyUA=MyUA=12. (8)

For fast and slow mode waves,

 UF,S=Kx+Kz+Mx+Mz+Θ (9)

and we follow Whang (1997) in expressing the partition fractions as follows,

 KxUF,S=fnsin2θ+ftcos2θ,KzUF,S=fncos2θ+ftsin2θ (10)
 fn=V2ph(V2ph−V2A)c2sΔ,ft=V2A(V2ph−c2s)cos2θV2phΔ (11)
 MxUF,S=(V2ph−c2s)cos2θΔ,MzUF,S=(V2ph−c2s)sin2θΔ (12)
 ΘUF,S=V2ph−V2AΔ (13)

where . The fast and slow velocity fluctuations () always occupy exactly half of the total energy density, and the combination of magnetic and thermal fluctuations () take up the other half.

The energy partition fractions given above are familiar components of plasma physics and MHD textbooks (e.g., Stix, 1992; Goedbloed & Poedts, 2004). However, it is difficult to see intuitively how these fractions vary throughout the heliosphere from Equations (10)–(13) alone. Thus, in Figure 2 we provide a schematic illustration of the energy partitioning for fast-mode waves. The three columns indicate the variation from low () to medium () and high () beta plasmas. The three rows show the results for purely parallel propagation (), an isotropic distribution of wavenumber vectors (see below), and purely perpendicular propagation (). In general, all five terms on the right-hand side of of Equation (9) are nonzero, but fractions less than 1% are not shown in Figure 2. This diagram can be transformed to show the properties of slow-mode waves by replacing with and interchanging the and subscripts with one another.

In order to determine how the total energy density of a given wave mode evolves with heliocentric distance, we solve equations of wave action conservation that contain multiple sources of wave damping. There have been many discussions of energy conservation for both pure acoustic waves and incompressible Alfvén waves (e.g., Dewar, 1970; Isenberg & Hollweg, 1982; Velli, 1993; Tu & Marsch, 1995; Verdini & Velli, 2007; Sokolov et al., 2009), but general derivations that can also be applied to fast and slow mode waves (for arbitrary ) are less frequently seen. We utilize the results of Jacques (1977) to write the damped wave action conservation equation as

 ∂∂t(UmΩ)+1A0∂∂r(⟨u0+Vgz,m⟩A0UmΩ)=−QmΩ (14)

where the subscript can be replaced by A, F, or S for the relevant mode, is the cross sectional area of the flux tube (i.e., ), and is the total dissipation rate for the mode in question. The dimensionless factor that takes account of the “stretching” effect of wavelengths in an accelerating reference frame is

 Ω=⟨Vph,mu0cosθ+Vph,m⟩ (15)

and the angle brackets denote a weighted average over all angles,

 ⟨f⟩=∫dθsinθf(θ)∫dθsinθ, (16)

where here we consider outward propagating waves with . The factor of in Equation (15) comes from the difference between the wave frequency in the Sun’s reference frame () and the comoving-frame frequency () that appears in the definition of the wave action; see Section III of Jacques (1977).111We note that the adopted form of Equations (14)–(15) is only one out of several possible ways of placing and grouping the angle brackets. For the fast and slow modes, there is also potential ambiguity about whether one should use Lagrangian or Eulerian averages for in the transport equation. In future work we will explore the consequences of different methods of averaging. Equation (14) implicitly assumes that remains constant, but it does not require the specification of any given value of .

Our use of weighted averages over is derived from the assumption that wave power is distributed isotropically in three-dimensional space. In Appendix A we discuss the motivations for assuming such an isotropic distribution of wavenumber vectors (specifically for the fast-mode waves). For Alfvén waves, this assumption has no impact on solving Equation (14), since the arguments of both angle-bracketed quantities given above are independent of (see also Hollweg, 1974). Thus, one obtains the same result for Alfvén waves whether one assumes a single value of or the isotropic distribution. For fast and slow mode MHD waves, some quantities depend strongly on and others do not. For example, slow-mode waves in low-beta plasmas have values of the angle-dependent quantity that are always nearly equal to . Figure 1 shows the radial dependence of for the isotropic distribution of fast-mode waves.

We solve Equation (14) for the energy densities of the three MHD modes (, , ), and we compute the dispersion and energy partition properties of all three wave types as given in Section 2.2. At this stage, we neglect couplings between multiple modes and other nonlinear effects. This is an approximation that is likely to break down wherever the wave amplitudes become large (e.g., Chin & Wentzel, 1972; Wentzel, 1974; Goldstein, 1978; Lacombe & Mangeney, 1980; Poedts et al., 1998; Vasquez & Hollweg, 1999; Del Zanna et al., 2001; Gogoberidze et al., 2007). In Section 4 we discuss the likelihood of rapid coupling between the high-wavenumber tails of the Alfvén and fast-mode power spectra. However, we continue to assume that the total energy densities are given by the solution of the individual transport equations.

To specify the dissipation rates , we include both linear collisional effects (e.g., viscosity, thermal conductivity, and electrical resistivity) for all three modes and nonlinear turbulent damping for the Alfvén and fast mode. Thus, we use

 QA=~QA+2γAUA,QF=~QF+2γFUF,QS=2γSUS. (17)

We give the amplitude damping rates , which include an approximation for the transition from strongly collisional to collisionless regimes, in Appendix B. The turbulent damping rates and are described in more detail below. In general, these rates depend on the parallel and perpendicular components of the wavenumber . For the purposes of evaluating these rates in the global wave transport equations, we assumed that for all three modes, where is the turbulent correlation length described below. For the fast and slow modes, our assumption of an isotropic distribution of wavenumbers is consistent with also assuming . For the Alfvén mode, we found that never depended on the assumed value of at all, but for completeness we used the critical balance condition (introduced in Appendix A) to specify .

We adopt phenomenological forms for the turbulent dissipation rates that are equivalent to the total energy fluxes that cascade from large to small scales. Thus, and are constrained only by the properties of fluctuations at the largest scales, and they do not specify the exact kinetic means of dissipation once the energy reaches the smallest scales (but see, however, Section 5). Dimensionally, these are similar to the rate of cascading energy flux derived by von Kármán & Howarth (1938) for isotropic hydrodynamic turbulence. For the nonlinear dissipation of Alfvénic fluctuations, we use

 ~QA=ρ0~αAEturbZ2−Z++Z2+Z−4λ⊥ (18)

(see also Hossain et al., 1995; Zhou & Matthaeus, 1990a; Matthaeus et al., 1999; Dmitruk et al., 2001, 2002; Breech et al., 2008). For the fast-mode waves, we use

 ~QF=ρ0~αF(v2x+v2z)2VAλ⊥ (19)

where the quantity collects together the total kinetic energy in fast-mode velocity fluctuations (Chandran, 2005; Suzuki et al., 2007). Many of the terms introduced in Equations (18)–(19) are defined throughout the remainder of this subsection.

Equation (18) depends on the magnitudes of the Elsasser (1950) variables, , which specify the power in outward () and inward () propagating Alfvénic fluctuations. Alfvénic turbulent heating occurs only when there is energy in both modes. In practice we compute an effective reflection coefficient whose magnitude is always less than unity, and thus we express the Elsasser variables in terms of the Alfvénic energy density as

 Z−=√4UAρ0(1+R2),Z+=RZ−. (20)

An accurate solution for requires the integration of non-WKB equations of Alfvén wave reflection (e.g., Heinemann & Olbert, 1980; Verdini & Velli, 2007). However, our assumption that the total power varies in accord with straightforward wave action conservation has been shown to be reasonable, even in environments where is not small such as the chromosphere (van Ballegooijen et al., 2011) and interplanetary space (Zank et al., 1996; Cranmer & van Ballegooijen, 2005).

We estimate the reflection coefficient using a modification of the low-frequency approximation of Chandran & Hollweg (2009). Specifically, we examine the magnitudes of terms in the transport equation for the inward Elsasser variable,

 ∂Z+∂t+(u0−VA)∂Z+∂r=(u0+VA)(Z+4HD+Z−2HA)−Z+Z−2λ⊥ (21)

where

 HA=VA∂VA/∂r,HD=ρ0∂ρ0/∂r. (22)

Chandran & Hollweg (2009) neglected both terms on the left-hand side of Equation (21) as well as the term containing , and thus were able to solve for straightforwardly. However, in cases of strong reflection, the term containing may have a magnitude comparable to the other dominant terms. Thus, we keep all three terms on the right-hand side and solve for

 R≈2h/|HA|1+(h/|HD|) (23)

where

 h=λ⊥(u0+VA)2Z−. (24)

Equation (22) of Chandran & Hollweg (2009) is recovered in the limit of , with . In the case of purely linear reflection, Cranmer (2010) found that the most accurate local estimates for were obtained when was replaced with the positive-definite quantity

 ~HA=VAtref=(r+R⊙)(1−R⊙r). (25)

We used instead of in Equations (23)–(24) to compute .

The definitions of the turbulent dissipation rates contain the perpendicular length scale , which is an effective transverse correlation length of the turbulence for the largest “outer scale” eddies. For simplicity we use the same correlation length for both the Alfvénic and fast-mode fluctuations, but this may not be universally valid (e.g., Suzuki et al., 2007). In previous papers we assumed that scales with the transverse width of the magnetic flux tube; i.e., that (Hollweg, 1986). Here we describe the evolution of the transverse correlation length with the following transport equation,

 ∂λ⊥∂r=λ⊥2A0∂A0∂r+~βAu0+VA(Z2−Z++Z2+Z−Z2−+Z2+), (26)

where is a dimensionless constant that is often assumed to be equal to (e.g., Hossain et al., 1995). The first term on the right-hand side of Equation (26) drives the correlation length to expand linearly with the perpendicular flux-tube cross section (Hollweg, 1986). The second term takes account of the nonlinear coupling between the fluctuations and the background plasma properties. It is given in a form suggested initially by Matthaeus et al. (1994) and later generalized to nonzero cross-helicity turbulence by Breech et al. (2008) and others. Our transport equation attempts to bridge together the effects of the two terms. In the lower solar atmosphere (between the photosphere and the chosen lower boundary of 0.01 for the wave transport models) we assumed that the first term in Equation (26) is dominant, and thus .

The turbulent dissipation rates also depend on dimensionless Kolmogorov-type constants and that are often assumed to have values of order unity. For example, Hossain et al. (1995) and Breech et al. (2009) found that gives rise to dissipation rates that agree well with both numerical simulations and heliospheric observations. In our case, we used this value as a starting point, but we also varied as a free parameter in order to produce the best match to the well-constrained Alfvénic fluctuations. On the other hand, the properties of heliospheric fast-mode turbulence are not known nearly as well as the Alfvén-wave turbulence. We thus relied on the independent wave-kinetic simulations of Pongkitiwanichakul & Chandran (2012) to fix at a value of 2.3.

The Alfvénic cascade rate contains an efficiency factor that attempts to account for regions where the turbulent cascade may not have time to develop before the fluctuations are carried away by the wind. Cranmer et al. (2007) estimated this efficiency factor to scale as

 Eturb=11+(teddy/tref), (27)

where the two timescales above are , a nonlinear eddy cascade time, and , a timescale for large-scale Alfvén wave reflection (see also Dmitruk & Matthaeus, 2003; Oughton et al., 2006). The reflection time is often defined as , but we solved Equation (25) for in order to remain consistent with the adopted model for . The eddy cascade time is given by

 teddy=λ⊥√3π(1+MA)vy, (28)

where the Alfvén Mach number and the numerical factor of comes from the normalization of an assumed shape of the turbulence spectrum (see Appendix C of Cranmer & van Ballegooijen, 2005). The two limiting cases of and are roughly equivalent to the “weak” and “strong” cascade phenomenologies discussed in Section 3, but they are not precisely the same.

### 2.4. Representative Solutions

We solved the transport equations given in Section 2.3 by numerically integrating upwards from a specified set of lower boundary conditions at and assuming time-steady conditions (i.e., ). We used a logarithmic grid of 500 radial zones in that expands out to a maximum distance of  AU. The transport equations were solved with straightforward first-order Euler steps. The values of the Elsasser variables in each zone were determined by iteration, since Equations (20) and (23) do not give a simple closed-form solution for and by themselves.

There are a number of free parameters in this model whose values were not easily obtained from either theoretical calculations or observations. In addition to the lower boundary conditions on the wave energy densities , , and , there is also the lower boundary condition on the correlation length and the values of the two von Kármán constants and . Initially, we varied these six parameters randomly in order to build up a large Monte Carlo ensemble of trial solutions. For each model, we synthesized the radial variation of observable plasma fluctuations such as the root mean squared (rms) parallel and perpendicular fluctuation speeds,

 v∥=vz,v⊥=(v2x+v2y)1/2, (29)

the Elsasser variables , and the rms fractional density fluctuation amplitude . The velocity amplitudes and contain contributions from all three MHD wave types. In nearly all models produced here, is dominated by the fast mode and is dominated by the Alfvén mode. Observations of these quantities are discussed below.

There was no single set of parameter values that gave rise to perfect agreement between all of the synthesized and observed fluctuation quantities. This is not surprising, since the models are certainly incomplete and there are significant uncertainties in the observations and their interpretation. Also, even though we aimed to restrict ourselves to measurements made in “quiet” high-latitude fast wind streams, sometimes only low-latitude data were available. Thus, in Table 1 we give a set of optimized parameters that were chosen because they produce adequate agreement with the full set of observed quantities. There were other combinations of the six parameters that gave better agreement on any single observation, but in most of these cases the agreement became worse for other observations. Although the ratio of the two von Kármán constants was allowed to vary freely, the optimal value was nonetheless found to be close to the commonly used value of 2 (Hossain et al., 1995). The best photospheric value of km is intermediate between the values of 75 km (Cranmer et al., 2007) and 300 km (Cranmer & van Ballegooijen, 2005) found from earlier models.

Figure 3 shows the comparison between synthesized and observed fluctuation quantities for the model parameters given in Table 1. The observational constraints on at are a combination of the off-limb nonthermal emission line widths given by Banerjee et al. (1998) and Landi & Cranmer (2009). The observations shown between 0.3 and 1 are from Esser et al. (1999). At larger heights, becomes approximately equal to , so we truncate the curve in favor of showing the radial dependence of and more clearly. The latter are compared directly with high-speed wind data from Helios and Ulysses (Bavassano et al., 2000). Observations of longitudinal velocity fluctuations are more difficult to find, and we show only the on-disk nonthermal line width velocities of Chae et al. (1998) as a way to compare with the modeled values of .

Figure 3(a) shows how the modeled density fluctuation amplitude is dominated by slow-mode waves in the low corona () and by fast-mode waves in the extended corona and solar wind (). The low-corona observations are drawn as an approximate boundary region around the polar plume data given by Ofman et al. (1999). The intermediate data point at is an empirical value of estimated from radio sounding data (Coles & Harmon, 1989; Spangler, 2002; Harmon & Coles, 2005; Chandran et al., 2009), but it is still unclear what fraction of the measured density fluctuations are due to anything even close to ideal MHD waves. At larger distances, we show approximate ranges of density fluctuations as reported by Marsch & Tu (1990) (blue rectangles at ), Tu & Marsch (1994) (open rectangle), and Issautier et al. (1998) (blue rectangle at ).

Figure 3(b) compares the result of solving Equation (26) for with the simpler approximation of . The plot also shows fast-wind estimates of between 1.4 and 5 AU from Ulysses (Breech et al., 2008). Figure 3(b) also compares the total heating rate with observational constraints and with the modeled coronal heating rate from Cranmer et al. (2007). The shaded area between 0.2 and 5 is an envelope surrounding a collection of empirical and theoretical heating curves from Wang (1994), Hansteen & Leer (1995), and Allen et al. (1998). These rates illustrate what is needed to produce the observed coronal heating and solar wind acceleration. The area shown at larger distances () is a representation of the range of total (proton and electron) empirical heating rates estimated by Cranmer et al. (2009). Note that the turbulent heating rates and dominate the total heating rate, with approximately 70% of the total coming from and 20% from . Less than 10% of comes from the linear damping terms.

There are additional measurement techniques that may be used to further constrain the model parameters, and in future work we will incorporate as many of these as possible. For example, Hollweg et al. (2010) argued that radio measurements of Faraday rotation fluctuations may put unique empirical constraints on the value of in the corona. Also, Sahraoui et al. (2010) used multi-spacecraft data to tease out new details of the wavenumber anisotropy of MHD fluctuations, which may lead to better limits on, e.g., in the heliosphere. Unfortunately, the vast majority of these measurements have been made for the slow solar wind and not the much less structured fast wind associated with polar coronal holes. Nearer to the Sun, Kitagawa et al. (2010) used the dispersive and energy partition properties of thin-tube MHD waves to diagnose the presence and strengths of various modes in active regions. These techniques may be useful in open-field regions as well.

Although we did not include any explicit multi-mode coupling in the transport equations of Section 2.3, there is some feedback between the modes. For example, the correlation length is used in both the Alfvénic and fast-mode turbulent heating expressions, and it is also used to set the wavenumbers and in the linear dissipation rates . Thus, the choice of the lower boundary condition on can have a significant impact on the radial evolution of all three wave types. Figure 4 illustrates this by varying the photospheric value of between 30 and 300 km and using the other standard parameters from Table 1. The integrated energy densities are plotted in velocity units as . The power in the Alfvén waves changes by only a small amount because the damping is never a strong contributor to the transport equation. However, damping is a major effect for the fast and slow modes, and thus small changes in the damping rate’s normalization can have large relative impacts on the resulting energy densities.

Figure 4 shows that, no matter the choice of normalization for , it seems unlikely for the slow-mode waves to have significantly large amplitudes anywhere but in the lowest few tenths of a solar radius. This appears to be consistent with models of slow-mode shock formation and dissipation in polar plumes (Cuntz & Suess, 2001). Therefore, in the remainder of this paper, our models of turbulence in the fast solar wind ignore the slow-mode waves altogether. We also note that Figure 4 suggests that the actual fast-mode wave properties in the high-speed solar wind may be more highly variable than the Alfvén wave properties. Our use of a “standard” model for the fast-mode waves (using the parameters given in Table 1) is thus presented as an example case and not a definitive prediction.

In this section we begin constructing a model of the wavenumber distribution of Alfvén and fast-mode fluctuation power at each radial distance. We make use of a general assumption of “scale separation;” i.e., we presume that the turbulence becomes fully developed on timescales short compared to the bulk solar wind outflow and the large-scale expansion of open flux tubes. This allows us to model the turbulence as spatially homogeneous in a small volume element with constant background plasma properties. This seems to be the general assumption made by the majority of MHD simulations of turbulence in the solar wind.222See, however, “expanding box” type simulations (Grappin & Velli, 1996; Liewer et al., 2001) that attempt to include some aspects of the large-scale radial evolution of the plasma parcel undergoing a turbulent cascade, and collisionless kinetic models that include expansion effects together with local diffusion in velocity space (Isenberg & Vasquez, 2009, 2011). Whether or not this approximation is valid, it is useful to begin studying the wavenumber dependence of the cascade in this manner.

We model the MHD fluctuations as time-steady Fourier distributions of wave power in three-dimensional wavenumber space. Although additional information about the physics of turbulence can be found in more complex statistical measures of the system (e.g., higher-order structure functions), we limit ourselves to describing the power spectrum because that is the basic quantity needed to compute the quasilinear particle heating rates.

Because of the simplified flux-tube geometry discussed in Section 2.1, we assume the background magnetic field is parallel to the bulk flow velocity, and thus the system has only one preferred spatial direction (see, however, Narita et al., 2010). The random turbulent motions create a statistical equivalence between the and directions transverse to the background field, so that we can describe the power spectra as two-dimensional functions of and only. By convention, we define the full three-dimensional power spectrum in effective velocity-squared units; i.e., when integrated over the full volume of wavenumber space, the spectrum gives the fluctuation energy density per unit mass, or

 Umρ0=∫d3kEm(k). (30)

In Appendix A we review some of the basic physical processes that determine the shape of the spectrum for Alfvénic () and fast-mode () fluctuations.

We describe the driven turbulent cascade as a combination of advection and diffusion in wavenumber space. At first, it may appear that a smooth and continuous description of the spectral “spreading” of a cascade ignores too much of the inherently stochastic and nonlocal nature of turbulence. However, Chandrasekhar (1943) showed that such a model can be made to capture the essential statistics of a large ensemble of random-walk-like (i.e., Brownian) processes. Specific models of turbulent wavenumber transport using diffusion or advection equations include those of Pao (1965), Leith (1967), Tu et al. (1984), Tu (1988), Zhou & Matthaeus (1990b), Miller et al. (1996), Stawicki et al. (2001), Chandran (2008b), Matthaeus et al. (2009), Jiang et al. (2009), and Galtier & Buchlin (2010). For the cascade of Alfvénic fluctuations, we generally follow the approach taken by Cranmer & van Ballegooijen (2003). The general forms of these equations are given as

 ∂EA∂t=1k⊥∂∂k⊥{DA⊥[α⊥k⊥∂∂k⊥(k2⊥EA)−μ⊥EA]}
 +α∥∂∂k∥(DA∥∂EA∂k∥)+SA−2γAEA+CAF (31)
 ∂EF∂t=αFk2∂∂k(k2DF∂EF∂k)+SF−2γFEF−CAF (32)

and the terms on the right-hand sides of Equations (31)–(32) are defined throughout the remainder of this subsection. The mode coupling term is described further in Section 4, and the dissipation rates and are described in Section 5.

The perpendicular Alfvénic cascade is described by the first term on the right-hand side of Equation (31), and we assume an arbitrary linear combination of advection and diffusion. Cranmer & van Ballegooijen (2003) found that many key properties of the turbulence do not depend on whether the cascade is modeled as advection, diffusion, or both, so we retain all terms for maximum generality. For both the parallel Alfvénic spectral transport and the isotropic fast-mode transport, a more standard diffusion coefficient is assumed. The dimensionless multipliers to the diffusion coefficients are denoted and , to correspond roughly to in Equation (18), and the dimensionless multiplier for the wavenumber advection coefficient is denoted .

For the Alfvénic cascade, the overall behavior of wavenumber transport in the perpendicular and parallel directions is specified by the diffusion-like coefficients

 DA⊥=k2⊥τA,DA∥=(v⊥VA)2DA⊥ (33)

where is the cascade timescale defined below, and is the -dependent velocity response of the waves. Note that is independent of , so it can be pulled out of the derivative in Equation (31). Cranmer & van Ballegooijen (2003) showed that the above form for the diffusion coefficients tends to reproduce the Goldreich & Sridhar (1995) critical balance, and Matthaeus et al. (2009) derived similar functional forms for the coefficients. When specifying the properties of the wavenumber cascade, we apply the scalings for “balanced” turbulence (i.e., zero cross helicity, or ), which is more straightforward to implement but is formally inconsistent with the large-scale transport model of Section 2.

For ideal MHD Alfvénic fluctuations, is equal to , the latter representing the transverse magnetic variance spectrum divided by to convert it to units of velocity squared. Following the usual convention, the power spectrum tracks the magnetic fluctuations, so the reduced spectra are defined formally as

 b2⊥=k2⊥∫dk∥EA,v2⊥=ϕb2⊥. (34)

The dimensionless factor describes the departure from ideal MHD energy equipartition. For small values of , we assume . However, as increases into the regime of kinetic Alfvén waves (KAWs), can become much larger than 1. Hollweg (1999) described how the main difference between and in the KAW regime comes from an enhanced response of the electron velocity distribution to the electric and magnetic fluctuations. For simplicity, we use an approximate analytic expression

 ϕ=ω2k2∥V2A≈1+k2⊥ρ2p1+k2⊥ρ2pme/(βmp), (35)

where is the proton thermal gyroradius, with the proton most-probable speed given by and the proton cyclotron frequency by . Our term is equivalent to as defined by Howes et al. (2008).

Inspired by Equation (A6), we define the Alfvénic spectral transport timescale as

 τA=1+χ0k⊥v⊥ (36)

where we chose to replace the general critical balance parameter by its value at the outer-scale parallel wavenumber . Thus,

 χ0=ω0k⊥v⊥≈k0∥VAk⊥b⊥, (37)

and is the appropriate critical balance parameter to use when solving for the properties of the dominant low-frequency cascade. From Equation (35) we see that KAW outer-scale frequency , so that a factor of cancels out of both the numerator and denominator to give the final approximate expression above. The wavenumber specifies the spatial scale along the field at which energy is injected in the source term (see below). Because is assumed to be constant (at a given heliocentric distance ), the parameters and are both functions of and not . The above form for Equation (36) was motivated by the analysis of Zhou & Matthaeus (1990b), Chandran (2008b), and Howes et al. (2012), who described how the cascade and wavenumber anisotropy change when the system transitions from weak () to strong () turbulence.

As mentioned above, our expressions for , , and assume zero cross helicity (i.e., ). There is still no agreement about how to generalize these terms when inefficient wave reflection gives rise to nonzero cross helicity. Lithwick et al. (2007) found that the cascade timescales for outward and inward wave modes are different from one another when , but their parallel spatial scales are the same. However, Beresnyak & Lazarian (2008, 2009) found that for the outward mode should be larger than for the inward mode, and thus the Goldreich & Sridhar (1995) critical balance must be modified (see Equation (45) below). Chandran (2008b) outlined a method for setting up the advection-diffusion equations in the case of , but we defer a full implementation of that approach to future work.

Putting aside the issue of imbalanced turbulence, the dominant perpendicular nature of the Alfvénic cascade allows us to define a reduced transport equation that follows the evolution of the spectrum as a function of only. If we ignore the mode coupling term for now, we can multiply Equation (31) by and integrate over to obtain

 ∂b2⊥∂t=k⊥∂∂k⊥[1τA(α⊥k⊥∂b2⊥∂k⊥−μ⊥b2⊥)]+~SA−2~γAb2⊥. (38)

This is essentially the same as Equation (11) of Cranmer & van Ballegooijen (2003). The reduced source term and dissipation rate are defined similarly to the corresponding terms in Equation (31), but they are weighted toward the low- regions of wavenumber space that are “filled” by the cascade. In Appendices C.1C.3 we derive analytic solutions for the time-steady Alfvén-wave power spectrum in various limiting cases.

The cascade of fast-mode waves, described by Equation (32), appears to be conceptually simpler than the strongly anisotropic Alfvén-wave cascade. The diffusion coefficient is given by , where is related to the IK-like cascade time given by Equation (A2) with . There is increasing evidence (e.g., Markovskii et al., 2010) that a fast-mode cascade is more rapid in the directions perpendicular to the field than along the field. However, the cascade does appear to proceed outward “radially” in the direction of increasing . Thus, it makes the most sense to use an isotropic diffusion formalism as in Equation (32), but scale the magnitude of the diffusion timescale with . Following the weak turbulence model of Chandran (2005), we adopt

 τF=VAkv2ksinθ, (39)

which implies that

 DF=k3v2ksinθVA=4πk6EFsinθVA. (40)

Chandran (2005) showed that the dependence in the denominator of is consistent with an isotropic energy flux for the cascade, but it does not guarantee an isotropic wavenumber spectrum . More information about how we chose to implement the fast-mode cascade is given in Appendices C.4 and C.5.

In order to fully describe the cascade in the advection-diffusion equations, four dimensionless spectral transport constants (, , , ) need to be specified. Matthaeus et al. (2009) summarized the results of many MHD turbulence models and found that often takes on values between 0.2 and 0.5, and seems to be a useful parameterization (see Equation 13 of Matthaeus et al., 2009). Zhou & Matthaeus (1990b) and Matthaeus et al. (2009) made a case for a classical form of the diffusion operator that implies . Alternately, van Ballegooijen (1986) found that a cascade of random-walk-like displacements of magnetic flux tubes is described well by . Howes et al. (2008) and Chandran (2008b) used a straightforward advection equation to model an Alfvénic cascade, which sets and assumes . For this type of model, Howes et al. (2008) derived .

In our models, we are constrained by the values of the cascade constants and used in the global transport equations of Section 2. We related these constants to the ones defined above by integrating the cascade advection-diffusion terms over wavenumber to find . By demanding this quantity be equal to the heating rate , we obtained

 2α⊥3+μ⊥=3√6π4~αA, (41)

which assumes that the perpendicular cascade is dominant and that , and

 αF=327π~αF. (42)

We keep the ratio as a free parameter and we explore the ramifications of varying it below. Note, however, that if we used (as assumed by Zhou & Matthaeus, 1990b), then Equation (41) gives and . These are roughly consistent with the constants given by Zhou & Matthaeus (1990b) and Matthaeus et al. (2009). To complete the system of cascade constants, we adopt the Matthaeus et al. (2009) choice for , but we compute this quantity using the Matthaeus et al. (2009) assumption of .

The source terms, in Equation (31) and in Equation (32), describe the outer-scale injection of fluctuation energy. The global energy balance of the waves is already described by the radial transport model of Section 2. Thus, we specify the magnitudes of and by demanding that the time-steady total energy densities and be maintained at their known values at a given distance . From a physical standpoint, however, it is unclear whether the passive propagation of waves dominates the source terms, or whether there is significant local “stirring” that converts large-scale dynamical motions (e.g., velocity shears in evolved corotating streams) into new fluctuations.

We adopt specific functional forms for and that are described in detail in Appendix C. Generally, the source terms are nonzero only at the lowest wavenumbers, at which the fluctuations are driven. For the Alfvén waves, we continue to use the assumption from Section 2.3 that the perpendicular driving scale is set by the turbulence correlation length; i.e., . For the fast-mode fluctuations, we assume their outer-scale wavenumber magnitude is also equal to , since the largest-scale transverse stirring motions are likely to be common to both Alfvénic and fast-mode waves. There are several ways that one could imagine defining the parallel outer-scale Alfvén wavenumber :

1. Monochromatic Alfvén waves that propagate up from the corona retain a constant frequency in the Sun’s inertial frame. However, because the phase speed varies with distance, the corresponding wavelength undergoes “stretching” commensurate with the dispersion relation

 k0∥=ω0u0+VA. (43)
2. The fluctuations propagating up from the Sun may already be fully turbulent (see, e.g., van Ballegooijen et al., 2011). Thus, the outer-scale parallel wavenumber may be coupled continuously to the perpendicular wavenumber via critical balance (Goldreich & Sridhar, 1995), with

 k0∥≈k0⊥VA√UAρ0. (44)
3. In flux tubes with nonzero cross helicity (i.e., ), Beresnyak & Lazarian (2008, 2009) found that the inward waves should obey the Goldreich & Sridhar (1995) critical balance, but the outward waves (which are generally what we intend to model) obey a modified version of critical balance, which we approximate as

 k0∥≈k0⊥VA√UAρ01R. (45)
4. In some cases we assume that the dimensionless ratio remains fixed at a constant specified value. Many studies of MHD turbulence assume isotropic forcing at the outer scale, which is consistent with the fixed ratio . The lack of a physical justification for this approximation is offset by its simplicity.

Figure 5 illustrates the ratio for several of the above methods of setting the parallel outer scale. For example, it shows the result of evaluating Equation (43) for a range of wave periods between 1 and 100 minutes. Constant assumed values of would correspond to horizontal lines in Figure 5.

### 3.2. Solutions in the Absence of Coupling

Here we present some example results for the power spectra and . These spectra are computed from Equations (31), (32), and (38) in the limiting cases of time independence and no mode coupling (). The Alfvénic spectrum was first computed in its reduced form using the solutions for given in Appendices C.1 and C.2, and then it was expanded into full wavenumber space by using the results of Appendix C.3. The shape of the fast-mode spectrum was determined from the analytic solutions given in Appendices C.4 and C.5.

To illustrate the wavenumber dependence of the power spectra, we chose a single coronal height at which . We typically plot the wavenumbers in terms of dimensionless quantities and . Dissipative wave-particle interactions tend to become important when these quantities reach order-unity values, and ideal MHD conditions apply when these quantities are small. Typically, the driving scale for Alfvénic turbulence occurs at to , with the larger values generally occurring at larger heliocentric distances.

In Figure 6 we show the time-steady dependence for the Alfvénic and fluctuations, both with and without KAW dissipation. To set the cascade properties, we utilized the values of the constants given in Section 3.1, and we also assumed and . The KAW damping ratio appropriate for the assumed value of , which was used in Equation (C11), is also shown in green (see also Section 5). At the outer scale, the peak value of is 42 km s. We caution that this value should not be assumed to be equivalent to the full rms velocity amplitude. In this case, km s, which is almost a factor of 5 larger than the maximum value of at this height.

The damped spectra shown in Figure 6 have several features that resemble those of measured KAWs in the solar wind. Using the conventional form of the reduced energy spectrum () we found that the magnetic fluctuation power made a transition from a Kolmogorov-like power law to a steeper spectrum with at . The spectrum becomes shallower again around because the wavenumber dependence of flattens out at low values of . This behavior is reminiscent of that predicted by Voitenko & De Keyser (2011). At larger radial distances where , the KAW dispersion relation (Equation (35)) gives rise to a more sustained increase in with increasing . This in turn produces spectra that remain steep, with persisting over several orders of magnitude of in agreement with both measurements (Smith et al., 2006; Sahraoui et al., 2010) and other models (Howes et al., 2008). We note that the predicted undamped KAW power-law decline of (see Appendix C.1) was not seen for any sustained range of .

Figure 7 shows the result of varying the normalization of the parallel outer scale wavenumber on the shape of . We kept the same value of that was used in Figure 6, but we varied the constant ratio over five orders of magnitude. For the lowest values of the outer-scale critical balance ratio always remains much smaller than unity. This means that the stirring or forcing takes place well within the “filled” region of wavenumber space, and thus strong turbulence occurs. In this case, and thus . The opposite extreme case of large corresponds to and weak turbulence with less anisotropic driving. In that limit, the inertial range spectra are given by and . Our model shows the gradual transition between these two extreme cases.

In Figure 8 we compare the Alfvén and fast-mode spectra with one another. As above, we used the background conditions at a coronal height of and we assumed . We illustrate the most extreme case of a lack of high-frequency Alfvénic power by showing the contours of for the case . In this limit, Equation (C16) describes an exponential decrease of power with increasing . Other comparable examples of this kind of spectrum can be found in Figure 4b of Cranmer & van Ballegooijen (2003) and Figures 1 and 2 of Jiang et al. (2009). We computed the Alfvénic and fast-mode spectra with the kinetic sources of damping that were described in Section 5. Note that experiences the strongest damping at intermediate values of . For or , the transit-time damping described by Equation (58) is relatively weak.

## 4. Coupling Between Alfvén and Fast-Mode Waves

### 4.1. Basic Physics and Phenomenological Rates

There are several ways that the ideal linear MHD wave modes can become coupled to one another in the corona and solar wind:

1. Inhomogeneities in the background plasma can blur the definitions of the individual modes. For example, linear reflection due to radial variations in (Ferraro & Plumpton, 1958; Heinemann & Olbert, 1980) may produce not only incoming Alfvén waves (i.e., ), but also fast and slow magnetosonic waves (e.g., Stein, 1971; McDougall & Hood, 2007). In addition, large-scale bends in the background magnetic field (Frisch, 1964; Wentzel, 1974), density variations between flux tubes (Valley, 1974; Markovskii, 2001; Mecheri & Marsch, 2008), or velocity shears (Poedts et al., 1998; Gogoberidze et al., 2007) can drive instabilities that partially convert Alfvén waves into other modes.

2. Even in a homogeneous medium, the MHD waves begin to lose their ideal linear character when their amplitudes become large. Nonlinear Alfvén waves naturally drive second order fluctuations in and that mimic the properties of both slow and fast magnetosonic waves (Hollweg, 1971; Spangler, 1989; Vasquez & Hollweg, 1999). Large-amplitude waves also excite a range of wave-wave interactions that can often be characterized either as two modes giving birth to a third, or one mode splitting into several others (e.g., Chin & Wentzel, 1972; Goldstein, 1978; Del Zanna et al., 2001; Sharma & Kumar, 2010). Models of weak turbulence, in which the wave-wave interactions describe the cascade process (Chandran, 2005, 2008a; Luo & Melrose, 2006; Yoon & Fang, 2008) also create this kind of coupling.

3. Although not strictly a multi-mode coupling, when the Alfvén mode begins to exhibit oscillations in density, parallel electron velocity, and the parallel electric and magnetic fields (Hasegawa & Chen, 1976; Hollweg, 1999). Observationally, it has proved difficult to separate such dispersive KAW density fluctuations from those arising from independent sources of fast or slow MHD waves (e.g., Harmon & Coles, 2005; Chandran et al., 2009).

In this paper we take account of one particular nonlinear effect from the second entry in the above list. Specifically, Chandran (2005) suggested that weak turbulence couplings between Alfvén and fast-mode fluctuations may provide enough power at high to induce substantial ion cyclotron heating. Suzuki et al. (2007) argued that this effect may be relatively unimportant because the fast-mode cascade timescale is long in comparison to the Alfvén cascade timescale . This may be the case in the low-frequency regime of wavenumber space where , but at the cyclotron resonant frequencies of interest () the Alfvénic cascade is quenched because . The fast-mode cascade may in fact even be faster than any intrinsic Alfvénic spectral transfer in this region of wavenumber space. Therefore, we proceed using the Chandran (2005) results for Alfvén/fast-mode coupling.

We express the coupling term in Equations (31)–(32) as

 CAF=EF−EAτAF (46)

such that, in the absence of other processes, the power spectra at a given wavenumber are driven toward a common value over a coupling timescale . The weak turbulence model of Chandran (2005) gave an approximate value for this timescale of

 τAF≈1523π2τFsin2θ (47)

which holds in the limiting cases of and nearly parallel propagation (). In the opposite case of , it’s likely that would no longer depend linearly on , and may scale instead with . However, the region of wavenumber space with which we are most concerned is the high-, low- ion cyclotron regime. At those wavenumbers, we know that in the absence of coupling the condition is likely to be satisfied, and the coupling will be a transfer of energy from the dominant fast-mode spectrum to the much less intense Alfvén mode.

The wave-wave conditions of frequency and wavenumber matching (e.g., Sagdeev & Galeev, 1969) confirm that the most rapid coupling should occur when the dispersive properties of the Alfvén and fast-mode waves are the most similar to one another; i.e., at . Note that Equation (39) gave , so the combined dependence for the coupling time is . In practice, however, we found that using this ideal expression for could lead to an unphysical singularity at . We removed this singularity by replacing in Equation (47) by . We set to a constant value of 0.01 to avoid having an infinitely fast coupling rate at parallel propagation.333Also note that the magnetic field in MHD turbulence undergoes a complex, multi-scale “wandering,” such that the direction corresponding to is continuously varying in time and space (see, e.g., Ragot, 2006; Shalchi & Kourakis, 2007). Thus, the plasma may seldom “see” exactly parallel wavenumber conditions. To retain the most generality, we chose to reparameterize the coupling timescale as

 τAF=1ΦτFsin2(θ+δθ) (48)

where we find it useful to vary the constant coupling strength up or down from the value of derived by Chandran (2005). The case corresponds to ignoring the coupling altogether.

Note that the above form for the coupling timescale implies that , so that the coupling is rapid at wavenumbers corresponding to ion cyclotron resonance (large , small ). The coupling is much slower at KAW wavenumbers favored by the pure Alfvénic cascade (small , large ). Thus, the bulk of the Alfvénic spectrum at is likely to be more or less unaffected by the coupling. This seems to be consistent with our assumption that the integrated energy densities and also remain uncoupled from one another. We realize that this may be a severe underestimate of the degree of energy transfer between Alfvén and magnetosonic modes in the corona and solar wind. However, one main purpose of this paper is to investigate how much can be accomplished with only this small degree of coupling in the high- tails of the power spectra.

### 4.2. Approximate Solutions for Coupled Spectra

The exact solutions to Equations (31) and (32) with coupling () must be found numerically. Here we present an approximate solution that is both (1) likely to reflect the proper behavior of more rigorous numerical solutions in many limiting regimes of parameter space, and (2) efficient to implement on a large grid of model spectra spanning a wide range of heliocentric distances. We begin by approaching the problem iteratively; i.e., we solve Equation (31) for under the assumption that is known, and we then solve Equation (32) for under the assumption that is known. The analytic solutions derived below suggest a natural way to terminate this iteration after only one round.

When solving the advection-diffusion equation for Alfvénic fluctuations, let us temporarily ignore the outer-scale source term and the dissipation-range damping term that depends on . Since we are most concerned with the generation and transport of wave power in the high- regions that undergo ion cyclotron resonance, we consider the weak turbulence regime of , in which the transport of energy is mainly from low to high and there is negligible parallel spreading (see also Oughton et al., 2006). Thus, we solve the advection-diffusion equation for discrete, non-interacting “strips” of wavenumber space each having constant . The nonlinear coupling supplies wave energy locally, and the Alfvénic cascade takes it from low to high . If we simplify further by assuming pure advection (i.e., ), the time-steady version of Equation (31) becomes

 μ⊥k⊥∂∂k