Particle-in-Cell Simulations of Sunward and Anti-sunward Whistler Waves in the Solar Wind

The early spacecraft measurements at 0.3–5 au (Feldman et al. 1975; Scime et al. 1994) and recent Parker Solar Probe (PSP) measurements at 0.1–0.3 au (Halekas et al. 2021) showed that electron heat conduction in the solar wind cannot be described by the Spitzer–Härm law (Spitzer & Härm 1953). The reason is that solar wind electrons are only weakly collisional because the collisional mean-free path typically exceeds the inverse gradient scale length of electron temperature in the heliosphere (Salem et al. 2003; Bale et al. 2013; Halekas et al. 2021). The observations of electron heat flux typically bounded by a threshold dependent on local electron beta indicate that wave-particle interactions are probably regulating electron heat conduction in the solar wind, and whistler-mode waves were suggested to be the most likely wave activity involved in the regulation process (Feldman et al. 1976; Gary et al. 1999; Halekas et al. 2021). Spacecraft measurements showed that whistler-mode waves (whistler waves is used further for simplicity) are indeed present in the solar wind (e.g., Tong et al. 2019a; Kretzschmar et al. 2021), but their properties and efficiency in electron heat flux regulation are still actively investigated.

Whistler waves in the solar wind can be produced by various electron-driven instabilities (e.g., Gary et al. 1975; Gary & Feldman 1977; Verscharen et al. 2022). The velocity distribution function (VDF) of electrons in the pristine solar wind consists of a dense thermal core population contributing about 90% of total electron density and tenuous superthermal halo and strahl populations carrying most of the electron heat flux (e.g., Maksimovic et al. 2005; Štverák et al. 2009; Halekas et al. 2020; Salem et al. 2023). The core and halo populations are relatively isotropic and can be described in the plasma rest frame by sunward-drifting Maxwell and anti-sunward-drifting κ distributions, respectively. In contrast, the strahl is a highly anisotropic population collimated around the local magnetic field and streaming anti-sunward.

It was recently suggested that oblique whistler waves driven by the strahl can potentially regulate the electron heat flux (Vasko et al. 2019; Verscharen et al. 2019). Numerical simulations showed this instability can indeed suppress the electron heat flux by pitch-angle scattering the strahl and converting it into more or less isotropic halo (Roberg-Clark et al. 2019; Micera et al. 2020, 2021). However, there are currently indications that this instability does not substantially regulate electron heat flux in the solar wind. First, PSP and Helios measurements at 0.1–1 au showed the strahl parameters are statistically well below the instability threshold (Jeong et al. 2022). Second, PSP measurements at 0.1–0.5 au revealed that the radial evolution of halo and strahl densities is not consistent with the halo produced by pitch-angle scattering of the strahl (Abraham et al. 2022; Horaites & Boldyrev 2022). Third, whistler waves observed in the pristine solar wind at 0.1–1 au are typically quasi-parallel and propagate within a few tens of degrees of the local magnetic field, while the oblique waves are observed significantly less frequently (Lacombe et al. 2014; Kajdič et al. 2016; Stansby et al. 2016; Tong et al. 2019a; Berčič et al. 2021; Kretzschmar et al. 2021; Cattell et al. 2022).

The early theoretical analysis by Gary et al. (1975, 1994) showed that quasi-parallel whistler waves in the solar wind can be produced by the whistler heat flux instability (WHFI). This instability operates when core and halo populations, isotropic or parallel-anisotropic in temperature, drift relative to each other parallel to the local magnetic field. There is a heat flux parallel to the halo drift (no net current in the plasma frame though) and the fastest-growing whistler waves propagate parallel to the heat flux. A strahl population that is drifting anti-sunward does not affect the WHFI because the unstable whistler waves are resonant only with a fraction of sunward-propagating halo electrons. The recent observations showed that the WHFI indeed operates in the solar wind and produces whistler waves propagating anti-sunward (Tong et al. 2019b). Particle-in-cell simulations showed that the WHFI can indeed produce whistler waves, whose properties are consistent with solar wind observations, but cannot regulate the electron heat flux (Kuzichev et al. 2019), contrary to previous speculations (Gary & Feldman 1977; Gary et al. 1999; Gary & Li 2000).

