Electron-ion heating partition in imbalanced solar-wind turbulence

A likely candidate mechanism to heat the solar corona and solar wind is low-frequency"Alfv\'enic"turbulence sourced by magnetic fluctuations near the solar surface. Depending on its properties, such turbulence can heat different species via different mechanisms, and the comparison of theoretical predictions to observed temperatures, wind speeds, anisotropies, and their variation with heliocentric radius provides a sensitive test of this physics. Here we explore the importance of normalized cross helicity, or imbalance, for controlling solar-wind heating, since it is a key parameter of magnetized turbulence and varies systematically with wind speed and radius. Based on a hybrid-kinetic simulation in which the forcing's imbalance decreases with time -- a crude model for a plasma parcel entrained in the outflowing wind -- we demonstrate how significant changes to the turbulence and heating result from the"helicity barrier"effect. Its dissolution at low imbalance causes its characteristic features -- strong perpendicular ion heating with a steep"transition-range"drop in electromagnetic fluctuation spectra -- to disappear, driving more energy into electrons and parallel ion heat, and halting the emission of ion-scale waves. These predictions seem to agree with a diverse array of solar-wind observations, offering to explain a variety of complex correlations and features within a single theoretical framework.


INTRODUCTION
The solar corona and its extended outflow, the solar wind, provide us with an unparalleled laboratory for studying the physics of magnetized collisionless plasmas.Decades of observations have revealed a highly complex system filled with electromagnetic fluctuations across a vast range of scales, whose properties (e.g., total power or spectra) are correlated in surprising and nontrivial ways with the plasma's flow speed, temperatures, anisotropies, and elemental abundances (Marsch 2006;Horbury et al. 2012;Bruno & Carbone 2013).These correlations, as well as the extended heating at large heliocentric distances needed to explain the high speed of fast-wind streams (Parker 1965), suggest that the solar wind is shaped by both properties of its low-coronal source and turbulent heating at larger altitudes.Of particular interest are the decades of observations that hint at the role played by the fluctuations' imbalance (i.e., normalized cross helicity), a key parameter in the theory of magnetized turbulence (Dobrowolny et al. 1980;Schekochihin 2022) that is observed to correlate with wind speed U and decrease with increasing heliocenjonathan.squire@otago.ac.nz tric radius R (e.g., Roberts et al. 1987;Marsch 2006;D'Amicis et al. 2021).
In this Letter, we argue for the importance of imbalance in shaping turbulence and heating in the low-β solar wind.We focus on the physics of the "helicity barrier" (Meyrand et al. 2021), which dramatically alters imbalanced turbulence in β ≲ 1 plasmas by restricting the turbulent energy cascade at perpendicular scales λ ⊥ below the ion gyroradius ρ i .Previous kinetic simulations have shown how numerous features of helicitybarrier-mediated turbulence match measurements of the low-β solar wind well (Squire et al. 2022;hereafter S+22), including properties of the steep "transitionrange" drop in electromagnetic field spectra around ρ i scales and the ion velocity distribution function (VDF; see, e.g., Duan et al. 2021;Bowen et al. 2022Bowen et al. , 2023)).
We study how the turbulence evolves as imbalance is slowly decreased in time, transitioning between the highly imbalanced and balanced regimes.This physics is of interest for two reasons: first, it is a crude model for the expected response of the turbulent heating as the plasma flies outwards from the Sun and its imbalance decreases; second, it captures how the helicity barrier becomes weaker and eventually disappears at low imbalance, thus probing its robustness.The results, which are based on a kinetic simulation of forced Alfvénic turbu-lence, can explain the observed correlation of wind speed with ion temperatures, as well as the switch from negative to positive correlation of wind speed with electron temperature at increasing R (e.g., Burlaga & Ogilvie 1973;Marsch et al. 1989;Shi et al. 2023;Bandyopadhyay et al. 2023).Correlations with other properties such as the spectral transition range, proton VDFs, and plasma-wave/instability activity seem to explain a diverse array of observations within a single theoretical framework.

