PIC Simulations of Velocity-space Instabilities in a Decreasing Magnetic Field: Viscosity and Thermal Conduction

In low-collisionality plasmas, the change in the magnitude of the local magnetic field ($B\equiv | {\boldsymbol{B}}|$) generically drives a pressure anisotropy with ${p}_{\parallel ,j}\ne {p}_{\perp ,j}$ (where ${p}_{\perp ,j}$ and ${p}_{\parallel ,j}$ correspond to the pressure of species j perpendicular and parallel to ${\boldsymbol{B}}$). This is a consequence of the adiabatic invariance of the magnetic moment of particles, ${\mu }_{j}\equiv {v}_{\perp ,j}^{2}/B$, in the absence of collisions (where ${v}_{\perp ,j}$ is the velocity of species j perpendicular to ${\boldsymbol{B}}$).

These pressure anisotropies can trigger various velocity-space instabilities, which are in principle expected to pitch-angle scatter the particles, to some extent mimicking the effect of collisions. The combined effect of pressure anisotropies and velocity-space instabilities can affect various large-scale properties of the plasma, including its effective viscosity (Sharma et al. 2006; Squire et al. 2017) and thermal conductivity (see, e.g., Komarov et al. 2016; Riquelme et al. 2016). This weakly collisional behavior is expected to be important in several astrophysical systems, including low-luminosity accretion flows around compact objects (Sharma et al. 2007), the intracluster medium (ICM) (Schekochihin et al. 2005; Lyutikov 2007), and the heliosphere (Maruca et al. 2011; Remya et al. 2013).

In a previous work (Riquelme et al. 2016) we studied how the plasma viscosity and thermal conductivity are affected by an increase in B, which naturally drives ${p}_{\perp ,j}\gt {p}_{| | ,j}$. In this paper we study the opposite case, where B decreases and ${p}_{\perp ,j}\lt {p}_{| | ,j}$. In this case several velocity-space instabilities can be excited, ultimately regulating the extent to which the ${p}_{\perp ,j}\lt {p}_{| | ,j}$ anisotropy can grow. When only the electron dynamics is considered, two types of plasma waves are expected to be driven unstable by the electron pressure anisotropy: (i) the oblique electron firehose (OEF) modes, which are purely growing modes, and (ii) the Alfvén/ion-cyclotron (A/IC) modes, which are quasi-parallel, propagating waves, driven unstable by cyclotron-resonant electrons (Li & Habbal 2000; Camporeale & Burgess 2008).⁵ Similarly, in the presence of an ion pressure anisotropy ${p}_{\perp ,i}\lt {p}_{| | ,i}$, there are also two types of modes that can grow unstable: (i) the oblique ion firehose (OIF) modes, which are purely growing modes, and (ii) the fast magnetosonic/whistler (FM/W) modes, which are quasi-parallel, propagating waves, excited by cyclotron-resonant ions (Quest & Shapiro 1996; Gary et al. 1998; Hellinger & Matsumoto 2000).⁶

In this work we studied the nonlinear, saturated properties of these instabilities, making use of particle-in-cell (PIC) simulations. This is achieved by continuously decreasing the strength of the background magnetic field by externally imposing a shear motion in the plasma. This setup is interesting since in realistic astrophysical scenarios the pressure anisotropies are expected to be driven (via decreasing B) for a time significantly longer than the initial regime where the instabilities grow exponentially.

Previous works have already studied this long-term regime by simulating an expanding (instead of shearing) plasma. These works have used both hybrid-PIC simulations, which focused on the evolution of the ion anisotropy-driven instabilities (Matteini et al. 2006; Hellinger & Travnicek 2008), and PIC simulations that mainly captured the role of electron anisotropy-driven modes (Camporeale & Burgess 2010). Thus, our work is intended to study the combined effect of the electron and ion pressure anisotropies on the nonlinear, saturated regime of the different unstable modes. This aspect of our study is motivated in part by previous linear dispersion analyses showing that the electron pressure anisotropy can significantly influence the evolution of both the FM/W and OIF modes (Michno et al. 2014; Maneva et al. 2016).

There are two important applications of our work. One is to quantify the so-called "anisotropic viscosity" of the plasma, which is controlled by the pressure anisotropies of the particles. This viscosity is believed to contribute significantly to the heating of electrons and ions in accretion disks and other low-collisionality plasmas (Sharma et al. 2006, 2007; Squire et al. 2017). Also, the nonlinear evolution of the different velocity-space instabilities sets the pitch-angle scattering rate of electrons, which is key to determining their mean free path and therefore the thermal conductivity of the plasma.

The paper is organized as follows. In Section 2 we describe our simulation setup and strategy. In Section 3 we determine the saturated pressure anisotropy ${\rm{\Delta }}{p}_{j}$ of ions and electrons, describing the physical mechanisms responsible. In Section 4 we quantify the ion and electron heating. In Section 5 we measure the mean free paths of electrons and ions, and determine their dependence on the physical parameters of the plasma. In Section 6 we summarize our results and discuss their implications for various low-collisionality astrophysical plasmas.

We use the electromagnetic, relativistic PIC code TRISTAN-MP (Buneman 1993; Spitkovsky 2005) in two dimensions. The simulation box consists of a square box in the x–y plane, containing an initially isotropic plasma with a homogeneous initial magnetic field ${{\boldsymbol{B}}}_{0}$. We simulate a decreasing magnetic field by imposing a velocity shear given by ${\boldsymbol{v}}=-{sx}\hat{y}$, where s is the shear rate of the plasma and x is the distance along $\hat{x}$ ($\hat{x}$ and $\hat{y}$ are the unit vectors parallel to the x and y axes, respectively). From flux conservation, the x and y components of the mean field evolve as $d\langle {B}_{x}\rangle /{dt}=0$ and $d\langle {B}_{y}\rangle /{dt}=-s\langle {B}_{x}\rangle$. Thus, if $\langle {B}_{x}\rangle$ and $\langle {B}_{y}\rangle$ are positive, there will be a decrease in $\langle {B}_{y}\rangle$ and therefore in $| \langle {\boldsymbol{B}}\rangle |$. Therefore, we initially choose ${{\boldsymbol{B}}}_{0}\propto \hat{x}+3.3\hat{y}$, which guarantees a decrease in $| \langle {\boldsymbol{B}}\rangle |$ and a ${p}_{\perp ,j}\lt {p}_{| | ,j}$ anisotropy during a simulation time $\sim 3{{\rm{s}}}^{-1}$.

By resolving the x–y plane, our simulations can capture the quasi-parallel A/IC and FM/W modes, as well as the oblique OEF and OIF modes with their wavevectors ${\boldsymbol{k}}$ forming any angle with the mean magnetic field $\langle {\boldsymbol{B}}\rangle$. The key parameters in our simulations are: the particle magnetization, quantified by the ratio between the initial cyclotron frequency of each species and the shear rate of the plasma, ${\omega }_{c,j}/s$ ($j=i,e$), and the ratio of ion to electron mass, ${m}_{i}/{m}_{e}$. In typical astrophysical cases, ${\omega }_{c,j}\gg s$ and ${m}_{i}/{m}_{e}=1836$. Due to computational constraints, however, we will use values of ${\omega }_{c,j}/s$ and ${m}_{i}/{m}_{e}$ much larger than unity, but still much smaller than expected in real environments. This limitation will be taken into account when applying our simulation results to relevant astrophysical cases.