There are currently strong indications that quasi-parallel whistler waves can be also produced by the instability associated with a perpendicular temperature anisotropy of the halo population (Jagarlamudi et al. 2020; Vasko et al. 2020). These indications consist in statistically significant observations of the halo population with perpendicular temperature anisotropy (Pierrard et al. 2016; Wilson et al. 2019; Jagarlamudi et al. 2020; Wilson et al. 2020; Salem et al. 2023) and preferential occurrence of whistler waves in association with isotropic or perpendicular anisotropic halo (Tong et al. 2019a, 2019b; Jagarlamudi et al. 2020). The recent reports of sunward and anti-sunward-propagating whistler waves in the near-Sun solar wind are also of relevance (Mozer et al. 2020; Froment et al. 2023).

In this paper, we present particle-in-cell simulations of a combined whistler heat flux and temperature anisotropy instability that is potentially operating in the solar wind and capable of producing both sunward and anti-sunward whistler waves. We address the evolution of this instability as well as its efficiency in regulating the electron heat flux. We also determine the saturated properties of the whistler waves along with their dependence on plasma parameters and demonstrate that saturated whistler wave amplitudes can be estimated using their initial linear growth rates. The applications of the presented results are discussed.

We use the particle-in-cell TRISTAN-MP code (Spitkovsky 2008) and perform 1D3V simulations restricted to whistler waves propagating parallel and antiparallel to the background magnetic field. Ions are assumed to be an immobile neutralizing background. Electrons are represented by core and halo populations, whose initial VDFs in the plasma frame in a nonrelativistic limit, which is the case in our simulations, are described by Maxwell distributions

$$\begin{eqnarray}&&{f}_{\alpha }({\boldsymbol{v}})={{ \mathcal N }}_{\alpha }\exp \left[-\displaystyle \frac{{m}_{e}{\left({v}_{| | }-{u}_{\alpha }\right)}^{2}}{2\ {T}_{\alpha }}-\displaystyle \frac{{m}_{e}{v}_{\perp }^{2}}{2\ {A}_{\alpha }\ {T}_{\alpha }}\right],\end{eqnarray} \tag{ 1 }$$

where α = c, h correspond to core and halo populations; v_∣∣ and v_⊥ are velocities parallel and perpendicular to background magnetic field; ${{ \mathcal N }}_{\alpha }={n}_{\alpha }{A}_{\alpha }^{-1}{\left({m}_{e}/2\pi {T}_{\alpha }\right)}^{3/2}$ is the normalization constant; and n_α , u_α , T_α , and A_α are, respectively, density, drift velocity, parallel temperature, and temperature anisotropy. The electron current is assumed to be zero, n_c u_c + n_h u_h = 0. The electron heat flux is parallel to the background magnetic field and carried predominantly by the halo population, q_e ≈ − n_c u_c T_h (3/2 + A_h ), because the halo is several times hotter than the core population, T_h /T_c ≈ 3–7 (e.g., Maksimovic et al. 2005; Salem et al. 2023).

The combination of core and halo populations relatively well describes the electron VDF beyond 0.3 au, where the halo density is more than 10 times larger than the strahl density (Maksimovic et al. 2005; Abraham et al. 2022; Salem et al. 2023). Although the halo population is better described by a κ distribution (e.g., Maksimovic et al. 2005), we consider a Maxwellian halo to reduce the number of free parameters. The use of a κ distribution would not affect the critical results of this study. The linear analysis of a combined whistler heat flux and temperature anisotropy instability shows that the growth rate normalized to electron cyclotron frequency ω_ce depends on the wavenumber normalized to electron inertial length c/ω_pe and the following parameters (Vasko et al. 2020):

• 1.  
  ${\beta }_{c}=8\pi {n}_{c}{T}_{c}/{B}_{0}^{2}$: core electron beta.
• 2.  
  n_c /n₀: core density relative to total electron density, n₀ = n_c + n_h .
• 3.  
  u_c /v_A : core drift velocity in units of Alfvén speed, ${v}_{A}={B}_{0}{(4\pi {n}_{0}{m}_{p})}^{-1/2}$.
• 4.  
  T_h /T_c : ratio of halo and core parallel temperatures.
• 5.  
  A_c and A_h : core and halo temperature anisotropies.

Note that linear stability as well as nonlinear evolution also depend on the ratio between electron plasma and cyclotron frequencies, but this dependence is negligible once ω_pe /ω_ce ≫ 1 (Kuzichev et al. 2019; Vasko et al. 2020). In this paper, we keep n_c /n₀ = 0.85, T_h /T_c = 6, and A_c = 1 and present numerical simulations for various combinations of β_c , u_c /v_A , and A_h . The typical values of these parameters in the solar wind are β_c = 0.1–10, ∣u_c ∣/v_A = 1–7, and A_h = 1.1–1.5 (Pierrard et al. 2016; Tong et al. 2018; Wilson et al. 2019; Jagarlamudi et al. 2020; Wilson et al. 2020; Salem et al. 2023). Table 1 presents values of these parameters used in three sets of our simulations (25 runs per set).

