Adaptive critical balance and firehose instability in an expanding, turbulent, collisionless plasma

Using hybrid-kinetic particle-in-cell simulation, we study the evolution of an expanding, collisionless, magnetized plasma in which strong Alfv\'enic turbulence is persistently driven. Temperature anisotropy generated adiabatically by the plasma expansion (and consequent decrease in the mean magnetic-field strength) gradually reduces the effective elasticity of the field lines, causing reductions in the linear frequency and residual energy of the Alfv\'{e}nic fluctuations. In response, these fluctuations modify their interactions and spatial anisotropy to maintain a scale-by-scale"critical balance"between their characteristic linear and nonlinear frequencies. Eventually the plasma becomes unstable to kinetic firehose instabilities, which excite rapidly growing magnetic fluctuations at ion-Larmor scales. The consequent pitch-angle scattering of particles maintains the temperature anisotropy near marginal stability, even as the turbulent plasma continues to expand. The resulting evolution of parallel and perpendicular temperatures does not satisfy double-adiabatic conservation laws, but is described accurately by a simple model that includes anomalous scattering. Our results have implications for understanding the complex interplay between macro- and micro-scale physics in various hot, dilute, astrophysical plasmas, and offer predictions concerning power spectra, residual energy, ion-Larmor-scale spectral breaks, and non-Maxwellian features in ion distribution functions that may be tested by measurements taken in high-beta regions of the solar wind.


Introduction
Many space and astrophysical plasmas are magnetized and weakly collisional, with the Larmor radii of the constituent particles being many orders of magnitude below their Coulomb mean free paths (e.g., Schekochihin & Cowley 2006). This feature results in a complex interplay between a plasma's macrophysical evolution (e.g., due to expansion, compression, or large-scale shear) and its microphysical response (e.g., departures from local thermodynamic equilibrium, triggering of kinetic instabilities) (Schekochihin et al. 2005;Kunz et al. 2014a;Hellinger & Trávníček 2015;Riquelme et al. 2015;Sironi & Narayan 2015;Squire et al. 2017;Kunz et al. 2020). This interplay becomes increasingly complex when that macrophysical evolution induces or accompanies a cascade of turbulent fluctuations down to microphysical scales, a situation thought to be ubiquitous in the solar wind, low-luminosity black hole accretion flows, and the intracluster medium (e.g., Alexandrova et al. 2013;Yuan & Narayan 2014;Simionescu et al. 2019).
In this Letter, we investigate to what extent the basic building blocks of strong, incompressible, Alfvénic turbulence -namely, the existence of a conservative cascade from large (injection) to small (dissipative) scales, the locality of interactions between turbulent fluctuations, and a scale-byscale balance between the characteristic linear oscillation time of the fluctuations and their nonlinear interaction time known as "critical balance" (Goldreich & Sridhar 1995;Mallet et al. 2015;Schekochihin 2020)-survive when subject to microphysical constraints dictated by the kinetic evolution of a collisionless plasma. Theoretical work describing magnetized turbulence in weakly collisional or collisionless plasma, but adopting a pressure-isotropic background, suggests that these organizing principles endure, with a local, conservative, Alfvénic cascade extending from macroscopic scales down to the ion-Larmor scale (Schekochihin et al. 2009). However, the assumption of an isotropic background pressure is not always justified; instead, the pressure tensor is more naturally anisotropic with respect to the magnetic field, with the evolution of field-parallel and perpendicular pressures influenced by approximate adiabatic invariance of the charged particles. How this pressure anisotropy alters the properties of Alfvénic turbulence has been a question of particular interest in recent years (e.g., Klein & Howes 2015;Kunz et al. 2015Kunz et al. , 2018Markovskii et al. 2019).
To address this question, we use results from a hybridkinetic simulation in which strong Alfvénic turbulence is driven in a collisionless, magnetized plasma undergoing steady expansion transverse to a mean magnetic field. This expansion drives pressure anisotropy in the plasma through approximate adiabatic invariance. We find that, despite the consequent decrease in the characteristic linear frequency and Alfvén ratio of the fluctuations, the Alfvénic cascade adapts to maintain critical balance. Eventually the plasma becomes unstable to kinetic firehose instabilities, which grow rapidly on ion-Larmor scales, scatter particles, and thereby impede the further production of pressure anisotropy. Even in this state, critical balance of the Alfvénic cascade persists, with the majority of the turbulent motions remaining stable.