The helicity barrier
The helicity barrier arises due to the conservation of a "generalized helicity" H in low-β gyrokinetics (Schekochihin et al. 2019;Meyrand et al. 2021).At large scales λ ⊥ ≫ ρ i , H is the cross helicity, which, along with the energy, cascades forward to smaller scales; at small scales λ ⊥ ≪ ρ i , H is a form of magnetic helicity, which cannot cascade to smaller λ ⊥ without violating energy conservation (Fjørtoft 1953).Thus, the sub-ρ i cascade must have H = 0, restricting its energy flux to be at most the "balanced portion" of that injected, ε − ε H , where ε and ε H are the large-scale injection rates of energy and helicity, respectively.The remaining energy input, ε − (ε − ε H ) = ε H , is "trapped" at λ ⊥ ≲ ρ i scales, causing the turbulence to reach large amplitudes and (through critical balance) small parallel scales approaching the ion inertial length d i .The fluctuations then take on the character of oblique ion-cyclotron waves (ICWs) with frequencies ω that approach the ion gyrofrequency Ω i .This activates strong wave-particle interactions and quasi-linear resonant heating (Kennel & Engelmann 1966), which can absorb the remaining energy input, allowing the turbulence to saturate (Li et al. 1999).It also causes the plasma to become unstable and emit small-scale parallel ICWs (Chandran et al. 2010b).By diverting a large fraction of the energy flux that would have reached electron scales, and thus electron heat, into ions instead, the helicity barrier causes large ion-to-electron heating ratios

Numerical method
We use the hybrid-kinetic method implemented in the Pegasus++ code (Kunz et al. 2014;Arzamasskiy et al. 2023).This approach treats the ion (proton) dynamics fully kinetically using a particle-in-cell (PIC) approach, while the electrons constitute a massless, neutralizing, isothermal fluid.The ion macro-particle positions (r) and velocities (v) are drawn from an initially Maxwellian f i (v) and evolved via where E and B are the electric and magnetic fields and F U ⊥ stirs turbulence by injecting incompressible motions perpendicular ("⊥") to a mean ("guide") magnetic field B 0 .The magnetic field satisfies a modified version of Faraday's law, where F B ⊥ forces solenoidal magnetic fluctuations perpendicular to B 0 and η 4 is a hyper-resistivity that absorbs small-scale magnetic energy.The electric field is where n is the ion (and electron) density and u is the ion flow velocity, both of which are computed via a weighted sum of the macro-particles in the relevant region of space; T e is the electron temperature, which is a parameter of the model; and e, m i , and c are the electron/ion charge, the ion mass, and the speed of light, respectively.We also define the Alfvén speed v A ≡ B 0 / √ 4πnm i , Ω i ≡ eB 0 /(m i c), the ion thermal speed v th , and the "Elsässer" fields Volume averages are denoted by ⟨ • ⟩ and f i (w ⊥ , w ∥ ) is the gyro-averaged VDF, with w ⊥ and w ∥ the particle velocities perpendicular and parallel to the local magnetic field in the frame of the plasma (i.e., w ≡ v − u).