All the dimensionless parameters used in our simulations are realistic of the solar wind, except for the ratio between electron plasma and cyclotron frequencies. We used ω_pe /ω_ce ∼ 10, which is about 10 times smaller than in the realistic solar wind but allows us to save computational resources. Since the core electron beta depends on the ratio between plasma and cyclotron frequencies, ${\beta }_{c}=\left(2{n}_{c}{T}_{c}/{m}_{e}{n}_{0}{c}^{2}\right){\left({\omega }_{{pe}}/{\omega }_{{ce}}\right)}^{2}$, we adjusted the core electron temperature to maintain realistic core electron betas. In all our simulations, we assumed core electron temperature of T_c = 2 keV and computed the frequency ratio ω_pe /ω_ce for a given realistic core electron beta. For β_c = 0.3–3, we obtained ω_pe /ω_ce ≈ 7–20. In all simulation runs, the length of the simulation box was L ≈ 105 c/ω_ce , which is about 1300 c/ω_pe in the case of β_c = 1. The temporal and spatial integration steps were 0.09 ${\omega }_{{pe}}^{-1}$ and 0.2 c/ω_pe , both adequate to resolve the expected whistler waves. The number of particles per cell for each population was 4·10⁴. We will preface the presentation of the simulation results by linear stability analysis.

Figure 1 presents the results of linear stability analysis of whistler waves at fixed values of core electron beta and halo temperature anisotropy (β_c = 1 and A_h = 1.3) but various values of electron heat flux determined by core drift velocity u_c /v_A . Panels (a) and (b) present the dispersion curves and growth rates of whistler waves propagating parallel and antiparallel to the electron heat flux. When the electron heat flux is absent (u_c /v_A = 0), identical parallel and antiparallel whistler waves are unstable due to the halo temperature anisotropy. The presence of electron heat flux breaks the symmetry, resulting in larger growth rates of whistler waves propagating parallel to the electron heat flux. Panels (c)–(e) present the maximum growth rates along with corresponding frequencies and wavenumbers of parallel and antiparallel whistler waves unstable at various values of u_c /v_A . In the considered range of u_c /v_A values, the parameters of the fastest-growing parallel whistler waves barely vary (γ₊/ω_ce ≈ 0.01, ω₊/ω_ce ≈ 0.1 and k₊ c/ω_pe ≈ 0.34). In contrast, the maximum growth rate of antiparallel whistler waves monotonously decreases from γ₋/ω_ce ≈ 0.01 to 10⁻³. The frequency and wavenumber also monotonously decrease by a factor of a few.

Zoom In Zoom Out Reset image size

Figure 2 presents the results of a simulation run performed at β_c = 1, A_h = 1.3, and u_c /v_A = − 3. We consider the dynamics of magnetic field $\delta {\boldsymbol{B}}(x,t)=\delta {B}_{y}(x,t)\hat{y}+\delta {B}_{z}(x,t)\hat{z}$ perpendicular to background magnetic field ${B}_{0}\hat{x}$. Panel (a) presents the magnetic field magnitude δ B(x, t)/B₀ and demonstrates the growth of magnetic field fluctuations propagating both parallel and antiparallel to the electron heat flux. Using Fourier transform, δ B (x, t) = ∫δ B _{k ω} e^{i(kx−ω t)} dkd ω, we decompose magnetic field fluctuations into those propagating parallel and antiparallel to the electron heat flux, δ B (x, t) = δ B ₊(x, t) + δ B ₋(x, t), where δ B ₊(x, t) = ∫_ω/k>0 δ B _{k ω} e^{i(kx−ω t)} dkd ω and δ B ₋(x, t) = ∫_ω/k<0 δ B _{k ω} e^{i(kx−ω t)} dkd ω. Both δ B ₊ and δ B ₋ have right-hand polarization (not shown here) and correspond to parallel and antiparallel whistler waves expected based on linear stability analysis (Figure 1). Panels (b) and (c) show that over the computation time, the parallel and antiparallel whistler waves reach peak amplitudes of about 0.1B₀ and 0.05B₀, respectively.

Zoom In Zoom Out Reset image size