Why Expansion?
Of the various types of macroscopic evolution that a turbulent, collisionless magnetized plasma can undergo, there are two compelling reasons to consider expansion.
First, plasma expansion on a timescale t exp much larger than the inverse cyclotron frequency Ws 1 of each particle species s (ä {e, i} for an electron-ion plasma) provides a natural way to drive temperature anisotropy, Δ s ≡ T ⊥s /T ∥s − 1 ≠ 0, where T ⊥s (T ∥s ) is the field-perpendicular (-parallel) component of the temperature of species s. For example, as plasma expands transversely to a mean magnetic ("guide") field, mass and magnetic-flux conservation dictate that the mean number density n s of each species s and the guide-field strength B g satisfy µn B L , s g 2 , where L ⊥ is the characteristic transverse size of the plasma (taken to be much larger than the thermal Larmor radius ρ s of each species; the characteristic parallel size L ∥ is held fixed). Combined with conservation of the first and second adiabatic invariants, viz. , increase. That the combination β ∥s Δ s grows increasingly negative has two important consequences. First, the effective Alfvén speed drops below the conventional Alfvén speed v A , tending toward zero as ∑ s β ∥s Δ s → −2 (at which point there is no energetic cost to bending the field). Thus, the effective tension in the magnetic-field lines is reduced, with Alfvén waves becoming unstable for ∑ s β ∥s Δ s < −2 (the "fluid firehose" threshold; Chandrasekhar et al. 1958;Parker 1958). Concurrently, when β ∥s Δ s  −1 the plasma becomes unstable to the kinetic parallel and oblique firehoses (Yoon et al. 1993). For plasma with β ∥i ≈ 2-4 and Maxwellian electrons, the oblique firehose operates when b D - Matteini et al. 2006). Both effects prompt several questions, including whether critical balance persists during the expansion, how the kinetic instabilities interact with the Alfvénic turbulence, and whether the turbulent motions themselves become unstable and disrupt the cascade.
The second reason to consider the problem of expanding Alfvénic turbulence is its relevance to the solar wind. A parcel of solar-wind plasma initially located at a large distance R ? L ⊥ , L ∥ from the Sun and moving radially outwards at speed v sw will undergo (approximately linear) expansion on a characteristic timescale t = R v exp sw (e.g., Matteini et al. 2012). Expansion is thought to play an important role in various key physical processes in the solar wind, including plasma heating, the generation of turbulence, and kinetic physics such as the production of temperature anisotropy (Velli et al. 1989;Verdini & Velli 2007;Chandran & Hollweg 2009;Matteini et al. 2013;Chandran & Perez 2019). There have therefore been many complementary investigations of expanding plasmas in the solar-wind context (e.g., Grappin et al. 1993;Liewer et al. 2001;Matteini et al. 2006;Camporeale & Burgess 2010;Hellinger 2017;Hellinger et al. 2019;Squire et al. 2020).