Problem set up
Our basic simulation set up follows S+22.The simulation represents a small co-moving patch of plasma, capturing realistic solar-wind turbulence amplitudes and anisotropies around k ⊥ ρ i ∼ 1 scales, where k ⊥ ∼ 1/λ ⊥ denotes the perpendicular wavenumber.The domain is Cartesian and periodic, with coordinates {x, y, z}, and is elongated along B 0 = −B 0 ẑ, with L z = 6L ⊥ and L ⊥ = 67.5di (d i = v A /Ω i ).With these parameters, the box shape approximately matches the (statistical) shape of turbulent eddies measured at similar scales in the solar wind (Chen et al. 2016;S+22).The grid resolution is N 2 ⊥ × N z = 392 2 × 2352, so that the smallest resolved scales are k ⊥,max d i ≃ πN ⊥ d i /L ⊥ ≈ 18.We use N ppc = 216 ion-macro-particles per cell.The hyperresistivity was increased from η 4 ≈ 2.4 × 10 −5 d 4 i Ω i to η 4 = 5 × 10 −5 d 4 i Ω i at t = 4.3τ A in response to the strengthening kinetic-range cascade.
The simulation is initialized at t = 0 using the final snapshot from S+22, and thus represents saturated highly imbalanced turbulence with strong perpendicular ion heating and a helicity barrier, similar to turbulence observed in the near-Sun solar wind by PSP (Bowen et al. 2023).The initial ion temperature is such that β i = 8π⟨nT i ⟩/⟨B 2 ⟩ ≈ 0.33; electrons are isothermal with β e0 = 8π⟨n⟩T e /B 2 0 = 0.3.The simulation is run for an additional ≈18τ A , where τ A = L z /v A is the outer-scale Alfvén time, which is also comparable to the turnover time of the turbulence τ turb ∼ L ⊥ /u rms due to our choice of forcing parameters (see below; u rms = ⟨u 2 ⟩ 1/2 is the root-mean-square velocity).Because the ions heat up, β i increases over the course of the simulation to ≃0.45 and ρ i = √ β i d i changes modestly.The domain resolves scales between k ⊥0 ρ i0 ≡ 2πρ i0 /L ⊥ ≈ 0.05 and k ⊥max ρ i0 ≈ 10, where ρ i0 = ρ i (t = 0).

Decreasing imbalance with distance from the Sun
A novel feature of this work is the forcing, which changes from imbalanced to balanced over the simulation, heuristically mimicking the radial evolution of the solar wind.We inject energy and cross helicity at the rates ε and ε H , respectively, with forcing functions F U ⊥ and F B ⊥ that are intended to capture the effect of stirring due to turbulent eddies above the box scale.These consist of random combinations of large-scale Fourier modes with wavenumbers k j satisfying 2π/L j ≤ k j ≤ 4π/L j for j = {x, y, z}.They are computed as where F 0 is divergence-free, perpendicular to B 0 , and evolved in time via an Ornstein-Uhlenbeck process with correlation time τ A /2.We fix ε and ε H at each time step by adjusting f U and f B so that nu • F U ⊥ and B • ∇ × F B ⊥ take the values needed to inject the required energies into z ± .F B ⊥ is computed from F 0 via a fast Fourier transform and used in the standard Pegasus++ constrained-transport algorithm to evolve B, ensuring that ∇ • B = 0 to machine precision. We , where C A = 0.29 is the Kolmogorov constant and V is the simulation volume; the factor (L ⊥ /L z ) 2 v 2 A /τ A guarantees critically balanced fluctuations at the outer scale with u rms ∼ (L ⊥ /L z )v A and τ turb ∼ τ A .In contrast, ε H decreases in time during the simulation so that the "injection imbalance" ε H /ε starts at 0.9 (as in S+22) then decreases linearly in time at a rate of 0.1 every 1.5τA , reaching ε H = 0 (balanced forcing) at t = 13.5τA , where it remains for the rest of the simulation (see Fig. 1).Because the cross helicity and energy are approximately conserved at k ⊥ ρ i ≪ 1 scales, this evolution is intended to mimic crudely the effect of larger scales becoming more balanced with radius, thereby driving the smaller scales with decreasing ε H /ε.However, in the actual solar wind, the timescale over which the imbalance decreases is comparable to both the turbulent decay timescale and the expansion timescale τ exp ≈ R/U (Roberts et al. 1987;Meyrand et al. 2023).By keeping ε fixed, we effectively assume that the turbulent decay and expansion are slow compared to the simulation's duration, which is appropriate since Bale et al. 2019).However, this also implies that our adopted imbalance-decrease timescale of ≈15τ A is too short.This trade-off is unavoidable given the extreme computational expense of kinetic-turbulence simulations, but it should be kept mind that the simulation can only qualitatively capture realistic features of this transition in the solar wind.