Figure 3 presents averaged amplitudes and growth rates of the parallel and antiparallel whistler waves. Panel (a) shows the temporal evolution of magnetic field amplitudes averaged over the simulation box, $\langle \delta {B}_{\pm }\rangle ={\left[{L}^{-1}{\int }_{0}^{L}| \delta {{\boldsymbol{B}}}_{\pm }{| }^{2}\ {dx}\right]}^{1/2}$, and shows that within the computation time, the parallel and antiparallel whistler waves saturate, and the saturated amplitudes are ${B}_{w}^{+}/{B}_{0}\approx 0.04$ and ${B}_{w}^{-}/{B}_{0}\approx 0.02$. Panel (b) presents the temporal evolution of the whistler wave growth rates computed as $d/{dt}\left[\mathrm{ln}\ \langle \delta {B}_{\pm }\rangle \right]$. The initial growth rates of the parallel and antiparallel whistler waves are, respectively, around 0.01 and 0.004 ω_ce , both consistent within a few tens of percent with the linear stability results (Figure 1(c)).

Zoom In Zoom Out Reset image size

Figure 4 presents the results of simulation runs performed at β_c = 1 and A_h = 1.3 but various values of u_c /v_A  are indicated in panels (c)–(e) in Figure 1. In all these simulation runs, parallel and antiparallel whistler waves saturated within the computation time, and we determined the averaged amplitudes ${B}_{w}^{+}$ and ${B}_{w}^{-}$ reached by the end of each simulation run. Panel (a) shows that the saturated amplitude ${B}_{w}^{+}$ of parallel whistler waves is around 0.04B₀ and varies by less than several tens of percent over the considered range of u_c /v_A values. In contrast, the saturated amplitude ${B}_{w}^{-}$ of antiparallel whistler waves monotonously decreases from 0.025 to 0.005B₀. Interestingly, according to panel (a), the dependencies of the saturated amplitudes ${B}_{w}^{+}$ and ${B}_{w}^{-}$ on u_c /v_A are almost identical with those of initial linear growth rates γ₊ and γ₋. Panel (b) demonstrates that the ratio ${B}_{w}^{+}/{B}_{w}^{-}$ is closely correlated with γ₊/γ₋ and the best power-law fit is ${B}_{w}^{+}/{B}_{w}^{-}\approx {({\gamma }_{+}/{\gamma }_{-})}^{0.73}$. This relation naturally predicts lower saturation amplitudes of antiparallel whistler waves compared to parallel whistler waves because the former always have lower initial linear growth rates (Figure 1).

Zoom In Zoom Out Reset image size

Figure 5 presents the temporal evolution of the electron heat flux in the considered simulation runs. We demonstrate the electron heat flux variation $\delta {q}_{e}(t)=\left[{q}_{e}(t)-{q}_{e}(0)\right]/{q}_{e}(0)$ in percent, where q_e (t) is the electron heat flux averaged over the simulation box and q_e (0) is the initial heat flux value. The electron heat flux suppression is most efficient, δ q_e ≈ −8%, in the simulation run with u_c /v_A = −1.5, while the efficiency drops to about 1% for u_c /v_A = −7.5. Note that electron heat flux suppression of 8% is relatively large compared to less than 2% suppression observed in the simulations of a pure WHFI (Kuzichev et al. 2019). The individual contributions of parallel or antiparallel whistler waves to the observed electron heat flux suppression can be revealed by the analysis of electron VDF evolution.

Zoom In Zoom Out Reset image size

