{'fraction_non_alphanumeric': 0.0404356031308975, 'fraction_numerical': 0.018933886087306253, 'mean_word_length': 4.432521008403361, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 12, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
A Hamiltonian gyrofluid model based on a quasi-static closure

A Hamiltonian six-field gyrofluid model is constructed, based on closure relations derived from the so-called"quasi-static"gyrokinetic linear theory where the fields are assumed to propagate with a parallel phase velocity much smaller than the parallel particle thermal velocities. The main properties of this model, primarily aimed at exploring basic phenomena of interest for space plasmas such as the solar wind, are its ability to provide a reasonable agreement with kinetic theory for linear low-frequency modes, and at the same time to ensure a Hamiltonian structure in the absence of explicit dissipation. The model accounts for equilibrium temperature anisotropy, ion and electron finite Larmor radius corrections, electron inertia, magnetic fluctuations along the direction of a strong guide field, and parallel Landau damping. Remarkably, the quasi-static closure leads to exact expressions for the nonlinear terms involving gyroaveraged electromagnetic perturbations. One of the consequences is that a rather natural identification of the Hamiltonian structure of the model becomes possible when Landau damping is neglected. A slight variant of the model consists of a four-field Hamiltonian reduction of the original six-field model, which we use for the subsequent linear analysis. In the latter, the dispersion relations of kinetic Alfv\'en waves and the firehose instability are shown to be correctly reproduced, relatively far in the sub-ion range (depending on the plasma parameters), while the spectral range where the slow-wave dispersion relation and the field-swelling instabilities are precisely described is less extended. This loss of accuracy originates from the breaking of the condition of small phase velocity, relative to the parallel thermal velocity of the electrons (for kinetic Alfv\'en waves and firehose instability) or of the ions (in the case of the field-swelling instabilities).


Introduction
Modeling the dynamics of collisionless (or weakly collisional) plasmas at scales comparable to or smaller than the ion Larmor radius is an important issue both for laboratory and astrophysical plasmas. At the level of a kinetic description, a valuable tool is given by the gyrokinetic theory which provides a reduction of the Vlasov-Maxwell (VM) equations by focusing on phenomena with a characteristic time scale large compared with the ion gyro-period. This approach, which eliminates the dependency on the gyration angle, typically adopts, as dynamical variables, the distribution functions of the gyrocenters rather than those of the particles. Within the gyrokinetic framework, a subset of models consists of the so-called δf -gyrokinetic models, which assume the gyrocenter distribution functions of the various particle species, to be close to those of an equilibrium state. In spite of this reduction, numerical simulations of three-dimensional gyrokinetic equations in a turbulent regime (even in the δf framework) require huge computational resources, which justifies the development of simpler (although less complete) descriptions based on gyrofluid equations governing the evolution of a finite number of moments of the gyrocenter distribution functions. The relation between gyrocenter and particle moments is well-defined and can usually be computed perturbatively. As in the case of the fluid hierarchy derived from the VM equations, a gyrofluid hierarchy of equations needs to be closed, in order to obtain a gyrofluid model with a finite number of dynamical variables. An important condition to prescribe at the level of closure assumptions is the preservation, in the absence of dissipation, of the Hamiltonian character of the parent gyrokinetic equations. This guarantees that in the reduction from a gyrokinetic to a gyrofluid system, not only no uncontrolled dissipation of the total energy is introduced but also that further invariants (Casimir invariants) of the system exist, and that the dynamics takes place on hyper-surfaces in phase space where the values of these invariants are constant. Another constraint is the consistency with the linear gyrokinetic theory. In particular, this requires retaining the influence of resonant effects such as Landau damping. A closure accounting for Landau damping in the δf approach typically introduces dissipation and thus prevents the model from being Hamiltonian. Such form of dissipation is, however, voluntarily added and the main requirement is that the model possesses a Hamiltonian structure when the dissipative terms are removed.
In this spirit, the main goal of the present paper is to construct a gyrofluid model possessing the above mentioned properties, and primarily addressed to study phenomena relevant for collisionless space plasmas. Motivated by measurements of sub-proton fluctuations in the solar wind (see Sahraoui et al. (2010) and Alexandrova et al. (2013) or Bruno & Carbone (2013) for reviews), reduced fluid models have already been derived and numerically integrated to explore the dynamics of space plasmas. Kinetic Alfvén wave (KAW) turbulence was for example addressed in Boldyrev & Perez (2012). A more general Hamiltonian reduced gyrofluid model  which, in the appropriate asymptotic limit yields the model of Boldyrev & Perez (2012), was recently developed to simultaneously capture the three regimes of: dispersive Alfvén waves (at scales larger than the sonic and/or ion Larmor radius), of KAWs at sub-ion scales, and also of inertial kinetic Alfvén waves (at scales comparable to the electron inertial length (Chen & Boldyrev 2017;Passot et al. 2017)). It has been used to derive weak turbulence kinetic equations (Passot & Sulem 2019) and the properties of imbalanced KAW turbulence in the framework of a reduction to nonlinear diffusion equations in spectral space (Miloshevich et al. 2020). Nevertheless, to the best of our knowledge, existing reduced gyrofluid models can capture parallel magnetic field fluctuations, electron inertia and ion finite Larmor radius effects, but (with the exception of the recent models of Tassi (2019), which will be discussed in Sec. 3.2.2) do not take into account a possible temperature anisotropy of the equilibrium state. It however turns out that collisionless space plasmas, such as the solar wind, often exhibit anisotropic distribution functions (Marsch 2012) that can result from various heating effects such as Landau damping (Chen et al. 2019) or stochastic heating (Bourouaine & Chandran 2013;Hoppock et al. 2018) or from mechanical effects such as the action of a shear flow (De Camillis et al. 2016). Proton temperature anisotropies play an important role in the solar wind (Hellinger et al. 2006), and possibly even more so at closer distance from the Sun (Huang et al. 2020). Temperature anisotropies are usually constrained by the micro-instabilities they trigger (e.g. mirror and firehose instabilities), both for ions (Bale et al. 2009) and electrons (Štverák et al. 2008). The range of these accessible temperature anisotropies increases as the beta parameter (ratio of thermal to magnetic presure) decreases. In addition, temperature anisotropies are known to affect the development of the tearing instability (Shi et al. 1987), and thus the stability of current sheets (Matteini et al. 2013), which plays a major role in the turbulence evolution (Franci et al. 2017). Another important feature to account for is the coupling to ion acoustic waves, which permits the development of the parametric decay of Alfvén waves at the MHD scales at small beta (Del Zanna et al. 2001), an effect that contributes to the generation of counter propagating waves required for the development of a turbulent dynamics at these scales. Interestingly, this parametric instability, proposed by Bowen et al. (2018) for the generation of the solar wind compressive fluctuations, turns out to occur in a wider range of beta parameters in the presence of temperature anisotropy (Tenerani et al. 2017). A useful feature of a new gyrofluid model for space plasma studies would thus be its capability to account for equilibrium temperature anisotropies and the coupling to ion acoustic waves.
As starting point for the derivation of the model, we choose the δf gyrokinetic equations presented in Kunz et al. (2015) where, for the sake of simplicity, we assume an electron-proton plasma with an equilibrium state described by bi-Maxwellian distribution functions with no mean drift velocity. Such a system fully satisfies our requirements. Indeed, it is a δf gyrokinetic model mainly conceived for pressure-anisotropic astrophysical plasmas and, as such, it specifically accounts for equilibrium temperature anisotropy (unlike most of gyrokinetic models which consider a generic equilibrium distribution function or specialize to the case of a Maxwellian equilibrium). Also, it accounts for parallel magnetic perturbations and it has been shown to possess (at least in the limit of interest for our derivation) a Hamiltonian structure (Tassi 2019). As far as the number of moments to be retained in the gyrofluid model is concerned, the inclusion of Landau damping requires to retain at least the first three moments for each particle species. This is why, under the assumption of a two-species plasma, we opt for the derivation of a six-field gyrofluid model evolving three moments, including parallel temperature fluctuations, for each species. Nevertheless, a four-field Hamiltonian reduced version will also be presented and applied. Another novelty with respect to already existing Hamiltonian gyrofluid models, is the adoption of a closure relation, referred to as quasistatic, derived from linear gyrokinetic theory in the limit of slowly-evolving fields. More precisely, according to such closure, all gyrofluid moments that are not determined by gyrofluid evolution equations, are fixed according to their expression obtained from the gyrokinetic linear theory in the limit |ω{pk z v th s q| ! 1, where ω is the frequency of a mode, k z the component of its wave vector along the direction of the guide field, and v th s is the thermal speed of the species s, associated with the equilibrium temperature along the direction of the guide field. As such, this closure is suitable for fields slowlyevolving (i.e. quasi-static) with respect to particles travelling at the parallel equilibrium thermal speed v th s . The derivation of this closure relation will be presented in Appendix A (see in particular Eq. (A 21) to find the expressions for the various gyrofluid moments according to the quasi-static closure). We anticipate, however, two remarkable properties that this closure possesses. The first one is that the quasi-static closure relation turns out to be compatible with a Hamiltonian structure. The second one is that it allows for exact expressions, in terms of canonical Poisson brackets, for all the nonlinear terms in the gyrofluid equations, and in particular for those involving only gyroaveraged electromagnetic fields or potentials. This is not the case, to the best of our knowledge, for the previously derived reduced gyrofluid models.
The model assumes the presence of a strong magnetic guide field and evolves, for both electrons and ions, gyrocenter density fluctuations as well as velocity and temperature fluctuations referred to the direction parallel to the guide field. A dissipative variant of the model accounting for parallel Landau damping is then formulated through a Landau-fluid modeling of the parallel heat fluxes (Hammett & Perkins 1990;Hammett et al. 1992). A four-field reduction assuming isothermality is presented as well. The resulting gyrofluid model retains ion and electron finite Larmor radius (FLR) corrections, electron inertia and parallel magnetic fluctuations, enables anisotropic equilibrium temperatures and does not prescribe special restrictions on the ion or electron β e,i parameters, where β e,i indicates, for each species (e for electrons and i for ions), the ratio between the equilibrium kinetic pressure and the magnetic pressure exerted by the guide field. In this respect, this gyrofluid model differs from most of the presently available gyrofluid models, which require β e,i ! 1.
In addition to the derivation of the model, we also carry out a detailed analysis of its linearized version. Predictions of the six-field model extended with Landau damping and of its four-field Hamiltonian reduction are compared at the linear level with those of the parent gyrokinetic model, by considering the dispersion relations of KAWs and slow waves (SWs), and also by analyzing the firehose or the field-swelling (Basu & Coppi 1982, 1984 instabilities. This latter instability, which requires the electron temperature transverse to the magnetic field larger than the longitudinal one, leads to a local increase in the transverse pressure, which tends to make the magnetic field "swell" further locally, an effect which can be important in producing magnetic reconnection. Its properties are rather subtle, especially when considering its effect on fast modes when the disturbance propagation is nearly perpendicular to the ambient magnetic field. We thus include an appendix summarizing the results on this instability that are relevant for the discussion of the present model. As will be remarked in Sec. 4, for the non-dissipative case, the choice of the four-field model instead of the six-field model, for the comparison with the linear gyrokinetic theory, is due to the fact that a closure at an even order as in the four-field model, provides a better agreement with the Padé approximant chosen for the electron response function in the quasi-static limit.
The investigation of the linear dispersion relations including temperature anisotropies is a first application of the model, mainly devoted to test its capability of reproducing results of the linear gyrokinetic theory in the appropriate regimes. We mention that the derivation of the model was motivated also by further physical applications of relevance to space plasmas and which will be part of subsequent works. These include the investigation of tearing instability in the presence of temperature anisotropy in a strong guide field regime, the influence of electron FLR effects and temperature anisotropy on inertial reconnection, or the effect of the coupling of Alfvén and compressible modes on the turbulence development.
The paper is organized as follows. In Sec. 2 the gyrokinetic parent model is reviewed. In Sec. 3 the six-field gyrofluid model is introduced and its Hamiltonian structure is presented. The two variants, corresponding to the Landau gyrofluid extension and to the Hamiltonian four-field reduction are also described. Sec. 4 is devoted to the comparison of the linearized versions of the two variants of the model with other linear theories. We conclude in Sec. 5 where we also mention the interest of the present gyrofluid model for space plasma applications. At the end of the paper three Appendices are provided, presenting the derivation of the quasi-static closure relations from gyrokinetic theory, the derivation of the six-field gyrofluid model equations and a discussion of field-swelling instabilities, respectively.

The gyrokinetic parent model
In order to derive a gyrofluid model only based on a quasi-static closure assumption, we consider as starting point the following set of gyrokinetics equations, which corresponds to the system provided by Eqs. (C58), (C60), (C66)-(C68) of Kunz et al. (2015) when collisions and equilibrium velocities are neglected and bi-Maxwellian distribution functions are chosen as equilibrium distribution functions for all the particle species. For simplicity, we specialize to the case of a plasma consisting of two species: electrons and one species of single ionized particles. The equations of the resulting gyrokinetic model are given by: The index s P te, iu adopted above indicates the particle species, so that quantities labelled with s refer to the electron species when s " e and to the ion species when s " i.
In the system (2.1)-(2.4), the function r g s is defined by where r f s is the perturbation of the gyrocenter distribution function for particles of species s, q s is the charge of these particles (so that q e "´e and q i " e, with e indicating the proton charge) and m s their mass. Furthermore, c denotes the speed of light. The bi-Maxwellian equilibrium distribution function is given by where n 0 is the uniform and constant equilibrium density, T 0 s and T 0Ks are respectively the equilibrium temperatures of the s-th particle species parallel and in a plane perpendicular to an equilibrium magnetic guide field of amplitude B 0 , directed along the z direction of a Cartesian coordinate frame tx, y, zu. We suppose that the spatial domain of the system corresponds to a box D " tpx, y, zq :´L x ď x ď L x ,´L y ď y ď L y ,´L z ď z ď L z u, with L x , L y and L z positive constants. All quantities of the system which depend on the spatial variables x, y and z are supposed to satisfy periodic boundary conditions on the domain D, so that they can be expanded in Fourier series. We indicated with v } P R the velocity coordinate parallel to the guide field and with µ s " m s v 2 K {p2B 0 q P r0,`8q the magnetic moment of the particle of species s in the unperturbed guide field, where v K corresponds to the velocity coordinate perpendicular to the guide field. We assume that all functions depending on v } decay to zero as v } Ñ˘8. Functions depending on µ s are assumed to tend to zero as µ s Ñ`8 and to be bounded at µ s " 0. The coordinate t P r0,`8q refers to time. The expressions J 0s and J 1s are related to the standard gyroaverage operators for the species s in Fourier space. Their definition can be introduced explicitly in the following way: adopting the notation x to indicate a point of coordinates px, y, zq P D and, similarly, k to indicate a point pk x , k y , k z q P D, where D is the lattice defined by D " tp2πl{p2L x q, 2πm{p2L y q, 2πn{p2L z qq : pl, m, nq P Z 3 u, we can consider a function f : Dˆr0,`8q Ñ R, periodic over D, so that it admits the Fourier representation f px, tq " ř kPD f k ptq exppik¨xq. The action of the operators J 0s and J 1s on the function f is defined by where J 0 and J 1 indicate the zeroth and first order J Bessel functions, respectively, a s " k K ρ Ks is the perpendicular Larmor radius associated with the species s, with k K and ρ Ks corresponding to the perpendicular wave number and the Larmor radius of the particle of species s. The former is defined as k K " b k 2 x`k 2 y with k x " 2πl{p2L x q, k y " 2πm{p2L y q for l, m P Z, while the latter is given by ρ Ks " v K {ω cs , where ω cs " eB 0 {pm s cq is the cyclotron frequency referred to the guide field and related to the particle of species s.
The leading order expression (up to second order terms in the perturbations) for the magnetic field is given by Bpx, y, z, tq " ∇ r A px, y, z, tqˆẑ`p r B px, y, z, tq`B 0 qẑ, (2.9) whereẑ is the unit vector along the z direction, r A (referred to as magnetic flux function) corresponds to the z component of the magnetic vector potential and r B is the perturbation of the magnetic guide field, also referred to as parallel magnetic perturbation or parallel magnetic fluctuations. We remark that the guide field is assumed to be spatially homogeneous. This assumption is valid for a local description of space plasmas such as the solar wind, where the background magnetic field varies on scales so large that, in the local description, it can be assumed to be homogeneous. The situation would be different, for instance, in the case of tokamak plasmas. Indeed, gyrofluid models more oriented towards tokamak applications (as, for instance those of Snyder & Hammett (2001); Madsen (2013); Brizard (1992); Scott (2010); Waelbroeck & Tassi (2012);Keramidas Charidakos et al. (2015)) take into account background magnetic inhomogeneities. The set of electromagnetic quantities involved in the system is completed by the electrostatic potential r φ " r φpx, y, z, tq. In Eqs. (2.2)-(2.4) we adopted the symbol dW s " p2πB 0 {m s qdµ s dv } to indicate the volume element in space velocity.
The parameters Θ s and β Ks are defined by (2.10) and measure, for each species s, the equilibrium temperature anisotropy and the ratio between equilibrium kinetic and magnetic pressure, respectively. Equation (2.1) is the gyrokinetic equation related to the species s, whereas Eqs. (2.2)-(2.4) relate the gyrocenter distribution functions to electromagnetic quantities in the non-relativistic limit. In particular, Eq. (2.2) corresponds to the quasi-neutrality relation, whereas Eqs. (2.3) and (2.4) descend form Ampère's law projected along directions parallel and perpendicular to the guide field, respectively.
The gyrokinetic model (2.1)-(2.4) is valid for small perturbation of the equilibrium distribution function ( δf approximation) and weak variations of the fields along the direction of the guide field, the equilibrium temperature anisotropy Θ s and the parameter β Ks of all the species, being kept finite in this asymptotics. Further details about the derivation the regime of validity of the model and its derivation can be found in Kunz et al. (2015). Its Hamiltonian structure, on the other hand, is presented in Tassi (2019).

The gyrofluid model
We define the following gyrofluid moments: For each particle species, the fields r N s and r U s represent the fluctuations of the gyrocenter densities and parallel fluid velocities, respectively. On the other hand, r T s and r T Ks correspond to the fluctuations of the gyrofluid temperatures defined with respect to the parallel and perpendicular gyrocenter velocities, respectively, whereas r Q s indicates the gyrocenter parallel heat flux fluctuations. In defining the parallel temperature and heat flux fluctuations, we introduced the constant v th s " b T 0 s {m s , corresponding to the parallel thermal velocity associated with the species s.
We intend to derive a gyrofluid model by taking moments of the gyrokinetic equation (2.1) and by imposing a closure relation derived from the quasi-static linear theory. In particular, the gyrofluid model, accounting for equilibrium temperature anisotropy, should be able, in the limit of vanishing finite Larmor radius effects, to reproduce the field-swelling instability criterion of Basu & Coppi (1984). We restrict to the evolution of the first three moments referring to the parallel direction. Therefore, the resulting gyrofluid model should evolve the following six fields: r N e , r N i , r U e , r U i , r T e and r T i . Also, we show that the model conserves the total energy and, moreover, that it possesses a noncanonical Hamiltonian structure, as is the case for the parent gyrokinetic model (Tassi 2019). We refer to such model as to the 6-field gyrofluid (GF6) model.
For the sake of the comparison, carried out in Sec. 4, of the linear gyrofluid theory with the linear gyrokinetic theory we also consider an extension of GF6 accounting for a Landau-fluid closure, analogously to that discussed in Hammett & Perkins (1990); Hammett et al. (1992); Snyder et al. (1997); Passot & Sulem (2007, 2015; Tassi et al. (2018). This variant of the model, which we denote as GF6L, differs from GF6 for the expression of the parallel heat flux fluctuations r Q s . Therefore, in order to avoid some redundancy in the exposition, in the following we present the model equations leaving r Q s unspecified and we will subsequently indicate the corresponding expressions for r Q s leading to GF6 and GF6L, respectively. The closure leading to GF6L, in particular, will be given in Sec. 3.2. A further variant of the model, denoted as GF4, will also be considered in Sec. 4. This model evolves only the four fields r N e , r N i , r U e and r U i and represents a minimal Hamiltonian model, derived from the quasi-static closure, capable to reproduce the fieldswelling instability criterion. This variant of the model will be introduced in Sec. 3.2.
The six-field system, both in the Hamiltonian (GF6) and dissipative (GF6L) versions, can be written, in a dimensionless form, as In Eqs. (3.3)-(3.8), we adopted the following normalized quantities where the quantities with overbars in Eq. (3.12) are the dimensional spatial and time coordinates. In Eqs. (3.9)-(3.12) the quantities: were introduced, which indicate the sound speed and the sonic Larmor radius, respectively, based on the perpendicular equilibrium temperature. The parameter for s " e, i, on the other hand, measures the ratio of the equilibrium perpendicular temperatures. In Eqs. (3.3)-(3.8) we introduced the operator b s "´∇ 2 K ρ 2 thKs , with ∇ 2 K denoting the Laplacian relatively to the transverse variables and indicating the perpendicular thermal Larmor radius associated with the species s.
The canonical bracket r , s, on the other hand, is defined as rf, gs " B x f B y g´B y f B x g, for two functions f and g.
The gyroaverage operators G 10s , G 20s , Γ 0s , Γ 1s in Eqs. (3.3)-(3.8) can in turn be expressed in terms of the operators which have to be intended as Fourier multipliers whose symbols are obtained by replacing b s by k 2 K ρ 2 thKs . We make this statement more precise, as an example, in the case of the operator G 10i , referred to the ion species. The expression G 10i f is defined by G 10i f px, tq " ř kPD G 10 pb i qf k ptq exppik¨xq " for a function f periodic in space. Analogous expressions are valid for the other gyroaverage operators. In Eq. (3.17), in particular, the symbols I 0 and I 1 indicate the modified Bessel functions I of order zero and one, respectively.
The expressions for the operators G 10s and G 20s given in Eq. (3.16) correspond to those present in Brizard (1992) and follow from assuming that the perturbation of the distribution function can be written as (3.18) where H m and L n indicate the Hermite and Laguerre polynomials, respectively, of order m and n, with m and n non-negative integers. The functions f mn s are coefficients of the expansion and are related to the moments of r f s , with respect to Hermite polynomials in v } {v th s and Laguerre polynomials in µ s B 0 {T 0Ks . Indeed, from the orthogonality properties of Hermite and Laguerre polynomials, the following relation holds: : Note that, according to a customary notation, in the symbols cs K and ρs, the subscript s is to indicate sonic quantities and not the particle or gyrocenter species.
It can be useful, in particular, to write explicitly the following relations between the lowest order normalized moments and the quantities f mns : where, in Eq. (3.21), we introduced the normalized gyrocenter perpendicular temperature fluctuations T Ks " r T Ks {T 0Ks . As above anticipated, in the system (3.3)-(3.8) we temporarily left unspecified the expression for the parallel heat flux fluctuations Q }s . In order to obtain the model GF6, the infinite hierarchy of gyrofluid equations following from the parent gyrokinetic model (2.1)-(2.4) is closed by imposing relations obtained by computing the gyrocenter moments other than N s , U s and T }s from a linearization of the parent gyrokinetic system about a homogeneous equilibrium, in the quasi-static limit. In the case of Q s , this leads (consider Eqs. (3.22) and (A 21)) to The model GF6 is thus obtained by inserting the relation (3.23) into the system (3.3)-(3.8). Details on the derivation of the closure relations originated from the quasi-static assumption can be found in Appendix A. A remarkable property of this closure is that it leads to the annihilation of all the contribution of the higher order moments in the gyrofluid equations. We find it useful to provide also a reformulation of the evolution equations (3.3)-(3.5), which should help putting in evidence the physical nature of the terms contributing to the evolution of the various fields. Eqs. (3.3)-(3.5) can indeed be rewritten as where the parallel gradient operator ∇ s is defined, for each species s, by for a function f . From the formulation (3.24)-(3.26) it emerges that the gyrocenter density, parallel momentum and parallel temperature fluctuations are all advected, in the perpendicular plane, by the incompressible velocity field (3.28) Such velocity field includes a first contribution, associated with G 10s φ, which corresponds to the usual EˆB drift (based on the gyroaveraged electrostatic potential), ubiquitous in low-β gyrofluid models such as those discussed in Waelbroeck & Tassi (2012), Snyder & Hammett (2001), Keramidas Charidakos et al. (2015) and Waelbroeck et al. (2009). When higher β values are allowed, however, the perpendicular advection acquires a further contribution due to the parallel magnetic perturbations, as it transpires from Eq.
(3.28). We remark that, similarly to the EˆB contribution, also the latter contribution does not vanish in the limit b s Ñ 0 of negligible FLR corrections. This is a consequence of the fact that, in the presence of parallel magnetic perturbations, gyrocenter density fluctuations differ from particle density fluctuations, even in the absence of FLR corrections (Brizard 1992). From the continuity equation (3.24), we see that the evolution of gyrocenter density fluctuations has also a source in the last term on the left-hand side, which is due to the compressibility of the gyrocenter mean velocity along the direction of the magnetic field.
From Eq. (3.25) we see that parallel momentum, in addition to be advected by u Ks , evolves due to the term´sgnpq s qE s , where E s is defined by (3.29) Such term represents the force exerted by the gyroaveraged electric field, along the direction of the magnetic field. A further source for the parallel momentum is due to the term ∇ s P s , where This terms is associated with the parallel component of the divergence of an anisotropic pressure tensor. The parallel temperature equation (3.26) has, on its left-hand side, the same structure of the continuity equation (3.24). We just remark the presence of the coefficient 2 multiplying ∇ s U s . This coefficient, of course, follows directly from taking the second order moment, in Hermite polynomials for the normalized parallel velocity, of the gyrokinetic equation (2.1), as discussed in Appendix B. However, as pointed out by Keramidas Charidakos et al. (2015), in the presence of background magnetic curvature such coefficient, in general, has to be adjusted in order to obtain a Hamiltonian structure. Finally, we consider the term on the right-hand side of Eq. (3.26), associated with the parallel heat flux. If the expression for Q s is chosen according to the quasi-static closure, i.e. imposing the relation (3.23), this term vanishes and the system, as will be shown in Sec. 3.1, is Hamiltonian. On the other hand, if the Landau fluid closure (3.35) is chosen, this term acts as a sink and the system is not energy-conserving.

Hamiltonian structure of GF6
The quasi-static closure (3.23) leading to GF6, allows the resulting model to be cast in Hamiltonian form. In particular, it can be verified by direct computation that, when the electromagnetic fluctuations φ, A and B can be expressed in terms of the dynamical variables N e , N i , M e and M i (where we introduced the short-hand notation M s " pm s {m i qU s`s gnpq s qG 10s A to indicate the parallel canonical momenta), by making use of the relations (3.6)-(3.8), the evolution equations (3.3)-(3.5) complemented by Eq.
and t , u is a noncanonical Poisson bracket given by for two functionals F and G. For details about the noncanonical Hamiltonian formulation of fluid models one can refer, for instance, to Morrison (1998). In Eq. (3.33) the subscripts on the functionals indicate functional derivatives. In order to verify the formulation (3.31) it is convenient to remark that, from Eq. (3.32), one obtains (3.34) In order to derive the relations (3.34) we made use of the formal symmetry of the operators G 10s and G 20s , i.e.
ş d 3 x f G 10s g " ş d 3 x gG 10s f and ş d 3 x f G 20s g " ş d 3 x gG 20s f for two functions f and g, as well as of the formal symmetry of the linear operators in terms of which one can express φ, B and A in terms of N s and M s through Eqs. (3.6)-(3.8). We did not provide the explicit expression for such operators which can, however, be obtained considering the representation in Fourier series of the fields involved, following the procedure discussed in Tassi (2019).
The Hamiltonian functional (3.32) is a conserved quantity for the dynamics and corresponds to the total energy. In Eq. (3.33), the sum of all the terms with s " e is a Poisson bracket in its own right. Similarly, all the terms with s " i form a Poisson bracket. The sum of these two contributions is a direct sum of Poisson brackets which is in turn a Poisson bracket verifying in particular the Jacobi identity. The Poisson brackets referring to the electron and ion quantities correspond to those already discussed in other Hamiltonian reduced fluid models. The reader can in particular refer to Tassi (2015) and Keramidas Charidakos et al. (2015) for the verification of the Jacobi identity for brackets of such form and for a discussion of the corresponding Casimir invariants. The model GF6, in the two-dimensional limit when the dependence on the z coordinate is suppressed, can also be cast in the form of a system of advection of equations for Lagrangian invariants, as is the case for several other reduced fluid and gyrofluid models (see, e.g. Keramidas Charidakos et al.  Schep et al. (1994)).
Remark: The gyroaverage operators G 10 with a form different from (3.16) have been proposed in the literature (see, e.g. Dorland & Hammett (1993); Snyder & Hammett (2001); Scott (2010)), when the expansion (3.18) is not assumed, and are frequently adopted. In particular, G 10 pb s q " Γ 1{2 0 pb s q was shown to provide better agreement with the linear theory at large b s (Dorland & Hammett 1993). We point out, however, that many important features of the model (3.3)-(3.8), such as the total energy conservation and the Hamiltonian structure, are guaranteed whatever the form of the operators G 10 and G 20 of Eq. (3.16) is, provided these operators are linear and formally symmetric, in the sense defined above. This is in particular the case with G 10 pb s q " Γ 1{2 0 pb s q. These issues are also discussed by Mandell et al. (2018).

Six-field model with Landau closure (GF6L)
The variant GF6L of the six-field gyrofluid model, accounting for Landau damping, corresponds to the system (3.3)-(3.8) with Q s given by In Eq. (3.35) we introduced the constant α s " p2{πq 1{2 pm i {m s q 1{2 . The operator L holds for the Landau damping operator. Its modeling in the nonlinear regime is discussed in Tassi et al. (2018). In the linear approximation, it reduces to the negative of the Hilbert transform in the direction of the ambient magnetic field (here taken in the z direction). The presence of this Landau operator in reduced fluid models breaks the Hamiltonian structure by violating energy and Casimir conservation ; Grasso et al. (2020)). Its purpose, on the other hand, is to introduce terms that allow the linear dispersion relation of the gyrofluid model to reproduce that of the parent gyrokinetic model.

Four-field model with quasi-static closure (GF4)
The second variant GF4 is obtained by retaining the evolution equations for N s and M s and imposing the quasi-static closure on the parallel temperature fluctuations T s . Considering Eqs. (3.21) and (A 21), this amounts to setting s " e, i. (3.36) The resulting model reads If one neglects electron FLR effects (i.e. b e Ñ 0), parallel magnetic perturbations (i.e. B " 0), equilibrium temperature anisotropies (i.e. Θ e " Θ i " 1) and sets G 10 pb i q " Γ 1{2 0 pb i q (i.e. takes the alternative form of the ion gyroaverage operator mentioned in the above remark in Sec. 3.1), GF4 reduces to the Hamiltonian gyrofluid model of Waelbroeck & Tassi (2012), the latter taken in the limit of vanishing magnetic curvature and equilibrium density gradients. We note, however, that, once that the quasistatic closure relations are determined, as in Eq. (A 21), all the terms in GF4 (and likewise for GF6), are determined exactly. In particular no approximations of the gyroaverage operators (unlike, for instance, in Waelbroeck & Tassi (2012) and Scott (2010)) are carried out. Terms involving gyroaverage operators are determined exactly also in Brizard (1992), but without making use of the quasi-static closure. As a result, in Brizard (1992), nonlinear terms involving more than one gyroaverage operator do not result in having the single canonical bracket structure (as is the case in GF4, for instance, with the term sgnpq s qrG 10s φ, G 10s A s appearing in the second line of Eq. (3.38)) and which is crucial for determining the Lie-Poisson Hamiltonian structure. This property, which follows from the quasi-static closure, differentiates the models presented in our paper also from its closest predecessors, i.e. the gyrofluid models constructed with the technique recently presented in Tassi (2019). The latter gyrofluid models, in fact, are also Hamiltonian and account for equilibrium temperature anisotropy, but adopt a different closure. Namely, all gyrofluid moments involving finite powers of the magnetic moment µ s (e.g. the perpendicular temperature fluctuations) are set equal to zero. This allows for a Hamiltonian structure but the terms involving gyroaverage operators are not all determined exactly from taking the moments of the gyrokinetic equations. For this reason, we think that, for situations where the quasi-static assumption |ω{pk z v th s q| ! 1 is satisfied, the models presented in this paper might be preferable to the models described in Tassi (2019).

KAWs dispersion relation
We first discuss the KAW dispersion relation as predicted by the gyrofluid models, with and without Landau damping (GF6L and GF4 respectively), in comparison with the predictions of the linearized parent model (Kunz et al. 2018) which identifies with the low-frequency kinetic theory described in Kuznetsov et al. (2012). This latter dispersion relation involves the plasma response function Rpζ s q of the particles of species s, which is related to the corresponding plasma dispersion function Zpζ s q by Rpζ s q " 1`ζ s Zpζ s q with ζ s " ω{p2 1{2 |k z |v th }s q. Different Padé approximants R ij pζ e q (for which we follow the notations of Hunana et al. (2019)) are used to estimate the electron response function.
In order to test the validity of the quasi-static assumption that affects the form of the FLR terms, independently from the effects resulting from Landau damping, we are led to compare the prediction of the 4-field gyrofluid model GF4 with the gyrokinetic dispersion relation obtained by choosing for the electron plasma response function, the function R 21 pζ e q " 1{p1´2ζ 2 e q (model denoted GKNL). This choice directly results from the assumption T }e " 0 (see Eq. (B7) of Passot & Sulem (2007)) and ignores electron Landau damping:. Predictions of the GF6L model will, on the other hand, be compared with the full gyrokinetic dispersion relation, referred to as GK. In the latter description, we use for the electron Padé approximant the function R 86 pζ e q, which shows an excellent agreement with the exact response function. For the ions, we use in all the cases the Padé approximant R 20 pζ i q " 1{p1´2ζ 2 i´i π 1{2 ζ i q (or R 21 pζ i q in the absence of Landau damping). Higher order Padé approximants give almost identical results.
Because we are making use of several acronyms to refer to the different adopted models, we summarized them in Table 1. In all the examples presented in this Section, Θ i " 1. Introducing the ion to electron parallel temperature ratio at equilibrium τ } " T 0 }i {T 0 }e " τ Ki Θ e {Θ i and denoting by α the angle between the wave vector and the ambient magnetic field (propagation angle), we compare, in Fig. 1 (left), the (real) normalized KAW frequency ω{pk z c A q for β Ke " 1, τ } " 1, Θ e " 1, α " 89˝as a function of kd i (where k " a k 2 K`k 2 z and d i is the ion inertial length defined by d i " c A {ω ci , with c A " a B 0 {p4πm i n 0 q the Alfvén velocity), calculated using GF4 (black diamond : Note that, choosing a closure that sets equal to zero a higher even order moment in the hierarchy of parallel moments would only improve the matching with kinetic theory at very large ζe but not in the zero-ζe limit. Differently, closures at an odd order are not consistent with the quasi-static assumption, as they rather correspond to an adiabatic regime.  , with the GKNL prediction (green solid line). An excellent agreement is found, the slight deviation appearing as kd i ą 40 probably resulting from the failure of the quasi-static approximation used in the calculation of all the moments starting from the temperature. For the same plasma parameters, the case with Landau damping is displayed as a red solid line for GK and blue diamond symbols for GF6L, with a good agreement up to scales kd i AE 20. The damping rate is displayed with the same symbols in Fig. 1 (middle). Interestingly, the prediction of GF6 (not shown) departs from GKNL significantly at all the scales. Indeed, in GF6, the parallel temperature fluctuations that obey a dynamical equations with zero heat flux do not approach a quasi-static dynamics. Such a non-dissipative odd-order closure rather fits an adiabatic regime (see footnote in this Section). Landau damping is requested to ensure convergence to a quasi-static regime. In that sense, GF4 is preferable to GF6 for addressing a non-dissipative problem, the interest of the latter model being mainly to provide a framework where Landau damping can be supplemented to an otherwise Hamiltonian description.

Firehose instability
Figure 1 (right) displays for β Ke " 0.5, τ } " 10´3, α " 89˝, Θ e " 0.16, the growth rate of the firehose instability as a function of kd i predicted by GF4 and GF6L (black and blue diamond symbols respectively), together with the GKNL and GK dispersion relations (green and red solid lines). In all the cases, the agreement is excellent, the quenching of the instability being in particular well reproduced.
The agreement found between the KAWs dispersion relations predicted by the gyrofluid models and the gyrokinetic theory, even when limited to scales such that kd i AE 20, and despite a large value of ζ i " b 2Θ e {pβ Ke τ } q ω{p|k z |c A q, can be attributed to the fact that, at least within the linear theory, ion acoustic and kinetic Alfvén waves remain Figure 1. Left: Normalized real part of the KAW frequency ω{pkzcAq for βKe " 1, τ } " 1, Θe " 1, α " 89˝, as a function of kdi, calculated using GF4 (black diamond symbols) and GF6L (blue diamond symbols) models (respectively without and with Landau damping) and with the gyrokinetic dispersion relations GKNL (green solid line) and GK (red solid line). Middle: Damping rate from the GK theory (red solid line) and from the GF6L model (blue diamond symbols). Right: Normalized growth rate of the firehose instability as a function of kdi, for βKe " 0.5, τ } " 10´3, Θe " 0.16, α " 89˝, as predicted by GF4 and GF6L models, together with GKNL and GK dispersion relations, using the same graphic conventions. essentially decoupled. The influence of the ion closure relation on the KAW properties thus remains limited. Deviations from the gyrokinetic theory at small scales are mostly due to the fact that when reaching these scales ζ e " a m e {m i a 2Θ e {β Ke ω{p|k z |c A q becomes non-negligible. Figure 2 concerns a similar comparison in the case of the field swelling instability discussed in Appendix C, for β Ke " 1, τ } " 1, α " 89˝, Θ e " 2.2 (left). For both the cases with and without Landau damping, the relatively good agreement found at large scale between the gyrofluid theories and the gyrokinetic ones deteriorates at smaller scales. The stabilization scale is however correctly captured. This discrepancy is associated with a value of ζ i , which is not small enough (in this case, ζ e remains reasonably small). We show in Fig. 2 (right) with Θ e " 2.41, again with β Ke " 1 and α " 89˝, that a much better agreement can be found when ions are hotter (τ } " 50), which prescribes a small ζ i . The case with cold ions (τ } " 10´5), for which the ion dynamics is decoupled, is displayed in Fig. 3, showing an even better agreement. The left panel displays the real part of the slow wave for Θ e " 1, obtained with the GKNL model (solid green line) or the GF4 model (diamond symbols). The right panel shows the growth rate of the field-swelling instability for Θ e " 2.01 (keeping unchanged the other parameters).

Conclusion
We derived a 6-field Hamiltonian gyrofluid model, referred to as GF6, retaining the gyrocenter density, the parallel velocity and temperature fluctuations for each species, under the sole assumption that all the other gyrocenter moments are calculated from the quasi-static linear kinetic theory. Such an assumption on the closure turns out to yield exact expressions for all the terms of the model, without, in particular, requiring approximated expressions for the terms involving gyroaverage operators. Nonlinear terms involving more than one gyroaverage operator, in particular, appear in the form of a single canonical bracket, which naturally lets the model fit in the class of Hamiltonian models with a Lie-Poisson structure. The model accounts for equilibrium temperature Figure 2. Normalized growth rate of the field swelling instability versus kdi for α " 89a nd βKe " 1, Left: Predictions of GF4 and GF6L, compared with those of GKNL and GK respectively, for τ } " 1 and Θe " 2.2. Right: Prediction of GF4 compared with that of GKNL, for τ } " 50 and Θe " 2.41 . Same graphic conventions as in Fig. 1 are used. Figure 3. Left: Normalized slow-wave frequency versus kdi in regimes where the ions are cold (τ } " 10´5), with βKe " 1 and α " 89˝, as predicted by the GF4 model and the GKNL dispersion relation. Right: normalized growth rate of the swelling instability for τ } " 10´5 and Θe " 2.01. Same conventions as in Fig. 1 are used. anisotropy and also retains both ion and electron FLR corrections, electron inertia and parallel magnetic fluctuations. In a variant of the model (GF6L) parallel Landau damping is retained through a Landau-fluid modelization of the gyrocenter parallel heat fluxes. A second variant of the model (GF4) is obtained by prescribing parallel isothermality, which still falls in the frame of the quasi-static closure and allows for a Hamiltonian formulation. The comparison of the dispersion relations of KAWs and SWs predicted from GF6L or GF4, with those derived from the parent gyrokinetic theory where the plasma response function is replaced by a Padé approximant, provides an estimate of the maximal transverse wavenumber beyond which the phase velocity of the corresponding wave is too large compared with the electron (in the former case) or ion (in the latter case) parallel thermal velocity for consistency with a closure condition based on a quasi-static assumption. It turns out that the agreement extends to transverse scales significantly smaller than the ion Larmor radius in the case of KAWs, mostly because, at least at the linear level, SWs and KAWs are essentially decoupled, making the influence of the ion closure relation on the KAW properties relatively limited. This situation contrasts with the case of the SWs for which the dispersion relation is accurately reproduced only at scales larger than a significant fraction of the ion Larmor radius. Under these conditions, the model reproduces the instabilities induced by temperature anisotropy, such as firehose or field-swelling instabilities. It should nevertheless be noted that, as it assumes small perturbations of an equilibrium state, the model does not permit evolution of the mean temperatures, an effect usually considered as contributing efficiently to the saturation of these instabilities. The subcritical nonlinear regime is however expected to be accurately described. The model will in particular be most useful for studying the coupling of KAWs with SWs which can generate large-scale parametric decay instabilities at small β e , a regime especially relevant in the regions of the solar wind relatively close to the Sun explored by space missions such as Parker Solar Probe or Solar Orbiter.
In general, to the best of our knowledge, our gyrofluid model is the only one, at the present moment, possessing the following features, which could make it a valuable tool for local investigations of basic plasma phenomena of interest for space plasmas: it accounts for equilibrium temperature anisotropies as well as parallel magnetic perturbations; it reproduces, in a rather wide range of values of parameters, compatible with the quasistatic closure, quantitative features of known kinetic linear dispersion relations; the model equations, and in particular the terms involving FLR corrections, are calculated exactly, unlike other gyrofluid models which adopt truncations or approximations of such terms; it possesses a Hamiltonian structure.

Appendix A. Derivation of closure relations from the gyrokinetic linear theory
We consider the linearization of the gyrokinetic system (2.1)-(2.4) about the equilibrium state r g s " 0 (or, equivalently, r f s " 0), with r φ " r A " r B " 0. The resulting linear system can be written in the form where we adopted the same notation with the tilde symbol, that we used in Eqs. (2.1)-(2.4), to indicate the dynamical variables r f s of the linearized system and the field perturbations r φ, r A and r B . We introduce the following Fourier series representation: with ω P C indicating the complex frequency. For any given k P D, from Eq. (A 1) we obtain the relation where r ζ s " ω{pk z v th s q. We consider now the quasi-static limit | r ζ s | ! 1. In this limit, the relation (A 7) reduces to In the following, we make use of the relation (A 8), derived from the linear theory, in order to determine the closure relations to insert in the hierarchy of nonlinear gyrofluid equations. For this purpose, we adopt the Hermite-Laguerre expansion of Eq. (3.18) for the perturbation of the distribution function in the linearized system and we write in Fourier representation. From Eqs. (A 9), (A 10), using the orthogonality relations for Hermite and Laguerre polynomials, one obtains Inserting the relation (A 8) into Eq. (A 11) and using the orthogonality relations for Hermite polynomials, one has where the operators G 1n and G 2n are defined by In order to obtain the expressions (A 16)-(A 17) use was made of the orthogonality of Laguerre polynomials as well as of the relations (Szegö 1975) where L p1q n is a generalized Laguerre polynomial. Making use of the Fourier representations (A 10) , (A 5) and (A 6) for for f mns , r φ and r B , respectively, one can deduce from Eq. (A 12) the relation (A 20) or, equivalently, where we also made use of the normalization (3.11) for φ and B . The operators G 1ns and G 2ns are defined, consistently with the definition of G 10s and G 20s given in Sec. 3, by G 1ns f px, tq " ř kPD G 1n pb s qf k ptq exppik¨xq and G 2ns f px, tq " ř kPD G 2n pb s qf k ptq exppik¨xq, for a function f and n ě 0 (in the specific case of the linear dispersion relation, the dependence on t is provided by the factor e´i ωt , but when the relations (A 21) are used as closures for the nonlinear models, the dependence on t is of course left arbitrary).
The relations (A 21) descending from the quasi-static assumption, are adopted as closure relations in GF6 and GF6L for all the moments involved in the model, except for N s , M s , T s , which are derived by solving the evolution equations (3.3)-(3.5). The parallel heat flux fluctuations Q s , on the other hand, are determined, as already mentioned, again by a quasi-static closure for GF6 (Eq. (3.23) which follows from Eq. (A 21) for m " 3, n " 0) or by the Landau closure (3.35) for GF6L. The closure (3.36) adopted for GF4, is again a quasi-static closure following from Eq. (A 21) when m " 2 and n " 0.

Appendix B. Derivation of the model equations
The gyrofluid system (3.3)-(3.8) descends from the parent gyrokinetic system (2.1)-(2.4) upon applying to the perturbations of the distribution functions the expansion (3.18). In order to obtain a closed system with a finite number of equations, such expansion is constrained in the following way. The moments f 00s , f 10s and f 20s (or, equivalently, by virtue of Eqs. (3.20)-(3.21), the gyrofluid densities, parallel velocities and temperatures N s , U s and T s ), for each species s, get determined by evolution equations obtained by making the product of all the terms of the gyrokinetic equation (2.1) with the zero, first and second order Hermite polynomial in the variable v } {v th s and integrating over the velocity volume element dW s . For GF6L, the parallel heat flux Q s gets determined by the relation (3.35). All the other moments, on the other hand, are assumed to be given by the relations (A 20), or, in normalized form, by Eq. (A 21), obtained from the linear theory in the quasi-static limit. With this prescription, the expansion (3.18) becomes where we made use of the fact that H 0 pv } {v th s q " 1, Inserting the expansion (B 1) into Eqs. (2.5) and (2.1), and integrating over dW s one obtains where we made use of the definitions (A 13) and (A 14) and of the orthogonality of Hermite polynomials. We remark, at this point, that the sum of the last term in the first line of Eq. (B 2) with the first term on the second line of Eq. (B 2) yields zero, because of the antisymmetry of the canonical bracket r , s. We thus conclude that the quasi-static closure (A 20) has the remarkable property of annihilating, in the continuity equation, all the contributions associated with the moments f 0ns , for n ě 1. By virtue of this cancellation, from Eq. (B 2) we obtain B Bt Applying to Eq. (B 3) the normalization (3.9)-(3.12), one obtains Eq. (3.3).
In order to derive Eq. (3.4) we point out first that, with the help of the identities we obtain the relation (see also Brizard (1992)) and, by an analogous procedure, the relation Upon multiplying Eq. (2.1) by p1{n 0 qv } {v th s and integrating over dW s one obtains, adopting the expansion (B 1) as well as the relations (A 13), (A 14), (B 10) and (B 11), the following equation Also in this case, the quasi-static closure leads to a remarkable cancellation. Indeed, among the nonlinear terms involving only electromagnetic quantities (i.e. r φ, r A and r B ) all those containing gyroaverage operators G 1ns and G 2ns , with n ě 1, vanish. As a The transverse velocity can be decomposed into compressible and solenoidal parts by writing u K "´∇ K χ c`∇Kˆp χ s p zq.
(C 1) One immediately gets (see e.g. Eqs. (48), (52) and (56) of Passot & Sulem (2006) where FLR corrections are neglected) where p K and p } denote the total (ion plus electron) perpendicular and parallel pressure fluctuations, with the subscript 0 referring to the equilibrium values. Furthermore, ρ 0 denotes the equilibrium plasma density. In Eq. (C 2), the perpendicular pressure fluctuations are given by the drift-kinetic theory, as found e.g. in Eq. (27) of Snyder et al. (1997), in the form where q r " e for ions and´e for electrons, respectively. The potential ψ is defined in terms of the parallel electric field by E z "´B z ψ. These expressions for the perpendicular pressure and density perturbations identify with the large-scale limit of formulas given in Appendix B of Passot & Sulem (2007). Quasineutrality requires the equality of the electron (n e ) and ion (n i ) number-density fluctuations, which prescribes Plugging the expressions for the pressures in Eqs. (C 2)-(C 3) thus leads to 2p1´Θ e Rpζ e qq`τ Ki p2´Θ e Rpζ e q´Θ i Rpζ i qq We note that ζ 2 i " T 0 e m i is the sound speed based on the parallel electron temperature.
C.2. Link with the fluid theory In order to make a link with fluid theory as performed in Basu & Coppi (1984), we note that keeping all the terms in the ion density and parallel velocity equations (in particular the time derivatives), is equivalent to expanding Rpζ i q for ζ i large, i.e. in the adiabatic limit where Rpζ i q «´1{p2ζ 2 i q. It thus follows that the ions cannot be hot, at least outside an angular boundary layer near the transverse direction.
In addition to adiabatic ions, let us also assume a quasi-static limit for the electrons, which leads to assume ζ e Ñ 0 and thus Rpζ e q « 1. Equation (C 7) then rewrites For cold ions (τ } " τ Ki " 0), we get which also rewrites in the form of Eq. (44) of Basu & Coppi (1984) The slow mode is obtained when ω " ω s , while the fast mode corresponds to ω " ω s .
As Θ e is further increased, a pair of purely imaginary complex conjugate roots appears, leading to the so-called slow-mode swelling instability.

C.3.2. Fast mode
According to the theory of Basu & Coppi (1984) that assumes ζ e very small, the fast mode becomes unstable when The validity conditions require in particular that pk z {k K q 2 ą pm e {m i qp1{β }e q in order to ensure that ζ e is small enough.

C.4. Instability growth rate
Le us now consider the full dispersion relation, given by Eq. (C 7), with a general electron response function. The instabilities are illustrated in Fig. 4 which displays, as a function of the perpendicular to parallel electron temperature anisotropy Θ e , the  7), for propagation angles α " 89˝(left) and 80˝(middle) as a function of the perpendicular to parallel electron temperature anisotropy Θe in the case βKe " 1 and cold ion temperatures (leading to take Rpζiq "´1{p2ζ 2 i q) for various approximations of Rpζeq: BC84 (black solid line) uses Rpζeq " 1, DK (red solid line) uses R20pζeq and DKNL (green solid line) uses R20pζeq (no electron Landau damping). Superimposed diamond symbols refer to predictions of GF4 model taken in the large-scale limit. Right panel corresponds to DKNL and GF4 models with α " 89i n the case τKi " 1.
imaginary part of r " ω{kc A for the unstable mode, solution of Eq. (C 7) with β Ke " 1 and cold ion temperatures, for propagation angle α " 89˝(left) and α " 80˝(middle), when using four different approximations for the plasma response functions. The first one, referred to as BC84 (black solid line) uses Rpζ i q "´1{p2ζ 2 i q and Rpζ e q " 1 and corresponds to the asymptotic solution of Basu & Coppi (1984), the second one, denoted DK (red solid line) uses for the electrons the Padé approximant R 20 pζ e q " 1{p1´2ζ 2 e´i π 1{2 ζ e q (the use of R 42 leads to almost identical results) and for the ions the same large-ζ i limit as in the BC84 model, since the ions are cold. The third model, called DKNL (green solid line) differs from the previous one by the fact that electron Landau damping is suppressed (using R 21 pζ e q " 1{p1´2ζ 2 e )). Superimposed diamond symbols refer to the prediction of GF4, taken in the large-scale limit. It is legitimate to consider this model, derived in the limit ζ i Ñ 0, in the regime of cold ions, because in this case the dynamics is insensitive to the closure assumption. The right panel corresponds to the case with τ Ki " 1 for α " 89˝. It shows that the GF4 model is still approximately valid even in a situation where ζ i is not small. Figure 5 displays, as a function of Θ e , the positive real part of the three roots of Eq. (C 7) (referred to as min, int and max in increasing order of magnitude, displayed in turquoise, magenta and brown colors respectively), in the case of DKNL model (solid lines) or for the DK model (dash-dotted lines), with the predictions of the GF4 model superimposed as diamond symbols, again for 89˝(left) and 80˝(right) propagation angles and cold ions. For the 89˝angle (which falls outside the range of admissible angles for the BC84 model), the destabilization of the slow mode for Θ e ą 2 is almost similar for the four models. Within the BC84 model, the fast mode becomes unstable for Θ e ą 4: the real part of the associated root vanishes (not shown) and the instability growth rate (black solid line of the left panel of Fig. 4)   is here the only extra mode for the present choice of the Padé approximant (R 21 pζ e q, leading to what we called DKNL)). The instability that continues to exist for Θ e ą 4 (its growth rate corresponds to the green solid line in the left panel of Fig. 4) is associated to the destabilization of this extra plasma mode (while both the fast and slow modes continue to exist and remain stable). The only small difference observed in the presence of Landau damping is that the intermediate mode becomes purely imaginary in a range of values of Θ e between 2 and 4 (see Fig. 5 left, dash-dotted line). We also note that the use of the R 42 Padé does not change the roots associated with the slow and fast modes, but only those associated to the extra damped plasma modes (not shown). We conclude that in an angular boundary layer close to 90˝, the fast mode is always stable, and the instability that continues to exist for Θ e ą 4 is of the same nature as the slow mode swelling instability, i.e. its growth rate tends to zero as α approaches 90˝.
For α " 80˝, the value of ζ e is sufficiently small for the approximation Rpζ e q " 1 to be valid. In this case, the behavior of the slow mode swelling instability is similar to the previous case, the main difference affecting the fast mode. Its phase velocity for Θ e ă 4 now corresponds to the intermediate branch, the one with the largest real part corresponding to the extra plasma mode which, in the presence of Landau damping, is heavily damped (damping rate not shown). For Θ e ą 4, the real frequency of the fast mode vanishes and this mode becomes unstable, as predicted in Basu & Coppi (1984). The slow mode reappears on this intermediate branch with a very good match between the GF4 and DK as well as DKNL predictions.
Complementary information is presented in Fig. 6 which displays the growth rate of the unstable mode as a function of the angle α for the three models described above for Θ e " 3 (left) and Θ e " 4.5 (right). For the case corresponding to the slow mode swelling instability of Basu & Coppi (1984) (Θ e " 3), the three models are very similar (the case without Landau damping is actually almost identical to the BC84 model so that both curves are superimposed). The GF4 model gives very similar growth rates for angles between 80˝and 90˝but its predictions deviate for smaller angles. For Θ e " 4.5 the fast mode instability as predicted by BC84 displays a growth rate proportional to k Figure 6. Growth rate of the unstable mode as a function of the propagation angle α for the BC84 (black), DK (red), DKNL (green) and GF4 (diamond symbols) models in the case Θe " 3 (left) and Θe " 4.5 (right). For Θe " 3, BC84 and DKNL curves are superimposed.
up to angles α " 90˝. This fast-mode instability is only recovered with DKNL (and DK) at oblique angles, the deviation with BC84 starting to be significant for α ą 75˝. For this value of Θ e , the slow mode is always stable in the GF4 model.

Remarks:
‚ Influence of warm ions: If one considers the fast mode at an angle close (but not equal) to 90 degrees, one can assume warm ions and at the same time Rpζ i q " 0. Taking also Rpζ e q " 1, and using Eq. (C 7), one gets the dispersion relation (Eq. 20) obtained by Pokhotelov & Onishchenko (2014) who show that ions are stabilizing.
‚ The case with ion temperature anisotropy: With finite (but isotropic) ion temperatures, other modes are present, but the one which becomes first unstable when electron temperature anisotropy is increased is still the slow mode. The case where the instability is driven by ion temperature anisotropy (so-called classical mirror instability with isotropic electrons) is in contrast different since the instability originates from the extra mode associated with the finite ion temperature fluctuations and not from the slow mode which continues to exist and to be stable above the mirror threshold. An interesting point is that the mirror mode, usually thought of being non-propagating, originates from one of the damped propagating "ion temperature modes". These modes only become non-propagating for a large enough ion temperature anisotropy.