Hybrid-kinetic Description of Expanding Alfvénic Turbulence
We adopt a hybrid-kinetic approach to solve for the multiscale dynamics of Alfvénic turbulence in a collisionless, expanding plasma. A nonrelativistic, quasi-neutral (n ≡ n i = n e ) plasma with kinetic ions (mass m i , charge e) and massless, fluid electrons is threaded by a uniform magnetic field = B z B g gˆa nd subjected to a random, time-correlated, solenoidal driving force F(t, r) ⊥ B g . This driving is the same as described in Arzamasskiy et al. (2019); it is designed to mimic the action of random inertial forces arising from an anisotropic cascade of turbulent fluctuations at scales larger than the simulation domain. The electrons are assumed to be pressure-isotropic and isothermal with temperature T e = T i0 , the initial ion temperature. A fourth-order hyper-resistivity is used to remove magnetic energy at the smallest scales.
The subsequent evolution of this plasma is solved using the second-order-accurate, particle-in-cell code Pegasus++ (Arzamasskiy et al. 2021, in preparation), which is an optimized implementation of the algorithms detailed in Kunz et al. (2014b). Well-resolved 3D hybrid-kinetic simulations of Alfvénic turbulence are essential for modeling this problem, in particular for simultaneously capturing both the turbulent cascade above and below ion-Larmor scales and the physics of ion-firehose instabilities. That being said, our treatment of the electrons as an isothermal, isotropic fluid precludes any kinetic instabilities driven by electron temperature anisotropy (e.g., the electron firehose; Li & Habbal 2000). While the properties of inertial-range Alfvénic fluctuations and ion-scale firehose instabilities are not expected to be affected appreciably by electron kinetics, it remains an open question as to how electron anisotropy affects the sub-ion-Larmor cascade of kinetic Alfvén waves (KAWs; see Sections 3.6.2, 4.4, and 4.5 of Kunz et al. 2018). For now, we simply note that, in the near-Earth solar wind, the electrons' collisional age seems to control the electron temperature anisotropy (Salem et al. 2003) and the total temperature anisotropy at β  1 is dominated by protons (Chen et al. 2016). By modeling only a single ion species (protons), our simulations also preclude some other effects thought to be relevant in the solar wind, e.g., instabilities driven by drifting helium ions (Verscharen et al. 2019).
To model the expansion, Pegasus++ enacts a coordinate transform from a comoving, nonexpanding frame (position vector r) to the comoving expanding frame (position vector ¢ r ) using the time-dependent (diagonal) Jacobian transformation matrix L º ¶ ¶ ¢ r r t ( ) , as in the Hybrid Expanding Box (HEB) model of Hellinger & Trávníček (2005, Appendix A). Pegasus++ solves the following modified versions of Faraday's and Ohm's laws in the expanding frame for the magnetic field lL ¢ º -B B 1 and the electric field L ¢ º E E: where the primed-frame number density l ¢ º n n and ion-flow velocity L ¢ ºu u 1 , l L º det , and ¢ = t t. These fields are used to update the simulation ion-particle positions L ¢ =r r The final (velocity-dependent) term in Equation (5) is straightforwardly incorporated into the semi-implicit Boris algorithm for solving particle trajectories alongside the ¢´¢ v B p rotation. Quantities in the nonexpanding frame are easily obtained ex post facto.
The expansion is taken to be perpendicular to ẑ and linear in time: where t exp is the expansion time and I the unit dyadic. Thus the perpendicular size of the simulated plasma increases in time as ) , while the parallel size of the simulated plasma remains constant, L ∥ (t) = L ∥0 . (We denote any given quantity X evaluated at the start of the simulation by X 0 .) Magnetic-flux conservation ) . This prescription is physically relevant to the expanding solar wind at  0.1 au, on account of the solar wind's constant speed and radial direction at those distances (Verscharen et al. 2019), although our treatment of the mean magnetic field as radial is a simplifying assumption.