Figure 6 presents the analysis of electron VDF evolution for the simulation run with u_c /v_A = −1.5, exhibiting the largest electron heat flux suppression of δ q_e ≈ −8%. Panels (a) and (b) demonstrate the initial and saturated VDFs of halo electrons. The relative difference of these VDFs in panel (c) shows that the most substantial VDF evolution occurs around ${v}_{| | }\approx \left({\omega }_{+}-{\omega }_{c}\right)/{k}_{+}\lt 0$ and ${v}_{| | }\approx \left({\omega }_{-}-{\omega }_{c}\right)/{k}_{-}\gt 0$ corresponding to electrons resonant with the parallel and antiparallel whistler waves, respectively. Around these resonant velocities, the parallel and antiparallel whistler waves drive electron transport from parallel speeds smaller to those larger than the resonant speed. The transport process is reflected by regions shaded blue and yellow in panel (c). Panel (d) presents the linear density of the initial and saturated heat flux of halo electrons in the parallel velocity space, ${Q}_{h}({v}_{| | })=0.5\int {v}_{| | }\ {m}_{e}{{\boldsymbol{v}}}^{2}\cdot {f}_{h}\left({v}_{| | },{{\boldsymbol{v}}}_{\perp }\right){d}^{2}{v}_{\perp }$, while panel (e) shows their difference δ Q_h (v_∣∣). In accordance with the electron transport observed in panel (c), δ Q_h (v_∣∣) is mostly negative at v_∣∣ < 0 and mostly positive at v_∣∣ > 0. The corresponding variations of the partial heat fluxes are ${\int }_{{v}_{| | }\lt 0}\delta {Q}_{h}{{dv}}_{| | }\,\approx -0.41\ {q}_{e}(0)$ and ${\int }_{{v}_{| | }\gt 0}\delta {Q}_{h}{{dv}}_{| | }\approx 0.32\ {q}_{e}(0)$, where q_e (0) is the initial electron heat flux value. The net variation of the heat flux of halo electrons, ∫δ Q_h dv_∣∣ ≈ − 0.09 q_e (0), does not exactly coincide with δ q_e ≈ − 8% reported in Figure 5 because δ q_e also includes the heat flux variation of core electrons. Since electron transport at v_∣∣ > 0 and v_∣∣ < 0 is driven by, respectively, antiparallel and parallel whistler waves, we conclude that the electron heat flux is actually decreased by parallel whistler waves and increased by antiparallel whistler waves. The net effect is the electron heat flux suppression by a few to 10% reported in Figure 5.

Zoom In Zoom Out Reset image size

Figure 7 presents averaged amplitudes ${B}_{w}^{+}$ and ${B}_{w}^{-}$ of parallel and antiparallel whistler waves for all of the 75 simulation runs (Table 1). Note that we demonstrate a whistler wave amplitude only if the initial linear growth rate of the whistler wave is larger than 10⁻³ ω_ce ; otherwise, the computation time of $5000\ {\omega }_{{ce}}^{-1}$ is insufficient for whistler waves to saturate. For this reason, the number of points corresponding to parallel and antiparallel whistler waves in panels (a)–(c) can be different and also less than 25. Panels (a)–(c) show that for a fixed core electron beta β_c , the saturated amplitudes are larger for larger halo temperature anisotropies A_h and increase by a factor of a few between A_h = 1.1 and 1.5. Furthermore, both parallel and antiparallel whistler waves tend to saturate at larger amplitudes for larger core electron betas. For identical anisotropies, saturated amplitudes increase by a factor of a few between β_c = 0.3 and 3. The observed dependencies of the saturated amplitudes on the halo temperature anisotropy and other plasma parameters can be deduced from a more fundamental relation to be presented below.

Zoom In Zoom Out Reset image size

Figure 8 shows that for every fixed core electron beta, the saturated amplitudes of parallel and antiparallel whistler waves are correlated with their initial linear growth rates. The observed trends can be fitted to power-law functions

$$\begin{eqnarray}&&{B}_{w}^{\pm }/{B}_{0}={C}_{\pm }{\left({\gamma }_{\pm }/{\omega }_{{ce}}\right)}^{{\nu }_{\pm }},\end{eqnarray} \tag{ 2 }$$

where the best-fit parameters C_± and ν_± are indicated in the panels and presented in Table 2. The power-law indices and multipliers corresponding to parallel and antiparallel whistler waves are different and vary with core electron beta. For the considered core electron betas, the power-law indices vary in a relatively narrow range, 0.6 ≲ ν_± ≲ 0.9. The fundamental relation between the saturated amplitude and the initial linear growth rate naturally predicts larger whistler wave amplitudes for larger electron betas and temperature anisotropies because the increase of these plasma parameters results in larger linear growth rates (Vasko et al. 2020). This relation also indicates that antiparallel whistler waves are expected to saturate at lower amplitudes than parallel whistler waves because the presence of electron heat flux results in smaller growth rates of the former (Figure 1). Equation (2) also naturally explains the correlations reported in Section 3 (Figure 4).

Zoom In Zoom Out Reset image size

Figure 9 demonstrates the efficiency of electron heat flux suppression quantified by the relative heat flux variation δ q_e reached at the saturation stage. We present δ q_e only if both parallel and antiparallel whistler waves saturate within the computation time. Panel (a) shows that electron heat flux suppression is within about 10% for low-temperature anisotropies (A_h ≲ 1.1) but can be as large as 60% for A_h = 1.5. At a fixed halo temperature anisotropy, the electron heat flux suppression is more efficient for larger core electron betas. The efficiency of electron heat flux suppression is expected to correlate with the saturated whistler wave amplitudes since the latter are positively correlated with both electron beta and temperature anisotropy (Figure 5). Panel (b) shows that electron heat flux suppression δ q_e is positively correlated with the saturated amplitudes of parallel and antiparallel whistler waves. The heat flux suppression is within 10% for ${B}_{w}^{\pm }/{B}_{0}\lesssim 0.01$ but can be as large as 10%–60% for ${B}_{w}^{\pm }/{B}_{0}\gtrsim 0.01$.