Our simulations initially have ${\beta }_{i}={\beta }_{e}=2$ (${\beta }_{j}\equiv 8\pi {p}_{j}/| {\boldsymbol{B}}{| }^{2}$). In almost all of our runs ${k}_{{\rm{B}}}{T}_{e}/{m}_{e}{c}^{2}=0.28$, which implies ${\omega }_{c,e}/{\omega }_{p,e}=0.53$ (where k_B, T_e, and ${\omega }_{p,e}$ are Boltzmann's constant, the electron temperature, and the electron plasma frequency, respectively). We will change our simulation conditions by varying ${\omega }_{c,e}/s$ and ${m}_{i}/{m}_{e}$ (which uniquely fix ${\omega }_{c,i}/s$ and ${k}_{{\rm{B}}}{T}_{i}/{m}_{i}{c}^{2}$). Some of our simulations use "infinite mass ions" (the ions are technically immobile, so they just provide a neutralizing charge), with the goal of focusing on the electron-scale physics. These provide a useful contrast with our finite ${m}_{i}/{m}_{e}$ runs and allow us to isolate the impact of ion physics on the electrons. The numerical parameters in our simulations will be ${N}_{{\rm{ppc}}}$ (number of particles per cell), $c/{\omega }_{p,e}/{{\rm{\Delta }}}_{x}$ (the electron skin depth in terms of grid size), $L/{R}_{L,i}$ (box size in terms of the initial ion Larmor radius for runs with finite ${m}_{i}/{m}_{e};$ ${R}_{L,i}={v}_{\mathrm{th},i}/{\omega }_{c,i}$, where ${v}_{\mathrm{th},i}\,=\sqrt{{k}_{{\rm{B}}}{T}_{i}/{m}_{i}}$), and $L/{R}_{L,e}$ (box size in terms of the initial electron Larmor radius for runs with infinite ${m}_{i}/{m}_{e}$). Table 1 shows a summary of our key runs. We ran a series of simulations ensuring that the numerical parameters (e.g., different ${N}_{{\rm{ppc}}}$) do not significantly affect our results. Note that most runs used just for numerical convergence are not in Table 1.

Table 1.  Physical and Numerical Parameters of the Simulations

Runs	${m}_{i}/{m}_{e}$	${\omega }_{c,e}/s$	${k}_{{\rm{B}}}{T}_{e}/{m}_{e}{c}^{2}$	${N}_{{\rm{ppc}}}$	$L/{R}_{L,e}$
I1	$\infty$	3600	0.28	40	210
I2	$\infty$	7200	0.28	40	210
I3	$\infty$	7200	0.1	40	210
F1	64	7200	0.28	40	640
F2	25	7200	0.28	40	400
F3	25	3600	0.28	40	400
F4	10	7200	0.28	40	250

Note. The parameters are the mass ratio ${m}_{i}/{m}_{e}$, the initial electron magnetization ${\omega }_{c,e}/s$, the ratio between electron temperature and rest mass energy, ${{kT}}_{e}/{m}_{e}{c}^{2}$, the number of particles per cell ${N}_{{\rm{ppc}}}$ (including ions and electrons), and the box size in units of the typical initial electron Larmor radius $L/{R}_{L,e}$ (${R}_{L,e}={v}_{\mathrm{th},e}/{\omega }_{c,e}$, where ${v}_{\mathrm{th},e}^{2}={k}_{{\rm{B}}}{T}_{e}/{m}_{e}$, k_B is Boltzmann's constant, and T_e is the electron temperature). All of the runs initially have ${\beta }_{i}={\beta }_{e}=2$ and $c=0.225{{\rm{\Delta }}}_{x}/{{\rm{\Delta }}}_{t}$, where ${{\rm{\Delta }}}_{t}$ is the simulation time step. The runs shown in this Table share the same electron skin depth $c/{\omega }_{p,e}/{{\rm{\Delta }}}_{x}=5$ (where ${{\rm{\Delta }}}_{x}$ is the separation of grid points), but we used several other simulations to confirm numerical convergence by varying $c/{\omega }_{c,e}/{{\rm{\Delta }}}_{x}$, N_ppc, and $L/{R}_{L,e}$.

Download table as:  ASCIITypeset image

In this section we focus on the nonlinear evolution of the ion and electron pressure anisotropies. As stated above, we will begin by showing simulations where ions have infinite mass.

3.1. Simulations with ${m}_{i}/{m}_{e}=\infty$

Figure 1 shows the early time evolution (until $t\cdot s\approx 1.3$) of the electron pressures perpendicular (${p}_{\perp ,e};$ black solid line) and parallel (${p}_{| | ,e};$ red solid line) to ${\boldsymbol{B}}$ for runs I1 and I2 of Table 1. These runs have ions of infinite mass so that the electrons can only be affected by the electron anisotropy-driven OEF and A/IC instabilities, and their magnetizations are ${\omega }_{c,e}/s=3600$ and 7200, respectively. The black dotted and red dotted lines show the expected evolutions of ${p}_{\perp ,e}$ and ${p}_{| | ,e}$ from the Chew–Goldberger–Low (CGL) or double adiabatic limit (Chew et al. 1956). We see that this adiabatic limit is reasonably well satisfied in the early stage of the simulations (until $t\cdot s\sim 0.7$), regardless of the magnetization ${\omega }_{c,e}/s$. After that, the growth of the electron anisotropy-driven instabilities provides enough pitch-angle scattering to stop the adiabatic evolution of the electron pressure.

Zoom In Zoom Out Reset image size

The presence of the OEF and A/IC instabilities can be seen from Figure 2, which shows the magnetic field fluctuations and plasma density in simulation I1. The upper row in Figure 2 corresponds to $t\cdot s=1$, i.e., after one shear time, while the lower row corresponds to $t\cdot s=2$. At all times the magnetic fluctuations are dominated by their $\delta {B}_{z}$ component, with wavenumbers k ($\equiv | {\boldsymbol{k}}|$, where ${\boldsymbol{k}}$ is the wavevector) satisfying ${{kR}}_{L,e}\sim 0.2$, and with ${\boldsymbol{k}}$ being mainly oblique to the mean direction of ${\boldsymbol{B}}$. The presence of these oblique modes can also be seen from Figure 3(b), which shows the Fourier transform of $\delta {B}_{z}$ at $t\cdot s=1$ as a function of the wavenumbers parallel and perpendicular to $\langle {\boldsymbol{B}}\rangle$ (${k}_{| | }$ and k_⊥, respectively). This suggests that the OEF modes contribute the most to the amplitude of the magnetic fluctuations. However, although smaller in amplitude, quasi-parallel A/IC modes can also be seen especially in the $\delta {B}_{x}$ component (see Figure 2(a)). This is also seen from Figure 3(a), which shows that the Fourier transform of $\delta {B}_{x}$ (as a function of ${k}_{| | }$ and k_⊥) is dominated by quasi-parallel modes.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Figure 4 shows the time evolution of the magnetic energy in $\langle {\boldsymbol{B}}\rangle$, the volume-averaged pressure anisotropy, and the electron magnetic moment for two runs, one with ${\omega }_{c,e}/s=3600$ and the other with ${\omega }_{c,e}/s=7200$ (runs I1 and I2 in Table 1, respectively). Panels (c) and (d) show the volume-averaged pressure anisotropy $-\langle {\rm{\Delta }}{p}_{e}\rangle /\langle {p}_{| | ,e}\rangle$ for these two runs, where ${\rm{\Delta }}{p}_{e}={p}_{\perp ,e}-{p}_{| | ,e}$. For comparison, in both cases we plot the anisotropy threshold that would make the OEF and A/IC modes grow at a rate equal to the shearing rate, s. These thresholds were calculated using the linear Vlasov solver developed by Verscharen et al. (2013). The calculations use mass ratio ${m}_{i}/{m}_{e}=1836$ and assume very cold ions (${\beta }_{i}={10}^{-4}$), which seeks to resemble our simulated ${m}_{i}/{m}_{e}=\infty$ situation where the ions only provide a neutralizing charge.⁷ We see that the OEF and A/IC thresholds are quite similar, although the OEF mode has a slightly lower threshold, especially in the case ${\omega }_{c,e}/s=3600$. This implies that both modes should play some role in regulating the electron anisotropy, with their relative importance depending weakly on the ratio ${\omega }_{c,e}/s$. Also, for both values of ${\omega }_{c,e}/s$ there is a reasonably good agreement between the electron anisotropy obtained from the simulations and the linear OFE and A/IC instability thresholds. This is thus consistent with the electron anisotropy being maintained at the level for the OEF and A/IC modes to grow at a rate close to s.

Zoom In Zoom Out Reset image size

The contribution of the different components of $\delta {\boldsymbol{B}}$ can be seen from Figures 4(a) and (b), which show the magnetic energy along different axes as a function of time, normalized by the average magnetic energy in the simulation, $\langle {B}^{2}\rangle /8\pi$. $\delta {\boldsymbol{B}}$ is decomposed in terms of $\delta {B}_{z}$ (component perpendicular to the simulation plane), $\delta {B}_{{xy},\perp }$ (component parallel to the simulation plane but perpendicular to $\langle {\boldsymbol{B}}\rangle$), and $\delta {B}_{| | }$ (component parallel to $\langle {\boldsymbol{B}}\rangle$). Clearly, $\delta {\boldsymbol{B}}$ is dominated by its z component (as already seen in Figure 2). This shows that, although the OEF and A/IC modes are expected to contribute to limiting the electron anisotropy, their contribution to the magnetic energy in $\delta {\boldsymbol{B}}$ is quite different. Indeed, our linear calculations show that, for the plasma parameters of runs I1 and I2, the OEF modes should satisfy $| \delta {B}_{z}| /(| \delta {B}_{| | }| +| \delta {B}_{{xy},\perp }| )\sim 4$. Thus the fact that in our runs $\delta {B}_{z}^{2}\sim 10\delta {B}_{{xy},\perp }^{2}$ implies that most of the magnetic energy is being contributed by the OEF modes. The quasi-parallel A/IC modes, which are most visible in the $\delta {B}_{x}$ component as can be seen from Figure 2(a), make a significantly smaller contribution to $\delta {\boldsymbol{B}}$. Another sign of the dominance of the OEF modes is the growing and damping phases of $\delta {B}_{z}$ observed in Figures 4(a) and (b), which are likely related to the conversion of the saturated OEF modes into propagating waves that are rapidly damped through scattering with electrons, as has been observed in previous initial-value PIC simulations (e.g. Hellinger et al. 2014). The different contributions to $\delta {\boldsymbol{B}}$ are likely due to the slightly different anisotropy thresholds of the OEF and A/IC instabilities, as well as to the different amplitude at saturation expected for these two modes. Another possible factor is that the growth rate of the A/IC instability is very sensitive to the orientation of the pitch-angle gradients of the distribution function. Therefore, the A/IC instability can relax the distribution faster toward a stable configuration through pitch-angle scattering than the OEF instability.

Finally, Figures 4(e) and (f) show the volume-averaged magnetic moment of the electrons, $\langle {\mu }_{e}\rangle$ ($\equiv \langle {p}_{\perp ,e}/B\rangle ;$ black solid line), for the same runs I1 and I2, respectively. It can be seen that until the onset of the exponential growth of OEF and A/IC ($t\cdot s\approx 0.7$), $\langle {\mu }_{e}\rangle$ is fairly constant, implying the lack of efficient pitch-angle scattering. After that, $\langle {\mu }_{e}\rangle$ tends to increase at a rate close to the shear rate s. This implies the appearance of an effective pitch-angle scattering rate for the electrons, ${\nu }_{\mathrm{eff},e}$, due to their interaction with the OEF and A/IC modes.

In order to help us to understand the way ${\rm{\Delta }}{p}_{j}$ is regulated by the different velocity-space instabilities, we propose a second way to calculate the average magnetic moment of species j:

$$\begin{eqnarray}&&{\mu }_{j,\mathrm{eff}}\equiv \displaystyle \frac{\langle {p}_{\perp ,j}\rangle }{\langle B\rangle }.\end{eqnarray} \tag{ 1 }$$

This definition is useful because there can be cases where ${\mu }_{j,\mathrm{eff}}\ne \langle {\mu }_{j}\rangle$. This is expected when, besides pitch-angle scattering, ${\rm{\Delta }}{p}_{j}$ is partly regulated by relatively large fluctuations in B, which may spatially correlate with ${p}_{\perp ,j}$ in a ${\mu }_{j}$-conserving way. This occurs, for instance, in the presence of large-amplitude mirror modes (Kunz et al. 2014; Riquelme et al. 2015, 2016). Figures 4(e) and (f) show that, after $\langle {\mu }_{e}\rangle$ conservation is broken, ${\mu }_{e,\mathrm{eff}}$ (dotted black) and $\langle {\mu }_{e}\rangle$ (solid black) are almost indistinguishable. This confirms that ${\rm{\Delta }}{p}_{e}$ is regulated by an effective pitch-angle scattering provided by the OEF and A/IC modes, with no significant contribution from fluctuations in B.

3.2. Simulations with Finite ${m}_{i}/{m}_{e}$

In order to study the effect of ions in regulating both the ion and electron pressure anisotropies, we now focus on simulations with finite ratios of ion to electron mass, ${m}_{i}/{m}_{e}$. Since using ${m}_{i}/{m}_{e}\simeq 1836$ is computationally infeasible with our current resources, we have instead tried to ensure that our simulation results do not depend significantly on ${m}_{i}/{m}_{e}$, which is reasonably well achieved for ${m}_{i}/{m}_{e}=25$ and 64.

As an example, Figure 5 shows the three components of $\delta {\boldsymbol{B}}$ for run F1 of Table 1 (${m}_{i}/{m}_{e}=64$ and ${\omega }_{c,e}/s=7200$). The upper and lower rows correspond to $t\cdot s=1$ and $t\cdot s=2$, respectively. At both times, quasi-parallel and oblique modes are present with similar amplitudes, and with wavenumbers satisfying ${{kR}}_{L,i}\sim 0.4$ (where ${R}_{L,i}$ is the ion Larmor radius). While the quasi-parallel modes are apparent in the three components of $\delta {\boldsymbol{B}}$, the oblique modes appear mainly in $\delta {B}_{z}$. This can be seen more clearly from Figure 6, which shows the Fourier transform of $\delta {B}_{x}$ (Figure 6(a)) and $\delta {B}_{z}$ (Figure 6(b)) at $t\cdot s=1$ as a function of ${k}_{| | }$ and k_⊥. The presence of quasi-parallel modes is clear in both panels, while the oblique modes mainly appear in $\delta {B}_{z}$. These features are consistent with the simultaneous presence of both OIF and FM/W modes.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The presence or absence of fluctuations on the scale of the electron Larmor radius is less clear from Figures 5 and 6. We will come back to this question below.

3.2.1. The Role of ${m}_{i}/{m}_{e}$

Figure 7 compares the evolution of the energy in $\delta {\boldsymbol{B}}$, the ion and electron anisotropies, and ${\mu }_{i}$ and ${\mu }_{e}$ for simulations with different mass ratios and electron magnetization. The first and second columns compare simulations with the same electron conditions (${\omega }_{c,e}/s=7200$, ${k}_{{\rm{B}}}{T}_{e}=0.28{m}_{e}{c}^{2}$, and ${\beta }_{e}={\beta }_{i}=2$) but with different mass ratios: ${m}_{i}/{m}_{e}=25$ (run F2) and ${m}_{i}/{m}_{e}=64$ (run F1), respectively.

Zoom In Zoom Out Reset image size

Figures 7(d) and (e) show the volume-averaged electron and ion pressure anisotropies as a function of time (green and black lines, respectively) for runs F2 and F1, respectively. These figures also show the anisotropy threshold for the growth of different instabilities, using a growth rate of s. As dotted lines we show thresholds that assume ${\rm{\Delta }}{p}_{e}={\rm{\Delta }}{p}_{i}$ and correspond to the instabilities: OEF (dotted green), OIF (dotted red), and FM/W (dotted blue). We do not include A/IC thresholds in this case. This is because our linear calculations show that the A/IC modes are subject to cyclotron-resonant ion damping when ${T}_{i}\sim {T}_{e}$, becoming stable in our runs with finite ${m}_{i}/{m}_{e}$ (this is true even for ${m}_{i}/{m}_{e}=1836$).⁸ We see that our simulation results are fairly independent of the mass ratio, and can be summarized as follows:

• 1.  
  The ion and electron anisotropies evolve quite similarly in both cases. (This justifies using instability thresholds that assume ${\rm{\Delta }}{p}_{i}={\rm{\Delta }}{p}_{e}$ for comparison.)
• 2.  
  The obtained electron anisotropy is a factor ∼1.5 smaller than the expected OEF threshold.
• 3.  
  The ion and electron anisotropies are best described by the thresholds of the OIF and FM/W instabilities with ${\rm{\Delta }}{p}_{i}={\rm{\Delta }}{p}_{e}$.
• 4.  
  The OIF and FM/W thresholds with ${\rm{\Delta }}{p}_{i}={\rm{\Delta }}{p}_{e}$ are quite similar, which is consistent with the simultaneous presence of these modes in Figure 5.

The fact that the electron anisotropy is close to the OIF and FM/W thresholds, and a factor ∼1.5 smaller than the expected OEF threshold, shows that the OIF and FM/W modes are the ones with the largest effect on the electron anisotropy. This can be understood as being due to the significant contribution of the electron pressure anisotropy to the growth of OIF and FM/W modes. Indeed, our linear calculations show that the ${\rm{\Delta }}{p}_{i}$ thresholds for the OIF and FM/W instabilities with ${\rm{\Delta }}{p}_{e}=0$ are a factor of ∼2 larger than in the ${\rm{\Delta }}{p}_{i}={\rm{\Delta }}{p}_{e}$ case. This conclusion is also supported by the fact that, for the obtained electron anisotropy, the OEF modes are stable, indicating that the contribution from the OEF modes to the scattering of electrons is not expected to be important.

Finally, we performed similar linear threshold calculations for the case ${m}_{i}/{m}_{e}=1836$ with ${\omega }_{c,e}/\gamma ={10}^{6}$ (where γ represents the growth rate of the different instabilities), which are shown in Figure 8(a). We find that the thresholds of the OIF (dotted red) and FM/W (dotted blue) modes with ${\rm{\Delta }}{p}_{e}={\rm{\Delta }}{p}_{i}$ continue to be similar and smaller than the OEF (dotted green) threshold by a factor ∼1.5 (the dotted green line almost coincides with the solid red line). Also, the OIF and FM/W thresholds with ${\rm{\Delta }}{p}_{e}={\rm{\Delta }}{p}_{i}$ are ∼1.5 times smaller than in the ${\rm{\Delta }}{p}_{e}=0$ case. This implies that, for realistic mass ratios and magnetizations, the dominant instability for the regulation of ion and electron anisotropies should continue to be the OIF and FM/W instabilities.⁹

Zoom In Zoom Out Reset image size

Figures 7(a) and (b) show the magnitude of the volume-averaged magnetic energy in the same three components of $\delta {\boldsymbol{B}}$ plotted in Figure 4: $\delta {B}_{z}$ (solid green), $\delta {B}_{| | }$ (solid black), and $\delta {B}_{{xy},\perp }$ (solid red), for ${m}_{i}/{m}_{e}=25$ and 64, respectively (for comparison the case ${m}_{i}/{m}_{e}=25$ is replicated in Figure 7(b) using dotted lines). For the two mass ratios, the two components perpendicular to $\langle \delta {\boldsymbol{B}}\rangle$ dominate, with $\langle \delta {B}_{z}^{2}\rangle$ being most of the time ∼3 times larger than $\langle \delta {B}_{{xy},\perp }^{2}\rangle$. This result implies that the OIF and FM/W modes are contributing comparable energy to $\delta {\boldsymbol{B}}$, as was also noticeable from Figure 5. Indeed, our linear calculations show that $\delta {B}_{z}^{2}\gg \delta {B}_{{xy},\perp }^{2}$ for the OIF modes, while $\delta {B}_{z}^{2}\sim \delta {B}_{{xy},\perp }^{2}$ for the FM/W modes,¹⁰ which implies that most of the $\delta {B}_{{xy},\perp }$ component is being produced by FM/W modes.

The amplitudes of the OIF and FM/W modes appear to depend on the mass ratio. Although time-dependent, the magnitude of $\delta {B}_{z}^{2}$ in the case ${m}_{i}/{m}_{e}=64$ is on average ∼1.5 times larger than in the case ${m}_{i}/{m}_{e}=25$. Since the magnetizations ${\omega }_{c,i}/s$ of the two runs differ by a factor ∼2.6 ($=64/25$), this is roughly consistent with previous studies of the OIF instability that show that $\delta {B}^{2}$ at saturation should scale as $\delta {B}^{2}/{B}^{2}\propto {(s/{\omega }_{c,i})}^{1/2}$ (Kunz et al. 2014). On the other hand, $\delta {B}_{{xy},\perp }^{2}$ in the case ${m}_{i}/{m}_{e}=64$ is about ∼2 times larger than for ${m}_{i}/{m}_{e}=25$, which is roughly consistent with the expectation for the FM/W modes to have a saturation amplitude that satisfies $\delta {B}^{2}/{B}^{2}\propto s/{\omega }_{c,i}$. This scaling can be obtained from the expected effective ion scattering frequency by resonant waves, ${\nu }_{\mathrm{eff},i}$, which scales as

$$\begin{eqnarray}&&{\nu }_{\mathrm{eff},i}\sim \displaystyle \frac{\delta {B}^{2}}{{B}^{2}}\displaystyle \frac{{\omega }_{c,i}^{2}}{{k}_{| | }{v}_{| | }},\end{eqnarray} \tag{ 2 }$$

where ${k}_{| | }$ and ${v}_{| | }$ are the wavevector component and the particle velocity component parallel to $\langle {\boldsymbol{B}}\rangle$ (Marsch 2006). For the case of the quasi-parallel FM/W waves, we obtained from the simulations that ${k}_{| | }{R}_{L,i}\sim 0.3$, meaning that ${k}_{| | }{v}_{| | }\propto {\omega }_{c,i}$ for most particles and that ${\nu }_{\mathrm{eff},i}\propto (\delta {B}^{2}/{B}^{2}){\omega }_{c,i}$. Thus, since one expects ${\nu }_{\mathrm{eff},i}\sim s\cdot {p}_{| | ,i}/{\rm{\Delta }}{p}_{i}$ at FM/W saturation (see Equation (4) below), this implies that $\delta {B}^{2}/{B}^{2}\propto s/{\omega }_{c,i}$ (considering that the change in ${\rm{\Delta }}{p}_{i}/{p}_{| | ,i}$ between the cases ${m}_{i}/{m}_{e}=25$ and 64 is small).

Finally, in panels 7(g) and (h) we compare $\langle {\mu }_{j}\rangle =\langle {p}_{\perp ,j}/B\rangle$ and ${\mu }_{j,\mathrm{eff}}=\langle {p}_{\perp ,j}\rangle /\langle B\rangle$ for both ions and electrons. We see that the change in ${\mu }_{j,\mathrm{eff}}$ tends to be somewhat larger than the one in $\langle {\mu }_{j}\rangle$ (by $\sim 20 \%$) for the two mass ratios tested. This implies that the combined effect of the OIF and FM/W modes reduces ${p}_{\perp ,j}$ in a way that mainly breaks the adiabatic invariance of ${\mu }_{j}$, with the preservation of ${\mu }_{j}$ due to changes in the field configuration playing a small role.

3.2.2. The Role of ${\omega }_{c,e}/s$

It is also important to understand the role of electron magnetization, ${\omega }_{c,e}/s$, while keeping the same mass ratio. This can be done by looking at the first and third columns in Figure 7, which compares simulations with ${m}_{i}/{m}_{e}=25$ but with electron magnetization ${\omega }_{c,e}/s=7200$ and 3600 (runs F2 and F3 in Table 1, respectively). We see that the two runs reproduce essentially the same results in terms of the evolution of ${\mu }_{j}$ and ${\rm{\Delta }}{p}_{j}$, with the only difference being in the amplitude of the magnetic fluctuations. Apart from some significant time variability, the amplitudes of the components $\delta {B}_{z}$ and $\delta {B}_{{xy},\perp }$ in the run with ${m}_{i}/{m}_{e}=25$ and ${\omega }_{c,e}/s=3600$ are quite similar to the case ${m}_{i}/{m}_{e}=64$, ${\omega }_{c,e}/s=7200$. Since these two runs have a very similar ratio ${\omega }_{c,i}/s$ ($=144$ and 113, respectively), this result is consistent with the dependence of the OIF and FM/W saturated amplitude on ${\omega }_{c,i}/s$ mentioned in Section 3.2.1.

3.2.3. Breaking of ${\mu }_{e}$ Adiabatic Invariance

An important question is whether the ion-scale instabilities alone are capable of explaining the break in the electron magnetic moment shown in Figures 7(g)–(i) (which starts at $t\cdot s\sim 0.7$).

We explore this issue by comparing the power spectra of the fluctuations in our finite ${m}_{i}/{m}_{e}$ runs with the power spectrum produced in the case with ${m}_{i}/{m}_{e}=\infty$. This is done in Figures 9(a) and (b), where the magnetic energy per logarithmic unit of ${k}_{| | }$ and k_⊥ is plotted at $t\cdot s=2$ for runs with ${m}_{i}/{m}_{e}=10,25,64,$ and $\infty$ (runs F4, F2, F1, and I2, respectively; ${k}_{| | }$ and k_⊥ are the magnitudes of the wavevector components parallel and perpendicular to $\langle {\boldsymbol{B}}\rangle$, respectively). The electrons in these simulations have the same conditions (${k}_{{\rm{B}}}{T}_{e}/{m}_{e}{c}^{2}=0.28$, ${\omega }_{c,e}/s=7200$, and initial ${\beta }_{e}=2$), so the different results are due only to the different ${m}_{i}/{m}_{e}$. We see that:

• 1.  
  In the cases with finite mass ratio, as ${m}_{i}/{m}_{e}$ increases, the peaks of the spectra shift to longer wavelengths (in units of ${R}_{L,e}$) in a way consistent with the growth of the ratio ${R}_{L,i}/{R}_{L,e}$.
• 2.  
  In the same way, as ${m}_{i}/{m}_{e}$ increases, there is a growth in the amplitude of the peak of the spectra, which accounts for the expected increase in the amplitude of the OIF and FM/W modes as ${\omega }_{c,i}/s$ decreases.
• 3.  
  The energy of the magnetic fluctuations on scales of ${k}_{\perp }{R}_{L,e}\sim {k}_{| | }{R}_{L,e}\sim 0.2$ is quite similar regardless of the mass ratio used.
• 4.  
  For finite mass ratios, the power spectra develop, via power cascade, a tail that behaves as¹¹ $d\delta {B}^{2}/d{\rm{ln}}({k}_{| | })\propto {k}_{| | }^{-2.8}$ and $d\delta {B}^{2}/d{\rm{ln}}({k}_{\perp })\propto {k}_{\perp }^{-2.8}$.

Zoom In Zoom Out Reset image size

The similar amplitude of the magnetic fluctuations on electron scales suggests that the break in the adiabatic invariance of ${\mu }_{e}$ can in principle be caused by the OIF and FM/W instabilities via a scenario of three steps: (i) ion and electron anisotropies create magnetic fluctuations through the OIF and FM/W instabilities, (ii) part of the energy in the magnetic fluctuations is transferred to electron scales via power cascade, and (iii) electrons are pitch-angle scattered by these fluctuations, producing the break in ${\mu }_{e}$ invariance. We can compare the contributions from the OIF and FM/W modes to the cascading process by looking at the power spectra by components: $\delta {B}_{z}^{2}$, $\delta {B}_{| | }^{2}$, and $\delta {B}_{{xy},\perp }^{2}$. This is done in Figures 9(c) and (d) for the cases ${m}_{i}/{m}_{e}=25$ and 64, respectively. These figures show that $\delta {B}_{z}^{2}$ and $\delta {B}_{{xy},\perp }^{2}$ are comparable within the power-law tails. Since $\delta {B}_{z}^{2}$ and $\delta {B}_{{xy},\perp }^{2}$ are expected to be dominated by the OIF and FM/W modes, respectively, this result suggests that these two modes contribute similar amounts of energy to the power-law tail.

A remaining question is whether the presented scenario is plausible in the more realistic case with ${m}_{i}/{m}_{e}=1836$, where a larger scale separation is expected between ${R}_{L,i}$ and ${R}_{L,e}$. We explore this question using Figure 8(b), where we plot the growth rate γ as a function of ${{kR}}_{L,e}$ for OIF and FM/W modes (green and black lines, respectively). This is done assuming ${\beta }_{i}={\beta }_{e}=10$, ${\rm{\Delta }}{p}_{e}={\rm{\Delta }}{p}_{i}$, and in two regimes: (i) ${m}_{i}/{m}_{e}=64$ with ${\rm{\Delta }}{p}_{j}$ such that the maximum growth rate ${\gamma }_{\max }={\omega }_{c,e}/7200$ (solid lines), and (ii) ${m}_{i}/{m}_{e}=1836$ with ${\rm{\Delta }}{p}_{j}$ such that ${\gamma }_{\max }={\omega }_{c,e}/{10}^{6}$ (dotted lines). Since for each value of γ there are multiple OIF and FM/W wavevectors ${\boldsymbol{k}}$, in Figure 9(b) we chose the maximum and minimum $| {\boldsymbol{k}}| =k$ for each γ. In this way we explore the possibility that the modes with a given γ could have wavelengths close to both the electron and ion Larmor radii. We obtained that:

• 1.  
  In the case with ${m}_{i}/{m}_{e}=64$ and ${\gamma }_{\max }={\omega }_{c,e}/7200$, the fastest growing OIF and FM/W modes have $0.03\lesssim {{kR}}_{L,e}\lesssim 0.06$. This range roughly coincides with the peak of the spectra shown in Figures 9(a) and (b).
• 2.  
  In the more realistic case with ${m}_{i}/{m}_{e}=1836$ and ${\gamma }_{\max }={\omega }_{c,e}/{10}^{6}$, the fastest growing OIF and FM/W modes appear at $0.006\lesssim {{kR}}_{L,e}\lesssim 0.01$. This implies that both modes have wavelengths a factor ∼6 larger than in the case with ${m}_{i}/{m}_{e}=64$ and ${\gamma }_{\max }={\omega }_{c,e}/7200$, which is expected because of the increase by a factor ∼6 in the ratio ${R}_{L,i}/{R}_{L,e}$.

Thus, for ${m}_{i}/{m}_{e}=64$ the scale separation between ${R}_{L,i}$ and ${R}_{L,e}$ allows the generation of magnetic fluctuations at electron scales with enough energy to pitch-angle scatter the electrons. This result relies on the existence of a power cascade with $d\delta {B}^{2}/d{\rm{ln}}(k)\propto {k}^{-2.8}$, which is observed in Figures 9(a) and (b). However, when ${m}_{i}/{m}_{e}=1836$ this scenario seems less likely. Indeed, at electron scales (${{kR}}_{L,e}\sim 0.2$) the cascade of OIF modes can produce an amount of energy $\sim {(0.2/0.01)}^{2.8}\sim 4400$ times smaller than at the ion scales (${{kR}}_{L,i}\sim 0.01;$ see dotted green line in Figure 8(b)). However, at saturation the OIF and OEF modes are expected to satisfy $\delta {B}^{2}/{B}^{2}\propto {(s/{\omega }_{c,j})}^{1/2}$ (j = i and e, respectively). Thus, in order to provide enough ion and electron pitch-angle scattering, the energy at electron scales should be a factor $\sim {({\omega }_{c,e}/{\omega }_{c,i})}^{1/2}\sim {1836}^{1/2}\sim 43$ smaller than at ion scales, which is ∼100 times larger than what can be produced through the cascade of OIF modes.¹²

This difficulty may get ameliorated if the power cascade process were further modified when ${m}_{i}/{m}_{e}=1836$, or if in 3D the spectral index of the cascade power-law tail were different from the one obtained in our 2D simulations. Unfortunately, our current simulations cannot clarify this aspect of the interplay between the electrons and the OIF and FM/W instabilities. It is important to point out, however, that in realistic settings we do not expect the OEF modes to produce the necessary electron-scale fluctuations either, since our linear calculations show that these modes are stable for the electron anisotropy set by the OIF and FM/W instabilities (with ${\rm{\Delta }}{p}_{e}={\rm{\Delta }}{p}_{i}$). Thus it seems likely that the electron anisotropy should continue to be determined by the OIF and FM/W marginal stability condition with ${\rm{\Delta }}{p}_{e}={\rm{\Delta }}{p}_{i}$.

The existence of electron and ion pressure anisotropies in general implies the presence of non-diagonal terms in the pressure tensor, which give rise to an effective viscosity for both species. In our case, particle velocities are nearly gyrotropic with respect to $\langle {\boldsymbol{B}}\rangle$, so the relevant component of the pressure tensor is ${p}_{{xy},j}\propto ({p}_{\perp ,j}-{p}_{\parallel ,j}){B}_{x}{B}_{y}/{B}^{2}$. It can be shown that this pressure component can tap into the velocity shear of the plasma, producing an increase in the internal energy of the particles. In our case, assuming no heat flux, the internal energy density of species j, U_j ($={p}_{\perp ,j}+{p}_{\parallel ,j}/2$), evolves as (Kulsrud 1983; Snyder et al. 1997; Sharma et al. 2007)

$$\begin{eqnarray}&&\displaystyle \frac{\partial {U}_{j}}{\partial t}=-s{\rm{\Delta }}{p}_{j}{B}_{x}{B}_{y}/{B}^{2}=q{\rm{\Delta }}{p}_{j},\end{eqnarray} \tag{ 3 }$$

where $q=-{{sB}}_{x}{B}_{y}/{B}^{2}$ corresponds to the growth rate of B. In the present context, both q and ${\rm{\Delta }}{p}_{j}$ are negative, which implies an increase in U_j. Before the onset of the instabilities, this process is adiabatic and therefore it is a reversible energy gain (in the sense that the increase in U_j would be reversed by reversing the direction of the plasma shear velocity). Indeed, as shown in Figures 1(a) and (b) for the case of electrons, the early evolution of ${p}_{\perp ,j}$ and ${p}_{\parallel ,j}$ follows the CGL or double adiabatic behavior with ${p}_{\perp ,j}\propto B$ and ${p}_{\parallel ,j}\propto 1/{B}^{2}$ (which gives rise to a net growth of ${U}_{j}={p}_{\perp ,j}+{p}_{\parallel ,j}/2$ since the growth in ${p}_{\parallel ,j}$ occurs faster than the decrease in ${p}_{\perp ,j}$). Thus, only after the instabilities start keeping ${\rm{\Delta }}{p}_{j}/{p}_{\parallel ,j}$ in a quasi-stationary regime by breaking ${\mu }_{j}$ invariance (after $t\cdot s\approx 0.7$ in Figures 1(a) and (b)) can the increase in U_j be considered as irreversible heating. Also it is important to point out that the role of the instabilities after $t\cdot s\approx 0.7$ is not the direct heating of the particles by wave–particle interactions. Instead, the role of the instabilities is to limit the pressure anisotropy and therefore to regulate the viscous heating provided by Equation (3).

Figure 10 quantifies the importance of this heating mechanism by showing the volume-averaged ion (solid black) and electron (solid green) heating rates for run F2. We also show the rate of heating by anisotropic viscosity predicted by Equation (3) for ions (dotted black) and electrons (dotted green). For both species there is reasonably good agreement between the particle heating in the simulation and the contribution from the anisotropic stress. This shows that anisotropic viscosity contributes most of the ion and electron heating in collisionless plasmas (with ${T}_{e}\sim {T}_{i}$), both in the case of decreasing magnetic field presented in this work and in the regime of growing field shown in Riquelme et al. (2016).

Zoom In Zoom Out Reset image size

Besides regulating the effective plasma viscosity, pitch-angle scattering by velocity-space instabilities is also expected to limit the mean free path of the particles. To quantify this effect we used the distance D_j(t) traveled along $\langle {\boldsymbol{B}}\rangle$ by $2\times {10}^{4}$ ions and electrons.¹³ We calculated their mean free paths assuming that $\langle {D}_{j}^{2}\rangle ={{tv}}_{\mathrm{th},j}\langle {\lambda }_{j}\rangle$ (where $\langle {\lambda }_{j}\rangle$ represents the average mean free path over species j, while ${v}_{\mathrm{th},j}={({k}_{{\rm{B}}}{T}_{j}/{m}_{j})}^{1/2}$ is their thermal speed), which is valid if the particles move diffusively. This allows us to estimate the average mean free path of species j as $\langle {\lambda }_{j} \rangle =(d \langle {D}_{j}^{2} \rangle /{dt})/{v}_{{\rm{th}},j}$.

Figure 11 shows our estimates of the electron (black) and ion (red) mean free paths for two simulations with ${m}_{i}/{m}_{e}=25$ and 64 (runs F2 and F1, respectively), normalized by ${v}_{\mathrm{th},j}/s$. We see that in both cases there is an initial period when $\langle {\lambda }_{j}\rangle /({v}_{\mathrm{th},j}/s)$ increases as $\sim 2t\cdot s$, which is consistent with the free-streaming behavior $\langle {D}_{j}\rangle \approx {{tv}}_{\mathrm{th},j}$. This is followed by the saturation of $\langle {\lambda }_{j}\rangle$, expected to start after a time of the order of the collision time of particles. The behaviors of $\langle {\lambda }_{e}\rangle$ and $\langle {\lambda }_{i}\rangle$ at this stage are quite similar, with $\langle {\lambda }_{e}\rangle /({v}_{\mathrm{th},e}/s)$ being somewhat smaller than $\langle {\lambda }_{i}\rangle /({v}_{\mathrm{th},i}/s)$ (by a factor ∼1.5).

Zoom In Zoom Out Reset image size

The evolutions of $\langle {\lambda }_{j}\rangle$ seen in Figure 11 are expected to be influenced by the pitch-angle scattering of the particles caused by the velocity-space instabilities. In Riquelme et al. (2016) we estimated this effect assuming an incompressible fluid with no heat flux, where the scattering produced by the instabilities on species j was incorporated using an effective scattering rate ${\nu }_{\mathrm{eff},j}$ (Kulsrud 1983; Snyder et al. 1997). In this model, valid for ${\rm{\Delta }}{p}_{j}/{p}_{| | ,j}\ll 1$, $\langle {\lambda }_{j}\rangle$ behaves as

$$\begin{eqnarray}&&\langle {\lambda }_{j}\rangle \approx \displaystyle \frac{{v}_{\mathrm{th},j}}{{\nu }_{\mathrm{eff},j}}\approx 0.3\displaystyle \frac{{v}_{\mathrm{th},j}}{q}\displaystyle \frac{{\rm{\Delta }}{p}_{j}}{{p}_{| | ,j}}.\end{eqnarray} \tag{ 4 }$$

Equation (4) provides a good approximation to the mean free path of particles in the case of growing magnetic fields (Riquelme et al. 2016), where $\langle {\lambda }_{j}\rangle$ is regulated by the whistler and mirror instabilities. In the cases of decreasing magnetic field, however, this is not the case. Using our measurements of ${\rm{\Delta }}{p}_{j}/{p}_{| | ,j}$ for runs F2 and F1 (see Figures 7(d) and (e)) one obtains that, at $t\cdot s=1.5$, $\langle {\lambda }_{j}\rangle \approx 0.15{v}_{\mathrm{th},j}/s$. However, Figure 11 shows that the average mean free paths of ions and electrons are ∼10 times larger than this simple estimate. On the other hand, considering the evolution of ${\rm{\Delta }}{p}_{j}/{p}_{| | ,j}$ and ${B}^{2}/{B}_{x}{B}_{y}$ (needed to determine q) from $t\cdot s=1.5$ to $t\cdot s=3$, $\langle {\lambda }_{i}\rangle$ should decrease by a factor ∼2, which is essentially what Figure 11 shows. Thus, putting aside the factor ∼10 difference, the scaling of $\langle {\lambda }_{j}\rangle$ on ${\rm{\Delta }}{p}_{j}/{p}_{| | ,j}$ and q presented in Equation (4) is reasonably well reproduced by the simulations with decreasing B.

The factor ∼10 discrepancy is interesting, partly because the behavior of ${\nu }_{\mathrm{eff},j}$ suggested by Equation (4) (${\nu }_{\mathrm{eff},j}\,\approx 3{{qp}}_{| | ,j}/{\rm{\Delta }}{p}_{j}$) is well reproduced in the case of ions in previous hybrid-PIC simulations that studied the saturated state of the firehose and mirror instabilities (Kunz et al. 2014). In that case, however, ${\nu }_{\mathrm{eff},i}$ is not measured using D_i. Instead, they constructed a distribution of the times taken by each ion to change its ${\mu }_{i}$ by a factor of e, and then approximated ${\nu }_{\mathrm{eff},i}^{-1}$ by the width of the distribution. Thus, in order to clarify this discrepancy, we compared two different measurements of the electron mean free path from run F1, one using the variations of ${\mu }_{e}$ (we will refer to this estimate as $\langle {\lambda }_{e,\mu }\rangle$) and the other one using D_e ($\langle {\lambda }_{e}\rangle$). These measurements of $\langle {\lambda }_{e}\rangle$ and $\langle {\lambda }_{e,\mu }\rangle$ are not defined for different times (as in Figure 11), but they correspond to averages between $t\cdot s=1$ and 3.¹⁴ The comparison was made for different groups of electrons, defined by the parameter ${A}_{e}\equiv \langle {v}_{\perp ,e}^{2}\rangle /\langle {v}_{\parallel ,e}^{2}\rangle$, where $\langle \rangle$ represents the average between $t\cdot s=1$ and 3 for each electron. We use A_e as a way to quantify the pitch angle of electrons, which, as we will see below, affects the behaviors of $\langle {\lambda }_{e}\rangle$ and $\langle {\lambda }_{e,\mu }\rangle$ in different ways.

Our results are shown in Figure 12(a). We see that $\langle {\lambda }_{e}\rangle$ and $\langle {\lambda }_{e,\mu }\rangle$ roughly coincide for ${A}_{e}\gtrsim 1$ (pitch angle $\gtrsim 45^\circ$). However, for ${A}_{e}\lesssim 1$ (pitch angle $\lesssim 45^\circ$), $\langle {\lambda }_{e}\rangle$ becomes significantly larger than $\langle {\lambda }_{e,\mu }\rangle$. This is consistent with the fact that, for electrons with small pitch angle, a variation in ${v}_{\perp ,e}^{2}$ of order unity due to scattering (which implies a variation in ${\mu }_{e}$ of order unity) does not imply that they reverse their velocity along ${\boldsymbol{B}}$. Also, when averaged over the entire A_e distribution (shown by the black line in Figure 12(a)), $\langle {\lambda }_{e}\rangle$ is a factor ∼10 larger than $\langle {\lambda }_{e,\mu }\rangle$, showing that using $\langle {\lambda }_{e,\mu }\rangle$ would essentially eliminate the discrepancy between our estimated mean free path and Equation (4).¹⁵

Zoom In Zoom Out Reset image size

Given this difference between $\langle {\lambda }_{e,\mu }\rangle$ and $\langle {\lambda }_{e}\rangle$, it is important to understand why, in the case of growing magnetic field studied by Riquelme et al. (2016), the behavior of $\langle {\lambda }_{e}\rangle$ is well reproduced by Equation (4). Figure 12(b) shows the same quantities as Figure 12(a) but for run MW1 of Riquelme et al. (2016). We see that the trend for $\langle {\lambda }_{e}\rangle$ to grow relative to $\langle {\lambda }_{e,\mu }\rangle$ as A_e decreases is maintained in this case. However, for all values of A_e, the case of growing B tends to have a smaller ratio $\langle {\lambda }_{e}\rangle /\langle {\lambda }_{e,\mu }\rangle$ than the case of decreasing B, making $\langle {\lambda }_{e}\rangle \sim \langle {\lambda }_{e,\mu }\rangle$ if the average over the whole distribution of A_e (black line) is considered.

This difference between the cases of growing and decreasing B appears to be due to the specific effect of the relevant instabilities on the electron velocities. This is suggested by Figure 13(a), which shows the electron velocity distribution $f({v}_{| | ,e},{v}_{\perp ,e}^{(z)})$ for run F1 at $t\cdot s=1.5$ (corresponding to the saturated stage of the FM/W and OIF instabilities), where ${v}_{| | ,e}$ and ${v}_{\perp ,e}^{(z)}$ are respectively the electron velocity parallel to $\langle {\boldsymbol{B}}\rangle$ and parallel to the z axis (which are mutually perpendicular). We see that, for ${v}_{e}\lesssim 0.6c$ (${v}_{e}\equiv {({v}_{| | ,e}^{2}+{v}_{\perp ,e}^{(z)2})}^{1/2}$), $f({v}_{| | ,e},{v}_{\perp ,e}^{(z)})$ is dominated by electrons with rather small pitch angle ($\lesssim 45^\circ$). This suggests that for ${v}_{e}\lesssim 0.6c$, the scattering process occurs in a way that disfavors the diffusion of electrons toward smaller values of $| {v}_{| | ,e}|$, which in turn precludes the reversal of ${v}_{| | ,e}$, contributing to increasing $\langle {\lambda }_{e}\rangle$. For comparison, in Figure 13(b) we show the analogous electron distribution for run MW1 of Riquelme et al. (2016) (growing B) at $t\cdot s=1.5$, where $f({v}_{| | ,e},{v}_{\perp ,e}^{(z)})$ appears more similar to a bi-Maxwellian distribution with ${p}_{\perp ,e}\gt {p}_{| | ,e}$. Notice that the dominance of electrons with small pitch angle for ${v}_{e}\lesssim 0.6c$ seen in the case of decreasing field is similar to the modification to the ion velocity distribution found by previous hybrid-PIC simulations of an expanding box, where the ion scattering is also provided by the FM/W and OIF instabilities (Matteini et al. 2006; Hellinger & Travnicek 2008).

Zoom In Zoom Out Reset image size

In summary, for both growing and decreasing B, the estimate of the electron mean free path provided by Equation (4) is fairly well reproduced by $\langle {\lambda }_{e,\mu }\rangle$ (as was shown by Kunz et al. 2014 in the case of ions). However, since $\langle {\lambda }_{e}\rangle$ is based on the direct calculation of the distance travelled by the electrons along ${\boldsymbol{B}}$, this quantity provides a more meaningful measurement of the electron mean free path for the purpose of quantifying the thermal conductivity of the plasma. In the case of decreasing field, $\langle {\lambda }_{e}\rangle$ is a factor ∼10 larger than $\langle {\lambda }_{e,\mu }\rangle$, most likely due to the specific electron scattering mechanism provided by the FW/W and OIF instabilities. Thus, in this case the electron mean free path is best described by the relation

$$\begin{eqnarray}&&\langle {\lambda }_{j}\rangle \approx 3\displaystyle \frac{{v}_{\mathrm{th},j}}{q}\displaystyle \frac{{\rm{\Delta }}{p}_{j}}{{p}_{| | ,j}},\end{eqnarray} \tag{ 5 }$$

which is valid for $2\lesssim {\beta }_{e}\lesssim 20$, and only differs from Equation (4) by its prefactor ∼3 instead of ∼0.3.

We used PIC plasma simulations to study the nonlinear, saturated stage of various ion and electron velocity-space instabilities relevant for collisionless plasmas. We focused on instabilities driven by pressure anisotropy with ${p}_{\perp ,j}\lt {p}_{| | ,j}$. To capture the nonlinear regime in a self-consistent way, we imposed a shear velocity in the plasma, which decreases the background magnetic field. This drives ${p}_{\perp ,j}\lt {p}_{| | ,j}$ due to the adiabatic invariance of the magnetic moment (the driving timescale is much longer than the gyroperiod of the particles). This, in turn, drives velocity-space instabilities, which inhibit the growth of pressure anisotropy. The relevant instabilities in this regime, as suggested by linear theory, are (i) the purely growing OIF and the resonant FM/W modes, which are mainly driven by the ions, and (ii) the purely growing OEF and the resonant A/IC modes, which are driven by the electrons. The nonlinear state of these instabilities is expected to be influenced by the simultaneous presence of ion and electron anisotropies on the different modes. In order to achieve reasonable scale separation between these modes, we mainly used ${m}_{i}/{m}_{e}=25$ and 64. Our results, valid for the regime $2\lesssim {\beta }_{i}\approx {\beta }_{e}\lesssim 20$, showed no significant difference between these two mass ratios.

We found that the mechanism for regulating the ion and electron anisotropies consists in the growth of OIF and FM/W modes, which affect the ions and electrons equally, rendering ${\rm{\Delta }}{p}_{e}\approx {\rm{\Delta }}{p}_{i}$. The numerically obtained ion and electron anisotropies are well approximated by the linear threshold for the growth of the OIF and FM/W modes with ${\rm{\Delta }}{p}_{e}={\rm{\Delta }}{p}_{i}$ and with growth rate $\sim s$. The electron pressure anisotropy in simulations with ions of infinite mass (where the ions only provide a neutralizing charge) is dominated by the OEF and A/IC modes, giving an anisotropy a factor ∼2 larger than in the cases with finite ${m}_{i}/{m}_{e}$. We attribute this result to the rather strong destabilizing effect of the electron pressure anisotropy, ${\rm{\Delta }}{p}_{e}$, on the OIF and FM/W modes (as already suggested by previous linear dispersion analyses, see Michno et al. 2014; Maneva et al. 2016), which in turn maintains ${\rm{\Delta }}{p}_{e}$ at a value significantly lower than the one necessary to make the OEF and A/IC modes grow at a rate $\sim s$.

Although the amplitude of the OIF and FM/W modes depends on the ratio ${\omega }_{c,i}/s$, the values of the parameters ${m}_{i}/{m}_{e}$ and ${\omega }_{c,e}/s$ used in the simulations do not affect our conclusions. Also, based on our linear Vlasov calculations (Verscharen et al. 2013), we infer that the presented scenario should hold in the highly magnetized (${\omega }_{c,i}/s\gg 1$) case with m_i/m_e = 1836 relevant for real astrophysical plasmas. However, an important point that our simulations could not completely clarify (due to the lack of sufficient separation of ion and electron scales) is the mechanism by which the electrons would be pitch-angle scattered in the saturated stage of the OIF and FM/W instabilities in the case of ${m}_{i}/{m}_{e}=1836$. Answering this question requires using significantly larger mass ratios and magnetizations (and possibly 3D runs), so we defer this aspect of the study to a future work.

We have also used our simulations to verify the expected viscous heating of particles, which is described in Equation (3), and arises due to pressure anisotropies tapping into the free energy in the shear motion of the plasma. Figure 10 shows good agreement between the heating of the particles in our simulations and the expectation from Equation (3). This result is valid for decreasing magnetic fields, as shown here, and for growing fields, as shown by Riquelme et al. (2016).

With the intention of quantifying the thermal conductivity in these plasmas, we have also computed the mean free path of species j, $\langle {\lambda }_{j}\rangle$, during the nonlinear stage of the OIF and FM/W instabilities. The average mean free path of both ions and electrons is reasonably well described by Equation (5). The scaling factors in this equation are the same as in Equation (4), which is based on a model where the mean free path of species j is determined by an effective scattering rate ${\nu }_{\mathrm{eff},j}$ that sets the rate at which the invariance of ${\mu }_{j}$ is broken (Kulsrud 1983; Snyder et al. 1997). However, the prefactor in Equation (5) is ∼10 times larger than in the case of Equation (4).

We explained this discrepancy for the case of the electrons by noticing that the variation of ${\mu }_{e}$ provides a good estimate of $\langle {\lambda }_{e}\rangle$ only for electrons with relatively large pitch angle ($\gtrsim 45^\circ$). For electrons with pitch angles smaller than $\sim 45^\circ$, $\langle {\lambda }_{e}\rangle$ can become a factor ∼10 larger. This is in contrast with the case of growing B (Riquelme et al. 2016), where Equation (4) provides a reasonably good estimate for $\langle {\lambda }_{j}\rangle$, despite the fact that $\langle {\lambda }_{j}\rangle$ also grows as the pitch angle decreases. This difference is likely due to the specific electron scattering mechanism provided by the FW/W and OIF instabilities, which tends to preclude the reversal of ${v}_{| | ,e}$, contributing to the increase in $\langle {\lambda }_{e}\rangle$.

The results shown in this work, as well as those presented in Riquelme et al. (2016) for the case of growing B, are relevant for quantifying the viscous heating and thermal conductivity in various low-collisionality astrophysical plasmas, including low-luminosity accretion flows around compact objects, the ICM, and the heliosphere.

This work was supported by NSF grants AST 13-33612 and 1715054, a Simons Investigator Award to E.Q. from the Simons Foundation and the David and Lucile Packard Foundation. D.V. also acknowledges support from NSF/SHINE grant AGS-1460190, NSF/SHINE grant AGS-1622498, NASA grant NNX16AG81G, and a UK STFC Ernest Rutherford Fellowship (ST/P003826/1). We are also grateful to the UC Berkeley–Chile Fund for support for collaborative trips that enabled this work. This work used the Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by National Science Foundation grant number ACI-1053575.