Physical Setup
At the start of the simulation (time t 0 ), N ppc = 10 3 simulation ion-particles per cell are drawn randomly from a stationary Maxwellian distribution with temperature T i0 and number density n 0 and placed uniformly in an elongated 3D computational domain of size r ŕ´=Ĺ L L 65 390 x y z i i 0 2 0 ( ) containing 256 2 × 1536 cells. At this box size and resolution, the captured wavenumbers are initially in the ranges k (x,y) ρ i0 ä [0.097, 12.37] and k z ρ i0 ä [0.016, 12.37]. The initial ion beta parameter is β ∥i0 = 2, representative of near-Earth conditions in the solar wind (Matteini et al. 2007). Prior to initiating expansion, steady-state Alfvénic turbulence is generated in the plasma by forcing the particles with an F(t,r) having the ) is the initial Alfvén speed, and Ω i0 ≡ eB g0 /m i c is the initial ion-cyclotron frequency. The magnitude of the force is such that critical balance is maintained for the box-scale fluctuations: u rms /v A0 ≈ L ⊥ /L ∥ , where u rms is the rms turbulent velocity. Assuming a −5/3 power-law scaling for turbulent fluctuations on scales larger than the box, the inferred perpendicular wavenumber at which the energy of the turbulent fluctuations becomes comparable to that of the guide magnetic field is r , a comparable degree of separation to that observed in the fast, β  1 solar wind (Wicks et al. 2010). This initial nonexpanding phase of the simulation lasts for five . The turbulent magnetic fields at t = 0 are visualized in Figure 1(a).
The plasma's expansion is then initiated as described in . This expansion time is comparable to the inferred Alfvén-crossing time at the outer scale~k 1 outer / of the turbulence. It also means that the turbulent heating time t t^ T L m u 3 2 5 i i heat rms 3 exp ( ) ( ) in our simulation; as a result, the thermodynamic evolution of the plasma is dominated not by turbulent heating but rather by the approximately double-adiabatic expansion and the feedback from firehose instabilities. As the plasma expands, strong Alfvénic turbulence is driven continuously such that In retrospect, our results suggest that a more realistic forcing prescription would maintain critical balance adaptively at the outer scale using v A,eff instead of v A . However, this prescription requires a priori knowledge of the temperature anisotropy's evolution to evolve v A,eff (t) and, moreover, becomes problematic if v A,eff were to approach 0. In practice, we find the only consequence of using v A (t) to determine the forcing amplitude to be a slight excess of energy in the turbulent fluctuations at the outer scale.