Zoom In Zoom Out Reset image size

We address spectral properties of the whistler waves by computing power spectral densities of parallel and antiparallel whistler waves over the saturation stage, ${\mathrm{PSD}}_{\omega }^{\pm }=\langle {\left|\int \delta {{\boldsymbol{B}}}_{\pm }(x,t){e}^{i\omega t}{dt}\right|}^{2}\rangle$, where the integration is over the saturation stage, while the brackets stand for spatial averaging over the simulation box. The Gaussian fitting of power spectral densities ${\mathrm{PSD}}_{\omega }^{\pm }$ was done to determine the central frequency and spectral width of parallel and antiparallel whistler waves. Figure 10(a) demonstrates that both parallel and antiparallel whistler waves have comparable relative spectral widths, Δω_±/ω_± ∼ 0.3–0.8. It is noteworthy that the relative spectral widths tend to be larger for larger core electron betas and positively correlated with an initial linear growth rate of the whistler waves. The frequency of the saturated whistler waves is consistent within a few tens of percent with the frequency of the initially fastest-growing whistler waves (not shown here).

Zoom In Zoom Out Reset image size

We presented the first particle-in-cell simulations of a combined whistler heat flux and temperature anisotropy instability, which generalize our previous simulations of a pure WHFI with an isotropic halo population (Kuzichev et al. 2019). In contrast to the pure WHFI capable of producing only whistler waves propagating parallel to the electron heat flux (anti-sunward), the combined whistler heat flux and temperature anisotropy instability produces both whistler waves propagating parallel (anti-sunward) and antiparallel (sunward) to the electron heat flux. We believe the combined whistler heat flux and temperature anisotropy instability may operate in the solar wind because there is statistically significant occurrence of the halo population with perpendicular temperature anisotropy (Pierrard et al. 2016; Jagarlamudi et al. 2020; Salem et al. 2023) and preferential occurrence of whistler waves in association with isotropic or perpendicular anisotropic halo (Tong et al. 2019a, 2019b; Jagarlamudi et al. 2020). The small and still uncertain occurrence of sunward whistler waves reported in the solar wind (Kretzschmar et al. 2021) may be caused by smaller saturated amplitudes and a higher probability of these waves being obscured by solar wind turbulence.

We showed that the saturated whistler wave amplitudes are correlated with their initial linear growth rates, ${B}_{w}/{B}_{0}\approx C{(\gamma /{\omega }_{{ce}})}^{\nu }$, where parameters C and ν are a bit different for sunward and anti-sunward whistler waves. For typical solar wind conditions considered in our simulations, the power-law index varies in a relatively narrow range, 0.6 ≲ ν ≲ 0.9 (Table 2). A similar scaling relation, though revealed using a few simulation runs, was previously reported for anti-sunward whistler waves produced by the pure WHFI (Kuzichev et al. 2019). Whistler waves in our simulations saturated at B_w /B₀ ∼ 0.01 because we required sufficiently high initial linear growth rates, γ/ω_ce ≳ 10⁻³, to save computational resources (Figure 8). Whistler waves with such high amplitudes rarely occur in the solar wind. The observed amplitudes are typically around 10⁻³ B₀ (Tong et al. 2019a; Jagarlamudi et al. 2021; Kretzschmar et al. 2021) and hence, correspond to initial growth rates of about 10⁻⁴ ω_ce . Nevertheless, we believe the revealed scaling relation is valid in a wide range of initial growth rates, and its predictions can be compared with solar wind observations.

The revealed scaling relation ${B}_{w}/{B}_{0}\approx C{(\gamma /{\omega }_{{ce}})}^{\nu }$ predicts a positive correlation of whistler wave amplitudes with electron beta and halo temperature anisotropy since the linear growth rates are larger for larger electron betas and temperature anisotropies (Vasko et al. 2020). The corresponding positive correlations were indeed observed in the solar wind (Tong et al. 2019a). We also showed that the whistler waves have relatively large spectral widths positively correlated with electron beta (Figure 10(a)). Similar values of the spectral width and its positive correlation with electron beta were reported for whistler waves observed in the solar wind (Tong et al. 2019a). Thus, the properties of whistler waves produced by the combined whistler heat flux and temperature anisotropy instability are basically consistent with those reported in the solar wind.