Evolution of the helicity barrier
The time evolution as the injection imbalance decreases is shown in Fig. 1.Properties of the saturated imbalanced turbulence, from which the system is initialized, are discussed in S+22; a portion of this phase is shown at t < 0 for comparison.The imbalance of the fluctuations, starts at σ c ≈ 0.98 for t < 0, then decreases as the system adjusts to the changing forcing (ε H /ε). We halt the simulation once σ c ≃ 0.4, a value similar to that in the solar wind around 1 au; as shown below, the heating properties have changed significantly by this point even though σ c is still nonzero.The blue line shows the fluctuation energy density which decreases over the simulation, as expected because the turbulence can dissipate more effectively at lower ε H /ε. The lower panel of Fig. 1 compares the measured heating rates, Q i and ε η , to helicity-barrier predictions.Here ε η , the hyper-resistive dissipation, is taken as a proxy for the electron-heating rate Q e , which assumes that ions cannot interact with k ⊥ ρ i ≫ 1 eddies and that there is minimal electron heating via ICWs or other k ⊥ ρ i ≲ 1 fields (Chandran et al. 2010b).For most of the simulation, ε η (green line) tracks well the helicitybarrier prediction Q e = ε − ε H (dotted-green line), aside from a delay of ∼τ A as the injected energy cascades towards smaller scales.Q i increases initially because the energy E ⊥ (t) decreases while the sub-ρ i energy flux is fixed by the helicity barrier.Its large value results from the unrealistically fast decrease of ε H /ε: in reality E ⊥ (t) would change more slowly, implying that ) and the helicity barrier erodes, the heating departs from Q e ≈ ε − ε H and approaches that found in kinetic simulations of balanced turbulence at β ≈ 0.3 (Kawazura et al. 2019;Cerri et al. 2021), in which ion heating via Landau damping and stochastic heating absorbs a portion of the energy flux before it reaches the smallest scales.At around the same time, Q i drops significantly.
Figure 2 shows the evolution of the perpendicular spectra E(k ⊥ ) of E ⊥ and B ⊥ .The flatter sub-ρ i range in the former helps to highlight the double-kinked "transition-range" power law; for k ⊥ ρ i ≲ 1.5, including in the transition range, E E ⊥ ≈ E B ⊥ .As the turbulence becomes balanced, the spectral break smoothly moves towards smaller scales, creating a less pronounced transition range that is narrower in k ⊥ and less steep, eventually transitioning directly into the kinetic range with E E ⊥ ∼ k −0.8 ⊥ .At late times, δB ⊥ is resistively damped by the larger resistivity and the magnetic spectrum does not exhibit a clear ∼k −2.8 ⊥ kinetic range.The lower panel of Fig. 2 provides the evolution of the break scale k * ⊥ , obtained by fitting a broken power-law to E B ⊥ (k ⊥ ). 1 We find a clear power-law dependence on the energy imbalance, , for σ c ≳ 0.8 (t ≲ 11τ A ) when the barrier is active (ε H /ε ≳ 0).Although this scaling remains unexplained theoretically, it matches the results of low-β gyrokinetic simulations for various choices of ε H /ε (Meyrand et al. 2021, figure 7), suggesting it is a robust consequence of the helicity barrier and independent of the energydissipation mechanism (gyrokinetics is ignorant of ICW kinetic physics).The persistence of the transition-range drop, as well as the evolution of ε η and Q i in Fig. 2, provide good evidence that a helicity barrier actively mediates the turbulence until t ≈ 11τ A , at which point it transitions towards the balanced regime.