Results
The overall plasma evolution is summarized in Figure 2. Panel (a) shows that the box-averaged magnetic-field strength B(t) and density n(t) decrease in tandem once the expansion begins, with The rms field strength actually decreases slightly slower due to the growth of the turbulent Alfvénic fluctuations, δB rms , relative to B g (t) as the plasma expands (Figure 2(b)). This growth, also evident in Figure 1, is caused by wave-action conservation and by the build-up of residual magnetic energy in the fluctuations from the reduced energetic cost of bending field lines in a plasma with Δ i < 0. Namely, the Alfvén ratio becomes smaller than unity as the expansion proceeds, an effect that may be compensated by instead using the "kinetic normalization" , when combined with critical balance of the box-scale fluctuations, viz.ũ rms  (Figure 2(b), dashed blue line), a manifestly good fit to the data.
Another key property of Alfvénic turbulence in an expanding collisionless plasma is the decreasing characteristic frequency of the fluctuations. This feature is demonstrated by Figure 2(c), in which the red pluses track the time evolution of the energetically dominant ("peak") oscillation frequency of the fluctuations, ω peak . 5 While some decrease in ω peak is caused by the decreasing Alfvén speed, A 0 e x p 1 ( ) (dashed blue line), it is mostly due to the reduction in the effective Alfvén speed caused by β ∥i Δ i becoming increasingly negative. Indeed, the effective Alfvén frequency of the box-scale fluctuations, 2πv A,eff /L ∥ (solid blue line), matches the data well.
The production of negative temperature anisotropy during the expansion is shown in Figure 2(d). During the initial phase, the parallel (blue line) and perpendicular (red line) ion temperatures evolve approximately double-adiabatically: 0.4 f e x p , an abrupt change in the evolution of T ⊥i (t) and T ∥i (t) occurs, and the double-adiabatic predictions no longer hold. This change is coincident with Δ i decreasing sufficiently (and β ∥i increasing sufficiently-see Figure 2(e)) that Δ i  −1.4/β ∥i (see Figure 2(f)), at which point the plasma is unstable to kinetic firehose instabilities. Such firehose fluctuations, visually evident near the ion-Larmor scale in Figure 1 5 ω peak is computed using time series of high-cadence magnetic-field data recorded during the simulation at 27 fixed points in space. These series are Fourier transformed and the frequencies corresponding to the peaks of their corresponding energy spectra are algebraically averaged. To isolate the peak frequency at a particular time, a Gaussian window function (FWHM t D = t 0.2 exp ) centered at that time is applied to each series before Fourier transforming. kinetic simulations of nonexpanding, β ∥i ∼ 1 turbulence (e.g.,

Arzamasskiy et al. 2019):
μ^- ( ) ( ) in the inertial ("MHD") range, before steepening at k ⊥ ρ i  1 due to finite-Larmor-radius effects. The "break point" at which this steepening occurs, r k i break ( ) (blue curve, Figure 3(b) inset), decreases at a rate quantitatively consistent with theoretical expectations (Kunz et al. 2018, Section 3.6.4) that Figure 3(b) inset). 6 These spectral features are maintained throughout the expansion, even for t  t f . Having provided evidence that various properties of the largescale fluctuations adapt to the changing background pressure anisotropy in a manner consistent with critical balance, we now utilize the spectra in Figure 3 to show that critical balance is in fact maintained adaptively, scale by scale, as the plasma expands. We do so by computing the spectral anisotropy of the fluctuations using an approach proposed by Cho & Lazarian (2009) in which the characteristic parallel wavenumber k ∥ (k ⊥ ) of magnetic-field fluctuations with perpendicular wavenumber k ⊥ is determined from their rms parallel lengthscale (see their Equation (34)). For fluctuations with a given k ⊥ , this measure is most sensitive to the energetically dominant fluctuations with the largest k ∥ , and so the approach can be used to determine the linear frequency ω A ≡ k ∥ v A,eff of these fluctuations and compare it with their nonlinear frequency w º . In critically balanced turbulence, the turbulent energy is concentrated in a cone satisfying ω A  ω nl , with the edge of the cone having µk k 2 3  (Goldreich & Sridhar 1995). The result of this calculation is shown at different times in Figure 3(c). At t = 0, the measured spectral anisotropy in the inertial range is consistent with the critical-balance scaling µk k 2 3  . As the expansion proceeds, this scaling is maintained as the overall degree of anisotropy decreases in tandem with the decreasing aspect ratio of the plasma. Furthermore, the inset shows that ω A ≈ ω nl scale by scale; thus critical balance holds adaptively. At t ≈ t f , firehose modes (which, unlike the Alfvénic fluctuations, are not highly elongated in the field-parallel direction) emerge and bias slightly the calculated scaling of k ∥ (k ⊥ ) in the inertial range. To mitigate this bias, a weight function is applied to the magnetic field that preferentially removes firehose-unstable regions before evaluating k ∥ . Using this weight function, adaptive critical balance of the Alfvénic cascade is seen to persist. 7 In summary, no dramatic alterations to the fundamental nature of the Alfvénic turbulence are observed during expansion, even when kinetic-scale firehose modes are present. Importantly, there is no noticeable destabilization of the inertial-range Alfvénic cascade. This result is due to the efficient regulation of the box-averaged temperature anisotropy, which (as shown in Figure 2(f)) barely drops below Δ i ≈ −1.4/β ∥i . While this value of Δ i is negative enough to destabilize the plasma to kinetic firehose instabilities, it is above the "fluid" firehose instability threshold and Alfvén waves cease to propagate. The character of the kinetic-scale firehose fluctuations can be ascertained by examining the 2D Fourier spectrum of the magnetic field E B in º (Figure 4(a), top), spectral power is concentrated in the region of (k r , k z ) space that satisfies k z = k r , affirming the quasi-perpendicular nature of the Alfvénic cascade. By t = t 0.4 exp (Figure 4(a), bottom), an additional region with spectral power is clearly visible, with its centroid located at (k r ρ i , k z ρ i ) ≈ (0.4, 0.3). We associate this power with growing oblique firehose fluctuations. 8 These fluctuations can be visualized by isolating the "firehose" part δB x,f of the magnetic field using a Fourier-space mask that filters out quasi-perpendicular ( ) falls below some fraction ofẼ B0 , denoted byẼ B,cut . We use k ⊥l ρ i = 0.4, k ⊥u ρ i = 0.8, and 0 the result is qualitatively insensitive to moderate variations in these parameters. 7 The weight function at a given time t is constructed by first identifying all cells in which, when time averaged over an interval of size τ A0 /2 prior to time t, the firehose instability parameter β ∥i Δ i −1.4. These regions are then masked, with the edges of the mask smoothed by a Gaussian filter of scale 4πρ i . 8 In principle, parallel firehose fluctuations sitting atop local field-line deformations caused by the Alfvénic turbulence could also appear as oblique modes in (k r , k z ) space. However, the characteristic angular deviation of the magnetic-field lines associated with the Alfvénic turbulence is relatively small (θ A ≈ 19°), while the observed modes have θ ≈ 53°. We thus conclude that the emergent region of spectral power seen in Figure 4(a) at t = t 0.40 exp is caused by the oblique firehose instability. modes; the region of (k r , k z ) space identified as the firehose part is indicated by the shaded region in Figure 4 Figure 4(c), which shows the evolution of the magnetic energy of ion-Larmor-scale modes at different angles to the guide field, confirms that oblique firehose modes are unstable, with maximum growth rate comparable to that predicted by linear theory at β ⊥i ≈ 3.6 and 3). Parallel firehose modes (measured in region III of Figure 4(a)) are also unstable, but they have a significantly smaller amplitude than the oblique modes.
The firehose fluctuations efficiently regulate the temperature anisotropy, even though their saturated rms magnetic-field strength is much smaller than that of the Alfvénic fluctuations at equivalent wavenumbers. They do so by pitch-angle scattering the ions so that the particles' first adiabatic invariants where v ⊥ is the peculiar perpendicular velocity) are no longer conserved (see Figure 4 Figure 4(d) indicates that t n 1 exp c  for t < t f (i.e., μ is approximately conserved prefirehose), while t n~1 exp c for t  t f (i.e., μ is significantly broken by the firehose fluctuations).
A simple model for ν c may be constructed by adopting three assumptions: (i) that n(t)/n 0 ≈ B(t)/B go ; (ii) that contributions from heat fluxes and turbulent heating to the temperature anisotropy are negligible (the latter being because t t ; The agreement between this model (Figure 4(d), blue line) and ν c evaluated directly from the simulation is good, although ν c fluctuates significantly. A direct calculation of the mean μ-breaking time of ∼ 10 4 tracked particles, following Kunz et al. (2014aKunz et al. ( , 2020 and Squire et al. (2017), yields an effective collisionality n CGL c  for t  t f . Setting n n = CGL c c in the above equations leads to a simple equation for the parallel temperature, . This model is plotted in Figure 2(d); given its simplicity, its agreement with the actual result is remarkable.
The regulation of temperature anisotropy can be elucidated further by considering PDFs of the simulation data in the (β ∥i , Δ i ) phase space (e.g., Bale et al. 2009). Figure 5(a) shows these PDFs at different stages: at the expansion's start (t = 0), at t = t f , and more than a full expansion time after t = 0; the phase-space trajectory of the PDF's average is indicated in the final panel by the black solid line. In all three cases, the relatively small dispersion in β ∥i and Δ i is consistent with the small rms amplitude of the turbulent fluctuations. The temperature anisotropy clearly approaches the oblique firehose instability threshold Δ i = −1.4/β ∥i (dashed line) and subsequently evolves along marginal instability.
Despite the success of our collisionality model, the compartmentalization of all of the kinetic physics into an effective collision frequency hides some interesting emergent features in the ion distribution function f (v ∥ , v ⊥ ). Figure 5  shows the difference between f and a Maxwellian distribution with the same temperature as f at three different times during the expansion (with all velocities normalized by the initial thermal speed . For comparison, the difference between f and a bi-Maxwellian distribution with the same values of T ∥i and T ⊥i as f is also shown. Prior to the start of the expansion, the slight deficit of particles with (peculiar) parallel velocities v ∥ just below the Alfvén velocity v A (Figure 5(b), left panel) is indicative of collisionless damping of the (kinetic) Alfvénic fluctuations. Once the expansion begins, these deviations are dwarfed by the expansion-driven temperature anisotropy ( Figure 5(b), middle panel), which, on account of approximate double-adiabaticity, causes f to look like a bi-Maxwellian. However, by late times in the simulation, significant deviations from a bi-Maxwellian are evident ( Figure 5(b), right panel), a finding seen in previous studies of the firehose instability (e.g., Hellinger 2017). In particular, the distribution function integrated over perpendicular velocities, , exhibits a flattened core (Figure 5(c)); the distribution function integrated over parallel velocities, shows that the anisotropy of the distribution function at subthermal velocities is much more pronounced than in a bi-Maxwellian. These features can be attributed to resonant interactions between ions and the oblique firehose modes (Bott et al. 2021, in preparation).

Discussion
That the nonlinear interactions between Alfvénic fluctuations adapt to satisfy critical balance, even as the characteristic linear frequency of those fluctuations is reduced by pressure anisotropy, is a vivid illustration of the complex interplay between velocity space and configuration space that is central to collisionless plasma physics. This interplay is made richer at β  1 by the emergence of ion-Larmor-scale firehose fluctuations, which establish a direct link between the microscales and and perpendicular ( f (v ⊥ )) distribution functions at the same times. Dashed lines denote the corresponding f biM . macroscales by regulating the pressure anisotropy and thereby controlling the effective tension of magnetic-field lines. Despite the small-scale injection of magnetic energy by the firehose, those fluctuations are not sufficient in amplitude to contribute significantly to the magnetic power spectrum (at least perpendicular to the guide field). This finding should ease the concern expressed in Bale et al. (2009) that "these local [kinetic] instabilities ... may confuse the interpretation of solar wind magnetic power spectra." From the standpoint of the Alfvénic cascade, the most important (and potentially observable) roles played by the firehose are as a direct regulator of pressure anisotropy and an indirect mediator of adaptive critical balance and the transition to the KAW range.
The evolution of purely decaying, magnetized turbulence in an expanding, collisionless plasma with β ∥i  1 was recently investigated by Hellinger et al. (2019) using HEB simulations. In their setup, an isotropic spectrum of Alfvénically polarized waves (amplitude δB rms /B g = 0.24) was initiated inside a cubic simulation domain with 512 2 × 256 cells spanning , before transverse expansion was introduced (t = W -10 i exp 4 0 1 ) and the system evolved. The initial ion distribution function had nonzero temperature anisotropy, Δ i0 = −0.25, with β ∥i0 = 2.4. Where there is overlap with their results, we find agreement: efficient regulation of the temperature anisotropy by kinetic firehose instabilities, persistence of a quasi-perpendicular Alfvénic cascade independent of firehose fluctuations, and distortion of the particle distribution function away from a bi-Maxwellian. There are, however, two important distinctions worth highlighting. First, because of the shape of the simulation domain (L ∥ L ⊥ ) in Hellinger et al. , we estimate λ CB ≈ 0.1L ∥ ∼ ρ i , placing the entire inertial range in the weak-turbulence regime. Our demonstration of adaptive critical balance of strong Alfvénic turbulence when the distribution function is anisotropic (even unstably so) is one of our key results. Second, we followed the evolution of the turbulence for well over an expansion time, and so could confirm that the temperature anisotropy remains pinned to the kinetic firehose instability threshold as the expansion proceeds. This is an important result for solar-wind applications, because the expansion time there is comparable to the turnover time (and thus the characteristic decay time) of the outer-scale turbulent eddies.
Our conclusions may not hold for plasmas with much higher β ∥i than have been considered here. First, it is possible to show using linear theory that, if t b b W - 10 ln , then Δ i would not be regulated fast enough by the oblique firehose to remain > −2/β ∥i . In this case, v A,eff 2 would pass through 0 and the entire inertial-range Alfvénic cascade would be destabilized. For the value of t exp used in our simulation, we expect this to occur for β ∥i  50. Second, negative pressure anisotropy driven by the Alfvénic fluctuations themselves can "interrupt" the fluctuations if d b - B B i rms g 1 2  , by nullifying the restoring tension force and exciting a sea of scattering firehose fluctuations (Squire et al. 2017). An investigation of strong Alfvénic turbulence at such high beta is already underway.