The revealed scaling relation (Equation (2)) does not correspond to whistler wave saturation via nonlinear trapping of cyclotron resonant electrons, which would occur once the bounce frequency of trapped electrons, ${\omega }_{T}\approx {\left({{kv}}_{\perp }{{eB}}_{w}/{m}_{e}c\right)}^{1/2}$, is comparable to the initial linear growth rate γ, where ${v}_{\perp }\sim {(2{T}_{h}/{m}_{e})}^{1/2}$ is a typical perpendicular speed of resonant electrons (Davidson et al. 1972; Karpman 1974; Shklyar 2011). The saturation via nonlinear trapping would result in ${B}_{w}/{B}_{0}\propto {(\gamma /{\omega }_{{ce}})}^{2}$ that is not consistent with the observed scaling relation. The nonlinear trapping does not occur when whistler waves have a sufficiently large spectral width or low amplitude such that the resonant velocity v_∣∣ = (ω − ω_ce )/k is distributed in a range much wider than the nonlinear resonance width ω_T /k (Sagdeev & Galeev 1969; Karpman 1974). In this case, the nonlinear evolution and saturation of whistler waves could be described within the quasi-linear theory. Thus, the quasi-linear description applies, provided that

$$\begin{eqnarray*}&&\displaystyle \frac{\partial }{\partial \omega }\left(\displaystyle \frac{\omega -{\omega }_{{ce}}}{k}\right){\rm{\Delta }}\omega \gg \displaystyle \frac{{\omega }_{T}}{k},\end{eqnarray*}$$

which can be rewritten as follows:

$$\begin{eqnarray}&&\displaystyle \frac{{\rm{\Delta }}\omega }{\omega }\gg {\left(\displaystyle \frac{{B}_{w}}{{B}_{0}}\right)}^{1/2}{\left(\displaystyle \frac{\beta \ \omega }{{\omega }_{{ce}}-\omega }\right)}^{1/4},\end{eqnarray} \tag{ 3 }$$

where β = β_h n₀/n_h and ${\beta }_{h}=8\pi {n}_{h}{T}_{h}/{B}_{0}^{2}$ is the halo electron beta. We computed the ratio between the left- and right-hand sides of Equation (3) for sunward and anti-sunward whistler waves observed in our simulations.

Figure 10(b) shows that both sunward and anti-sunward whistler waves have relatively large spectral widths or low amplitudes to make quasi-linear descriptions applicable for initial growth rates realistic of the solar wind. The quasi-linear description may fail for growth rates substantially exceeding 0.01ω_ce , which are, however, not realistic of solar wind plasma. The simulation results are consistent with spacecraft observations, which also demonstrated the applicability of quasi-linear descriptions for whistler waves in the solar wind (Tong et al. 2019a). The above quantitative estimates allow us to state that the phase-space transport of electrons around resonant velocities v_∣∣ ≈ (ω_± − ω_c )/k_± observed in Figure 6 is nothing but quasi-linear diffusion, which results in the instability saturation through the formation of plateaus around the resonant velocities (Sagdeev & Galeev 1969). Similar quasi-linear saturation via a plateau formation was previously observed for the pure WHFI (Kuzichev et al. 2019). We will present the derivation of the scaling relation (Equation (2)) using quasi-linear computations in a separate study, where its validity for both Maxwell and κ distributions of the halo population will be demonstrated.

The presented simulations demonstrate that the combined whistler heat flux and temperature anisotropy instability is much more efficient in the electron heat flux suppression than the pure WHFI considered previously (Kuzichev et al. 2019). We showed that the electron heat flux suppression occurs due to anti-sunward whistler waves, while sunward whistler waves actually increase the electron heat flux. These different effects of anti-sunward and sunward whistler waves could be foreseen based on the expected quasi-linear diffusion of electrons around the resonant velocities (Figure 7). We demonstrated that the net electron heat flux suppression is particularly efficient for large electron betas and halo temperature anisotropies. For example, the electron heat flux can be suppressed by up to 30%–60% for β_c ≳ 3 and A_h ≳ 1.3 (Figure 9(a)). Note that the more efficient electron heat flux suppression for higher electron betas is consistent with spacecraft observations (e.g., Gary et al. 1999; Tong et al. 2019a; Halekas et al. 2021), while the effect of halo temperature anisotropy has not been addressed experimentally yet.