Consequences for ion heating
This drop in ion heating associated with the balanced transition is diagnosed in Fig. 3.We show the perpendicular energy-diffusion coefficient D E ⊥⊥ (panel a), which is computed from the measured evolution of f i via where (Vasquez et al. 2020), averaging over the w ∥ variation of f i and assuming that the heating is predominantly perpendicular, as appropriate for t ≲ 12τ A (Fig. 3a inset).In quasilinear theory, waves strongly scatter "resonant" particles with w ∥ = w ∥res ≡ ω(k)/k ∥ − nΩ i /k ∥ for n ∈ Z, flattening f i along contours of constant energy in the wave's frame (the "resonance contours"; we consider 1 The functional form of the fit is , and n controls the break's sharpness; we fit E B ⊥ (k ⊥ ) in the range 0.14 < k ⊥ ρ i < 2. The measured proportionality between k * ⊥ ρ i and (1−σc) 1/4 varies with fitting choices, but the power-law exponent (1/4) is robust.
where there is a sharp dropoff in wave power (Fig. 4).By late times it drops significantly and an ill-defined peak around w ⊥ ≃ v th appears, reminiscent of stochastic heating.The Q ∥ profile (panels d and e) maintains its structure, although spreads out modestly (accounting partially for its lower values).
only n = 1 resonances).When parallel waves dominate the spectrum, D E ⊥⊥ ∝ w 2 ⊥ , because w ∥res and various wave-intensity/polarization factors are all independent of w ⊥ (Kennel & Engelmann 1966).Although the wave spectrum here is more complex, involving a mix of oblique and parallel ICW modes (see below), this prediction is nonetheless well satisfied until t ≈ 11τ A . 2 After this, D E ⊥⊥ drops and flattens, causing a large drop in the perpendicular ion heating Q ⊥,i , even while the parallel heating Q ∥,i remains almost constant (see inset).At later times, the form of D E ⊥⊥ does not clearly indicate a particular heating mechanism, but is plausibly consistent with stochastic heating (Chandran et al. 2010a): computing the D E ⊥⊥ predicted for stochastic heating using the method of Cerri et al. (2021) gives a similar shape (not shown).As a complementary analysis, panels (b) and (c) display ∂ 2 Q ⊥,i /∂w ⊥ ∂w ∥ , the differential perpendicular heating per unit velocity, which is computed directly from ⟨E ⊥ • w ⊥ ⟩ evaluated along particle trajectories (the so-called field-particle correlation technique; Klein & Howes 2016;Arzamasskiy et al. 2019).At early times a resonant feature is centered around the w ∥ that resonates with the smallest-(d i )-scale oblique ICWs with significant power (the shaded region shows w ∥res for d −1 i ≲ k ∥ ≲ 2d −1 i ; see Fig. 4); at later times, the 2 The effect of obliquity is to flatten D E ⊥⊥ at larger w ⊥ , particularly for shorter-wavelength waves (Isenberg & Vasquez 2011).The highest-k ∥ power here is predominantly parallel (Fig. 4) so these effects are likely minor.
magnitude drops significantly into a diffuse peak around w ⊥ ≃ v th , consistent with stochastic heating.In contrast, ∂ 2 Q ∥,i /∂w ⊥ ∂w ∥ (the differential parallel heating; panels d and e), maintains a form consistent with Landau damping of kinetic Alfvén waves (KAWs) throughout the simulation.
In Fig. 4 we provide evidence that the shut off in ion heating occurs because the turbulence amplitude decreases to the point where it can no longer drive quasilinear heating by oblique ICWs.In the top panels, we show 2D spectra of δB ⊥ fluctuations, E B ⊥ (k ⊥ , k ∥ ), which are computed by filtering each field into k ⊥ bins, then interpolating these onto the exact magnetic-field lines to take a k ∥ spectrum (see Methods in S+22).Energy is concentrated at k ⊥ > k ∥ (the turbulence) and in a bump at k ⊥ ≪ k ∥ ≃ 0.8d i (parallel ICWs).Examining the time evolution from left to right, we see that the energy migrates to lower k ∥ with time (it moves downwards), with the dotted white lines indicating how this follows the critical-balance scaling chosen to align with the "peak" on the cone and is consistent across time).A corollary is that largeramplitude turbulence feeds more power into smaller-k ∥ oblique ICWs with smaller w ∥res , driving more quasilinear heating.
We quantify this effect in panel (d) by computing the energy of ICWs (E ICW ), defined as fluctuations with 0.7 ≤ k ∥ d i ≤ 2, as a function of wavevector obliquity θ = tan −1 (k ⊥ /k z ).The minimum wavenumber k ∥ d i = 0.7 is chosen to capture fluctuations with (around the time of the heating transition; middle), and t ≈ 18τA (right).The dotted white lines show the critical-balance condition EICW for different ICW obliquities, which is the energy contained in the shaded regions of panel (b) (see text).(e)-(f) Frequency spectra at different times, using the same colors as Fig. 3a, showing (e) the magnetic-field spectrum ω 2 E δB ⊥ (ω) and (f) the ratio of density to magnetic-field spectra E δn /E δB ⊥ .Because oblique ICWs involve δn fluctuations, while parallel ICWs do not, the drop in the ratio E δn /E δB ⊥ with time at high frequencies is taken as a signature of a drop in the relative power in oblique ICWs compared to parallel ICWs.
w ∥res ≲ v th , which interact with the VDF core (Isenberg & Vasquez 2011), while the angle ranges capture quasi-parallel (0 < θ < 45 • ), moderately oblique (45 • < θ < 60 • ), and highly oblique (60 • < θ < 75 • ) populations; these ranges are indicated by the shaded regions in panel (b).Despite E ⊥ decreasing continuously, the energies of the oblique-ICW populations increase slightly for t ≲ 12τ A before dropping rapidly.The similarity of this functional form with that of Q ⊥,i (t) provides good evidence for oblique ICWs being the primary driver of heating.The energy of parallel ICWs (θ ≲ 45 • ) follows a similar trend, but they are not driven directly by the turbulence, which has little power at k ∥ ≳ k ⊥ .Instead they arise because oblique-ICW quasilinear heating causes f i to increase along the resonance contours of parallel ICWs, thus emitting waves (causing instability; Kennel & Wong 1967;Chandran et al. 2010b).Their energy is driven and undamped, building up substantially for t ≲ 12τ A , then decreasing at later times as f i changes shape and renders parallel ICWs stable.Further evidence for this scenario is presented in panels (e) and (f), which show frequency spectra E(ω) of δB ⊥ and density fluctuations.These are computed from high-cadence time data at 100 spatial points, segmented into time intervals of ±2τ A = ±809Ω −1 i around different times matching Fig. 3 (we highlight the ICW frequency range where k ∥ d i ≳ 0.7, 0.3 ≲ ω/Ω i ≲ 1).In E δB ⊥ (ω), the lower-frequency fluctuations approximately satisfy E δB ⊥ ∝ ω −2 , which signifies a constant-flux cascade (Corrsin 1963;Beresnyak 2015).The spectrum drops off dramatically for ω ≳ 0.5Ω i , which is presumably where ICW heating becomes dominant.In panel (f) the ratio of density to magnetic fluctuations, E δn /E δB ⊥ , provides a diagnostic of parallel and oblique ICWs: while both parallel and oblique ICWs involve δB ⊥ fluctuations, only oblique modes involve δn.We see that at ω ≳ 0.5Ω i , E δn /E δB ⊥ decreases at intermediate times, even though at low frequencies it increases over the same time period (presumably because larger density fluctuations are driven by the forcing in balanced turbulence).This high-frequency decrease is consistent with the ratio of oblique-ICW to parallel-ICW power seen in panel (d), which drops by a factor of ≃2 between t ≈ 5τ A and t ≈ 10τ A then remains approximately constant over the rest of the simulation.
The interpretation described above predicts that Q i scales with the oblique-ICW power, which itself scales with k ∥,max d i , where k ∥,max corresponds to the small- est parallel scale accessible by the imbalanced portion of the cascade.An independent estimate of k ∥,max d i can be obtained from the spectrum using k is the break scale (Fig. 2).
This gives k ∥,max In the simulation, (z + rms ) 2 z − rms increases modestly for t ≲ 11τ A before dropping significantly (not shown), remaining well correlated with Q i and E ICW at all times, including during pre-saturated phase of S+22.This agreement demonstrates the internal consistency of the results across diverse diagnostics (E ICW , Q i , and k * ⊥ via E B ⊥ (k ⊥ )), and could prove useful for developing a closure model for helicity-barrier-mediated turbulence.
The effect of these dynamics on f i is illustrated in Fig. 5, whose top panel superimposes isocontours of f i (w ⊥ , w ∥ ) at t ≈ 18τ A over f i (w ⊥ , w ∥ ) at t ≈ 5τ A .Particles are scattered to render f i constant along the resonant contours, which are computed from G res [f i ] = 0, where With heating being driven by oblique ICWs but exciting parallel ICWs for t ≲ 11τ A , we expect f i to decrease along the oblique-ICW contours (which approximately satisfy the sin θ ≈ 1 cold-plasma relation ω Stix 1992) and to increase along the parallel-ICW contours (ω = k ∥ v A 1 − ω/Ω i ; see Isenberg 2012 for further discussion of cold-plasma relations in this context).Its evolution towards a flatter core at t ≳ 11τ A is hard to discern on f i (w ⊥ , w ∥ ), so in the lower panel we plot f i along an oblique-ICW resonant contour (thin black line).It changes from decreasing slightly to increasing along the contour, with the latter a signature of f i (w ⊥ , w ∥ ) becoming flatter in w ⊥ , likely due to stochastic heating (a w ⊥ -independent f i is shown with the dashed line; Klein & Chandran 2016;Cerri et al. 2021).Various other modifications to f i over the simulation are clear in the top panel, including the smoothing of the sharp "ridge" bordering the quasilinearly heated region of phase space at w ∥ /v th0 ≈ −0.5, and the flattening of f i (w ⊥ ) at w ∥ > 0.
Another feature of f i is the strong parallel plateau, or beam, which results from the Landau damping of k ⊥ ≫ k ∥ fluctuations at k ⊥ ρ i ≲ 1 (see also S+22; Li et al. 2010).Such fluctuations, which lie in the k ⊥ range between Alfvén waves and KAWs, propagate with phase speed >v A (Howes et al. 2006), growing a modestly super-Alfvénic plateau in f i .As the turbulence becomes balanced, such fluctuations would likely flatten f i at w ∥ ≈ −v A also (see figure 7 of Arzamasskiy et al. 2019), but this feature is not yet observable.The total parallel heating Q ∥,i , which captures this beam formation, remains almost constant throughout the full simulation (Fig. 3a inset).This suggests that the k ⊥ ρ i ∼ 1 Alfvénic fluctuations are Landau damping some fixed portion of their energy before dissipating into oblique ICWs (at earlier times) or KAW turbulence (at later times).A corollary is that there is no fundamental difference between beam formation and standard resonant parallel heating, aside from the fluctuations' imbalance.
4. DISCUSSION A complete theory for the helicity barrier would be able to predict the parameters (ε H /ε, β, etc.) at which it is operative.Low-β gyrokinetics states only that the barrier occurs unless the generalized helicity can be destroyed at least as fast as ∼ε H before k ⊥ ρ i ∼ 1 scales are reached (Meyrand et al. 2021).Although generalized helicity is not a true invariant of hybrid kinetics, we can infer from our simulation, based on the agreement between Q e and ε − ε H (Fig. 1) and the continuous evolution of the transition-range spectrum (Fig. 2) for t ≲ 11τ A , that the critical ε H /ε at β ≈ 0.4 is ≈0.2.The changes that occur once ε H /ε ≲ 0.2 also coincide with the sharp drop in Q ⊥,i (Fig. 3), suggesting that breaking the helicity barrier curtails ion heating.These results highlight the surprising robustness of the helic-ity barrier in the face of the additional complexities of a true kinetic system.However, ε H /ε decreases unrealistically fast in our simulation; understanding how the barrier evolves when ε H /ε changes on timescales comparable to the turbulent decay time should be a priority for future work.
Our results provide a helpful roadmap for understanding the radially dependent interplay between turbulence, heating, and instabilities in the β ≲ 1 solar wind.At smaller R and in faster streams, which have higher observed imbalance, we predict strong perpendicular ion heating via quasi-linear resonance, continual emission (instability) of parallel ICWs, a steep and wide ion-Larmor-scale transition range, proton beam formation, and little electron heating.At larger R and in slower streams, which have lower imbalance (σ c ≲ 0.8 in our simulation), electron heating dominates, parallel ICWs are absorbed/damped, there is no transition-range spectrum, and sharp features in the VDF are smoothed out.These correlations and features match those measured in the low-β solar wind by in situ spacecraft (e.g., Marsch 2006;Bruno et al. 2014;Zhao et al. 2021;Shi et al. 2023).
Under the assumption that, far from the Sun, faster wind is heated more than slower wind (Hansteen & Leer 1995;Totten et al. 1995;Halekas et al. 2023), our results explain qualitatively various interesting properties of observed temperature profiles.Close to the Sun, proton and minor-ion temperatures correlate positively with wind speed (Burlaga & Ogilvie 1973), while the electron temperature is negatively correlated (Marsch et al. 1989); this is natural if the highly imbalanced turbulence of faster streams has low Q e /Q i (Shi et al. 2023).At larger distances (∼1 au), the electron-temperature correlation flips to positive for U ≲ 500 km s −1 (Shi et al. 2023), as would occur if Q e /Q i increased to ≳1 as the turbulence becomes balanced at larger radii.More directly, Abraham et al. (2022) measure Q e ∝ R −2 for R ≲ 0.3 au, followed by a rapid drop in Q e at larger R.This profile, Q e ∝ R −2 , is much flatter than the "standard" total-heating profile Q ∝ R −4 (e.g., Totten et al. 1995)-a natural explanation is that Q = Q i + Q e is steeper than Q e , because Q e /Q i increases as ε H /ε decreases with R.Then, once ε H /ε ≪ 1 (around R ≃ 0.3 in this scenario), Q e /Q i saturates and Q e drops rapidly.
Together, these observations suggest that imbalanced turbulence and the helicity barrier are actively shaping global coronal and solar-wind dynamics.More generally, they highlight the crucial role of imbalance in controlling collisionless plasma thermodynamics.