The efficiency of electron heat flux suppression is positively correlated with the saturated amplitudes of sunward and anti-sunward whistler waves. The electron heat flux can be suppressed by more than 10%–60% for B_w /B₀ ≳ 0.01 (Figure 10(b)). This correlation allows the observed amplitudes of sunward and anti-sunward whistler waves to serve as the indicator of the efficiency of local electron heat flux suppression. Note that the amplitude of anti-sunward whistler waves alone cannot indicate the efficiency of electron heat flux suppression since anti-sunward whistler waves are also produced by other instabilities proved to be inefficient in electron heat flux suppression (Tong et al. 2019b; Kuzichev et al. 2019). Since whistler waves in the solar wind have typical amplitudes below a few percent of the background magnetic field, we expect the typical efficiency of local electron heat flux suppression to be within about 10%. The recent reports of sunward whistler waves with amplitudes of 0.1B₀ (Mozer et al. 2020; Froment et al. 2023) indicate, however, that the local electron heat flux suppression may occasionally exceed 60%.

The presented results are directly applicable for the solar wind beyond about 0.3 au, where the halo density substantially dominates the density of the strahl population not considered in this study. The considered instability may therefore contribute to the electron heat flux suppression observed beyond 0.3 au (e.g., Feldman et al. 1975; Scime et al. 1994; Tong et al. 2019a). The results relevant to anti-sunward whistler waves may actually apply even for the solar wind within 0.3 au since strahl electrons do not resonantly interact with anti-sunward whistler waves propagating parallel to the local magnetic field. No properties of anti-sunward whistler waves would change in the presence of strahl electrons, and these waves would still scatter resonant halo electrons and suppress their heat flux. In turn, the presence of strahl electrons would definitely affect sunward whistler waves and could potentially stabilize these waves. The effect of strahl electrons on the electron heat flux suppression by the combined whistler heat flux and temperature anisotropy instability remains to be established.

We presented particle-in-cell simulations of the combined whistler heat flux and temperature anisotropy instability that is potentially operating in the solar wind and capable of producing whistler waves propagating parallel (anti-sunward) and antiparallel (sunward) to the electron heat flux. The results of these simulations can be summarized as follows

• 1.  
  The instability saturates through quasi-linear diffusion of electrons and the formation of plateaus around the resonant velocities. The quasi-linear description can be applied because both sunward and anti-sunward whistler waves have sufficiently large spectral widths or, equivalently, relatively low amplitudes.
• 2.  
  The saturated amplitudes of both sunward and anti-sunward whistler waves strongly correlate with their initial linear growth rates, ${B}_{w}/{B}_{0}\approx C{(\gamma /{\omega }_{{ce}})}^{\nu }$, where for typical electron betas, we have C ∼ 1 and 0.6 ≲ ν ≲ 0.9 (Table 2). The scaling relation explains the positive correlation of the saturated whistler wave amplitudes with electron beta and halo temperature anisotropy.
• 3.  
  The revealed correlations of saturated whistler wave amplitudes and spectral widths with electron parameters (electron beta and temperature anisotropy) are consistent with solar wind observations.
• 4.  
  The anti-sunward whistler waves result in an electron heat flux decrease, while the sunward whistler waves actually lead to an electron heat flux increase. The net effect is the electron heat flux suppression that is more efficient for higher core electron betas β_c and halo temperature anisotropies A_h . The electron heat flux can be suppressed by 30%–60% for β_c ≳ 3 and A_h ≳ 1.3.
• 5.  
  The efficiency of electron heat flux suppression is positively correlated with saturated amplitudes of parallel and antiparallel whistler waves. The electron heat flux can be suppressed by up to 10%–60%, provided that the saturated whistler wave amplitudes exceed about 1% of the background magnetic field.

We expect that the presented results will stimulate experimental analysis of the occurrence of the combined whistler heat flux and temperature anisotropy instability and allow establishing the contribution of this instability to the electron heat flux suppression in the solar wind.

The work of I.K., I.V., and A.A. was supported by NASA grants HGI 80NSSC21K0581 and HSR 80NSSC23K0100. We would like to acknowledge high-performance computing support from Cheyenne (doi:10.5065/D6RX99HX) provided by NCAR's Computational and Information Systems Laboratory, sponsored by NSF grant No. 1502923. I.V. thanks the International Space Science Institute, Bern, Switzerland for supporting the working group on "Heliospheric Energy Budget: From Kinetic Scales to Global Solar Wind Dynamics."