Figure 2 .
Figure 2. (Top) Perpendicular spectra of the electric field (EE ⊥ ; solid lines) and magnetic field (EB ⊥ ; dashed lines) at a selection of times.The sub-panel illustrates the local spectral slopes α, computed via a fit over range ± log 10 (k ⊥ ) ≈ 0.1 around each point; the dotted lines indicate k ⊥ ρi = 1 and the representative power laws k −5/3 ⊥ and k −0.8 ⊥ .(Bottom) Evolution of the break scale k * ⊥ with imbalance 1 − σc, with the line color indicating the time as in the top panel.The dotted line shows the empirical scaling k* ⊥ ρi = 2(1 − σc) 1/4, which provides a good fit to the simulation until t ≈ 11τA.

Figure 3 .
Figure 3. (Left) Time evolution of the perpendicular diffusion energy coefficient D E ⊥⊥ , showing the transition from the quasilinear expectation D E ⊥⊥ ∝ w 2 ⊥ to a flatter profile around t ≳ 11τA.The inset shows the normalized perpendicular (orange) and parallel (purple) ion heating rates, illustrating a sharp drop in Q ⊥,i associated with the change in D E ⊥⊥ , while Q ∥,i remains almost constant.(Right) Differential heating rates, v 2 th0 ε −1 ∂ 2 Q ⊥,∥,i /∂w ⊥ ∂w ∥ , computed from ⟨E ⊥,∥ • w ⊥,∥ ⟩.The Q ⊥ profile (panels b and c) exhibits resonant structure at early times (top); the shaded region shows w ∥res for oblique ICWs with d −1

Figure 5 .
Figure5.(Top) Contours of the ion VDF fi at t ≈ 5τA (white lines and colormap) and the t ≈ 18τA (red lines).Dotted vertical lines highlight w ∥ = 0 and w ∥ = vA.The red contour levels are chosen to emphasize changes in fi and do not correspond to the white contours (see text).(Bottom) fi versus distance s along the oblique-ICW resonant contour shown with the black line in the top panel (s = 0 at w ∥ = 0).Colors are as in Fig.3a.The thin-black line is flat (Gres[fi] = 0), while the dashed line shows an fi that is independent of w ⊥ (as expected from stochastic heating).