## Abstract

The fully nonlinear Grad's *N*×26-moment (*N*×*G*26) equations for a mixture of *N* monatomic-inert-ideal gases made up of Maxwell molecules are derived. The boundary conditions for these equations are derived by using Maxwell's accommodation model for each component in the mixture. The linear stability analysis is performed to show that the 2×G26 equations for a binary gas mixture of Maxwell molecules are linearly stable. The derived equations are used to study the heat flux problem for binary gas mixtures confined between parallel plates having different temperatures.

## 1. Introduction

Gases in nature are usually present in the form of mixtures and describing the flows of gas mixtures is a challenging problem, especially in rarefied processes. Rarefied or micro gas flows are characterized by a dimensionless flow parameter, the Knudsen number (*Kn*) which is defined as the ratio between the mean free path of a molecule and a macroscopic length scale in the problem. Based on the Knudsen number, the flow is classified into different regimes. The gas flows in all regimes are well described by the Boltzmann equation which, however, is extremely difficult to solve numerically or analytically. Furthermore, the direct solution of the Boltzmann equation or the solution obtained via the direct simulation Monte Carlo (DSMC) method [1] are computationally expensive, especially in the transition regime (intermediate Knudsen number flow) [2]. Also, DSMC cannot be run for dilute gases because there are not sufficient sample molecules. On the other hand, classical gas dynamic models like Fick's law are valid only for very small Knudsen numbers. Therefore, proper models to describe the flows in the transition regime are desirable and the models are developed through kinetic theory.

In kinetic theory, one looks for approximate solutions of the Boltzmann equation. The two well-known approximation techniques in kinetic theory are the Chapman–Enskog (CE) expansion method [3] and Grad's method of moments [2,4], and both techniques also work for gas mixtures. However, at least in the case of a single gas, the higher order approximations based on CE expansion suffer from instabilities [5] while moment equations obtained via Grad's method of moments are always stable and have been used successfully [6]. In order to avoid some drawbacks of classical moment methods, Struchtrup & Torrilhon [7] regularized the moment equations by performing an expansion on Grad's 13-moment (G13) equations, which led to the regularized 13-moment (R13) equations for a single gas. In the last decade, R13 equations have proved to be very successful in describing the rarefied gas flows (e.g. [8–12]). As low-dimensional continuum models, moment equations are an efficient alternative to direct numerical discretizations of the Boltzmann equation, especially at moderate Knudsen numbers. Unfortunately, they are not available for gas mixtures at present.

Moment equations for mixtures have been considered in the context of extended thermodynamics in [13,14], but typically simplified models for the collision operators have been used and specific equations have only been given up to 13 moments. The papers, e.g. [15,16], see also the text book [17], study Grad's moment method in a multi-component approach for plasma models. Complete 13-moment equations are given based on the linearized Boltzmann collision operator with arbitrary interaction potentials. However, boundary conditions are not derived and the equations are exclusively used to derive first-order constitutive laws based on CE expansion. Recent works [18,19] extend the moment approach to relativistic gas mixtures and polyatomic gases, respectively.

As a first step towards regularized moment equations for rarefied gas mixtures, this paper derives the fully nonlinear *N*×G26 moment equations from the complete Boltzmann equation for a gas mixture of *N* components made up of Maxwell molecules. To the best of our knowledge, the derivation of such a system has never been attempted because—among other reasons—the production terms which appear from the nonlinear collision term of the Boltzmann equation are very hard to evaluate explicitly. In this paper, we use a generalization of Grad's ansatz for the distribution function to mixtures and evaluate the higher moments as well as the collision integrals for each mixture component described with 26 moments. Additionally, we use Maxwell's accommodation model to derive appropriate wall boundary conditions for each component of the mixture. These boundary conditions are equivalent to the diffuse-specular scattering boundary condition used in [20,21] in the hydrodynamic setting.

Our equations are based on a multi-fluid approach which considers different field variables for the different components including multi-temperatures, that is, different temperatures for different species in the gas mixture. Single temperature theories which assume a common average temperature for the whole gas mixture is used more often than the multi-temperature theories in the literature for simplicity (e.g. [3,22–24]). However, multi-temperatures are expected to be more authentic in many physical situations. Especially in plasma physics these models become indispensable to use [17,25].

After deriving the equations, the wave modes of the final system of equation are studied in a linear wave analysis and the 2×G26 equations for binary mixtures are shown to be stable in time for spatial perturbations. Moreover, using the linearized 2×G26 equations and the boundary conditions, we study the classical heat conduction problem for a binary gas mixture confined between two parallel plates having different temperatures considered in [20,26–28] for instance.

The aim of this paper is to provide explicit moment equations as a possible tool for describing mixtures of rarefied gases. The equations are demonstrated to be stable, possess proper boundary conditions for walls and show agreement with kinetic solutions. Many open questions remain to be answered in future work, such as more comparison with DSMC/Boltzmann solutions in other geometries, a formal analysis of the asymptotic order of the system and the regularization of the system in order to simulate nonlinear phenomena like shock waves. Additionally, different molecule models like hard spheres or power potentials can be considered within the framework of this paper.

The structure of the paper is as follows. The key ingredients of kinetic theory for binary gas mixtures are introduced in §2. The fully nonlinear *N*×G26 equations and the boundary conditions are derived in §3. The linear stability analysis of 2×G26 equations is detailed in §4. The heat flux problem for a binary gas mixture confined between parallel plates having different temperatures is considered and the results are compared with those of [20] in §5. The paper ends with a conclusion and outlook in §6.

## 2. Kinetic theory for mixture of gases

In both the CE expansion method and Grad's method of moments, the field variables can be defined as averages over polynomials in self peculiar velocities, i.e. the peculiar velocity of the corresponding single gas (e.g. [22]). However, following Chapman & Cowling [3], we shall define the field variables as averages of polynomials in peculiar velocities with respect to the whole mixture. To be more precise, let us use the following notations for different velocities of the *α*-constituent in a gas mixture: the instantaneous velocity of an *α*-gas molecule *c*_{α}, the macroscopic velocity of the *α*-gas molecules *v*_{α}, the self-peculiar velocity of the *α*-gas molecule *c*_{α}−*v*_{α}, the macroscopic velocity of the mixture ** v**, the peculiar velocity of the

*α*-gas molecule with respect to the mixture

*C*_{α}=

*c*_{α}−

**and the diffusion velocity of the**

*v**α*-gas molecules

*u*_{α}=

*v*_{α}−

**. We shall define the moments as averages of polynomials in peculiar velocities with respect to the mixture,**

*v*

*C*_{α}.

### (a) Boltzmann equation

The state of a mixture consisting of *N* monatomic-inert-ideal gases in the phase space is characterized by the distribution functions (also termed as the phase densities) of individual gases, i.e. by *f*_{α}≡*f*_{α}(** x**,

*c*_{α},

*t*) for

*α*=1,2,…,

*N*in such a way that

*f*

_{α}d

**d**

*x*

*c*_{α}gives the number of molecules of the

*α*-constituent at time

*t*in an infinitesimal volume d

**centred around**

*x***whose velocities belong to an infinitesimal volume d**

*x*

*c*_{α}in velocity space located around

*c*_{α}. The distribution function of the

*α*-constituent,

*f*

_{α}is governed by the Boltzmann equation

*t*≡∂/∂

*t*+

**⋅**

*v***∇**is the material derivative,

*F*_{α}is the external force per unit mass acting on the

*α*-constituent, the phase densities with primes correspond to the phase densities with post-collisional velocities

*g*_{αβ}=

*C*_{α}−

*C*_{β}=

*c*_{α}−

*c*_{β}is the relative velocity of an

*α*-gas molecule with respect to a

*β*-gas molecule,

*b*is the collision parameter and the angle

*ϵ*describes the orientation of the collision plane. In (2.1) and throughout this paper, the integrations over any velocity space are the volume integrals over all components of that velocity and the limits of integrations are

*b*and

*ϵ*, to make the notations compact.

### (b) Moment approximations

Let *m*_{α} be the mass of an *α*-gas molecule for *α*=1,2,…,*N*. We introduce the general form of a moment for the *α*-constituent as
*ρ*_{α}(** x**,

*t*), diffusion velocity

*u*_{α}(

**,**

*x**t*), temperature

*T*

_{α}(

**,**

*x**t*), stress tensor

*k*is the Boltzmann constant,

*n*

_{α}=

*ρ*

_{α}/

*m*

_{α}is the number density of the

*α*-constituent and

*θ*

_{α}=

*kT*

_{α}/

*m*

_{α}is the temperature of the

*α*-constituent in energy units. Furthermore, the macroscopic velocity of the

*α*-constituent

*v*_{α}(

**,**

*x**t*) and that of the whole mixture

**(**

*v***,**

*x**t*) are given by

_{2},

In order to derive moment equations for the *α*-constituent, one first multiplies the Boltzmann equation (2.1) by an arbitrary function *ψ*_{α}≡*ψ*(** x**,

*c*_{α},

*t*) and integrates over the velocity space

*c*_{α}. After some algebra, one obtains the transfer equation

*production term*or the

*Boltzmann collision integral*corresponding to moment

Moment equations for the *α*-constituent in a gas mixture are obtained by substituting appropriate values for *ψ*_{α} in the transfer equation (2.5). The typical tensorial structure of *ψ*_{α} is: *ψ*_{α} can be chosen in infinitely many ways, the moment equations form a coupled system of infinitely many partial differential equations. However, in practice, only a finite set of *ψ*_{α} is chosen, which leads to a finite system of moment equations. Nevertheless, the system of moment equations is not closed because each moment equation contains a flux term comprising higher order moment(s) on the left-hand side as well as the unknown production term on the right-hand side.

The system of moment equations is closed with a certain approximation for the phase density. A standard choice for the approximation is an expansion of the phase density into Hermite polynomials, which yields so-called G13 or G26 distribution functions [2,4,22]. With this G13 or G26 closure, the unknown higher order fluxes and the production terms in the moment equations are expressed in terms of field variables considered and, eventually, 13- or 26-moment equations boil down to a closed system of G13 or G26 equations, respectively.

## 3. Moment equations for gas mixture

### (a) Full nonlinear equations

The first 26 moment equations for the *α*-constituent corresponding to 26 field variables
*ψ*_{α} in the transfer equation (2.5) with each element of
*α*-constituent (*α*=1,…,*N*) read
*α*-constituent. The mass, momentum and energy balance equations for the whole mixture can be obtained by summing (3.3)–(3.5) over *α*=1 to *N*. They turn out to be
*T* is the total temperature of the mixture. Note that the summations over production terms in (3.11) vanish owing to mass, momentum and energy conservation, respectively. Owing to (2.4)_{3}, it is possible to use the partial momentum balance equations (3.4) only for (*N*−1) components. In any case, the total momentum balance equation (3.11)_{2} needs to be added to the system (3.3)–(3.10) in order to evaluate the total velocity. The use of the total mass and energy balance equations is optional and could replace one of the *N* partial mass balance equations (3.3) and one of the *N* partial energy balance equations (3.5).

Clearly, the system of moment equations (3.3)–(3.10) is not a closed system, as it contains the unknown higher order moments—

*Closure*: In order to close the system of moment equations, we assume a Grad-type distribution function for each component in the mixture, e.g.
_{i1⋯ip} may depend on the non-equilibrium moments. Note that _{2} by choosing *T*_{α} as the temperature of the whole mixture, which is same as the temperature of each gas in equilibrium. The ansatz in (3.12) forms the basic ingredient of the moment theory. Many possibilities are available, each trying to optimize the flexibility and physicality of the distribution function while remaining easy to evaluate. In the single gas case, the Gaussian distribution function is chosen as local or global equilibrium, i.e. Maxwellian, and the non-equilibrium variables enter the ansatz through the linear coefficients of the sum. Maximum entropy closures [29] use fully nonlinear representations of the distribution function. In (3.12)_{2}, we choose to include the *α*-temperature into

For 26 field variables (3.1), the Grad-type distribution function (3.12)_{1} results into G26 distribution function for each component in the mixture, e.g.
*α*-constituent in a gas mixture of *N* gases made up of Maxwell molecules, using the computer algebra software Mathematica. The nonlinear production terms of the *α*-constituent for mass, momentum and energy balance equations are as follows:
*μ*_{γ}=*m*_{γ}/(*m*_{α}+*m*_{β}) for *γ*∈{*α*,*β*} is the mass ratio of the *γ*-constituent. The collision frequency for cross-collisions of the *α*-gas molecules with *β*-gas molecules is given by
*α*- and *β*-gas molecules is balanced, i.e. *n*_{α}*ν*_{αβ}=*n*_{β}*ν*_{βα} [31,32]; the constants *a*_{i} are defined as *a*_{i}=*A*_{i}/*A*_{2} with
*a*_{2}=1, *Ω*-integral, *α*- and a *β*-gas molecules and an explicit expression for *α*-constituent are quite involved and, therefore, they are given in appendix A.

Equations (3.3)–(3.10) along with (3.11)_{2}, (3.14), (3.15) and (A 2)–(A 5) form fully nonlinear G26 equations for the *α*-constituent in a mixture of *N* gases, made up of Maxwell molecules. The G26 equations for other constituents in the gas mixture can be written analogously. Thus, we have derived the fully nonlinear *N*×G26 equations for mixtures of *N* gases of Maxwell molecules.

### (b) Boundary conditions

Besides a large number of equations and production terms, another difficulty in dealing with moment equations is that they require proper boundary conditions. Maxwell's accommodation model [33] can be used to derive the boundary conditions associated with the moment equations in the case of a single gas [34–36]. In the case of a gas mixture, we consider Maxwell's accommodation model for each component separately.

#### (i) Maxwell's accommodation model

*Single gas case*. According to Maxwell's accommodation model, a fraction (1−*χ*) of molecules hitting the wall is specularly reflected into the gas and the remaining fraction *χ* of molecules is accommodated at the wall and they return into the gas with a Maxwell distribution function of the wall
*m* is the mass of a molecule, ** c** is the instantaneous velocity of the molecule,

*ρ*

_{w}is a parameter which is determined according to molecule conservation at the wall, and

*T*

_{w}and

*v*_{w}are the temperature and the velocity of the wall, respectively.

*χ*is termed as the accommodation coefficient for the gas.

According to Maxwell's accommodation model, the phase density in the infinitesimal neighbourhood of the wall is given by
** n** denotes unit normal from the wall pointing into the gas,

*f*

_{gas}is the distribution function of molecules hitting the wall,

*f**

_{gas}is the distribution function of specularly reflected molecules from the wall, and the subscript ‘

*nw*’ is used to denote ‘near wall’. The distribution function

*f**

_{gas}follows from

*f*

_{gas}analogously by reversing the sign of normal component of velocity

*c*

_{n}.

*Extension to gas mixture*. As in the case of a gas mixture, the accommodation coefficient as well as the parameter *ρ*_{w} for each component are not the same in general, we extend Maxwell’s accommodation model for a single gas to a mixture of gases by writing the Maxwell distribution function of the wall and the distribution function near the wall for each component in the mixture separately. More precisely, for any gas *α* in the gas mixture, the Maxwell distribution function of the wall and the distribution function near the wall are written by replacing *f*_{gas}, *f**_{gas}, *f*_{w}, *f*_{nw}, ** c**,

*ρ*

_{w},

*m*,

*θ*

_{w},

*χ*in (3.17) and (3.18) with

*f*

^{(α)}

_{gas},

*f*

^{*(α)}

_{gas},

*f*

^{(α)}

_{w},

*c*_{α},

*ρ*

^{(α)}

_{w},

*m*

_{α},

*θ*

^{(α)}

_{w},

*χ*

_{α}, respectively.

#### (ii) Boundary conditions for the *α*-component

The boundary condition for any function *ψ*_{α}≡*ψ*(** x**,

*c*_{α},

*t*) is obtained by assuming that the moment with respect to this function should be continuous near the wall, i.e.

*α*-constituent by choosing the appropriate values of

*ψ*

_{α}, but the problem now is that it provides too many boundary conditions than required. However, by considering the special case of

*χ*=0 for a single gas, Grad [4] concluded that only those moments which are odd functions in

*C*

_{n}need to be considered for the boundary conditions. We borrow Grad’s strategy for a single gas and in the case of gas mixture also, only those moments which are odd functions of

*ψ*

_{α}to

*i*,

*j*∈{1,2} and

*t*

_{1}and

*t*

_{2}are two orthonormal directions on the plane of the wall. By assuming that no

*α*-molecules are produced or eliminated at the wall, i.e.

*ψ*

_{α}provides the definition for the parameter

*ρ*

^{(α)}

_{w}

*ψ*

_{α}yield the following boundary conditions

**=**

*V***−**

*v*

*v*_{w}is the slip velocity. An additional boundary condition resulting from physics states that the normal component of the velocity of the gas mixture at the wall vanishes (

*v*

_{n}=0 at the wall) if the gas molecules cannot penetrate through the wall. The boundary conditions for the other constituents in the mixture can be written analogously. From here onwards, we shall focus on a binary gas mixture of

*α*and

*β*gases.

## 4. Linear stability analysis

### (a) Linearized equations

Consider a binary gas mixture of two gases *α* and *β* without any external forces acting on them. For linear stability analysis, we linearize the 2×G26 equations and the momentum balance equation for the mixture by perturbing the field variables around their respective ground states, i.e.
*ε* is a small parameter. Note that the temperature difference *α*-constituent read
*α* and *β*, respectively, in the ground state;
*δ*_{i}’s for *i*=1,2,…,15 depend only on the mass ratios of the constituents through
*δ*_{i}’s, the relation *μ*_{α}+*μ*_{β}=1 has been used.

The linear-dimensionless G26 equations for the *β*-constituent follow in a straightforward way by interchanging *α* and *β* in (4.3)–(4.10) and the linear-dimensionless momentum balance equation for the mixture reads

The velocity scale *v*_{°} is chosen as *β* or *α*, respectively.

### (b) Dispersion relation

In order to verify the stability of the 2×G26 equations, we consider them in a one-dimensional setting, such that all fields depend only on the *x*-coordinate. It suffices to consider the equations for the variables
_{3}.

Let the length scale be the inverse of the wavenumber *κ*. For the one-dimensional system of moment equations, we assume the plane wave solution of the form
*U*_{0} is the vector containing the complex amplitudes of the corresponding variables in ** U**,

*ω*as the usual complex frequency. In this case, the Knudsen number appearing in the equations above is

*Kn*=ℓ

*κ*, which now takes the role of a dimensionless wavenumber. The system of moment equations can be written as

*ω*and

*κ*(between

*Kn*here).

For spatial disturbances, the wavenumber *κ* is real and the frequency *ω* is complex. The corresponding wave has phase velocity *v*_{ph}=Re(*ω*)/*κ* and damping *ς*=Im(*ω*). For stability, *ς* must be non-positive.

In order to discuss the linear stability of the equations, we consider three binary mixtures of noble gases neon–argon (Ne–Ar), helium–argon (He–Ar) and helium–xenon (He–Xe) as our examples. The molecular masses of these gases are *m*_{He}=4.0026, *m*_{Ne}=20.1791, *m*_{Ar}=39.948 and *m*_{Xe}=131.293 in atomic units. The ratios of collisional cross-sections (4.13) for these mixtures are calculated through the values of *Ω*-integrals for hard spheres: *d*_{α} and *d*_{β} are the diameters of *α*- and *β*-gas molecules. The required diameters of the gases are calculated using the exact expression of viscosity for a single gas given in [26] and the experimental data on the viscosities of the single gases at temperature 300 K given in [37]. The diameters of these gas molecules turn out to be *d*_{He}=0.2166, *d*_{Ne}=0.2564, *d*_{Ar}=0.3606, *d*_{Xe}=0.4821 in nanometres.

Figure 1 illustrates the dispersion modes for the considered binary gas mixtures in *Kn*=*κ*=0. In each figure, two sound modes—the modes starting with zero damping—commence with non-zero velocities and one with zero velocity at *Kn*=*κ*=0. The speeds of sound modes starting with non-zero velocities in He–Ar and He–Xe mixtures with *d*,*f*) decrease significantly before increasing with increasing Knudsen number, which are the cases of mixtures having large difference in the mass ratios and high mole fractions of the heavier gases. There are two modes which begin with non-zero velocity and non-zero damping, only in the case of He–Ar with *c*). All other modes in each figure are pure diffusion modes starting with no initial velocity and non-zero damping. Some of these modes bifurcate into propagating waves with damping for larger *Kn*. As the mass ratios of the lighter gases in the gas mixtures of Ne–Ar, He–Ar and He–Xe is *μ*_{Ne}=0.3356, *μ*_{He}=0.0911 and *μ*_{He}=0.0296, we observe that the damping increases with increase in mass ratio of the lighter gas for any mole fraction. Moreover, for any mixture damping also increases with increase in the mole fraction of the lighter gas.

We have tested the value of damping *empirically*, that 2×G26 equations for a binary gas mixture of Maxwell molecules are linearly stable.

## 5. Heat flux between parallel plates

To validate the equations, we study the problem of heat flux in a binary gas mixture confined between two parallel plates having different temperatures and compare our results with those of [20] obtained from a direct discretization of the Boltzmann equation.

### (a) Problem description

Let us consider a binary gas mixture confined between two parallel plates at *x*=±*L*/2 and let the temperatures of the left plate (at *x*=−*L*/2) and the right plate (at *x*=*L*/2) be *T*_{°} so that the linear moment equations remain valid. The scales used for non-dimensionalization are same as those in (4.2).

For the problem under consideration, the system of linear-dimensionless 2×G26 can be obtained by substituting *b*. The relevant boundary conditions for the problem are (3.19)–(3.22) for the two species and an additional condition *v*_{n}=0 on the plates. The additional condition—*v*_{n}=0 on the plates—for the aforementioned problem means that the total number of the *α*- and *β*-gas molecules confined between the plates should remain equal to their given values, i.e.
** V**=

**0**. Conditions (3.19), (3.21) and (3.22) in linearized-dimensionless form read (

*γ*∈{

*α*,

*β*})

### (b) Method of solution

The system of moment equations for the problem along with the boundary conditions (5.1) and (5.2) is solved numerically using the finite difference method. The distance between the plates is discretized into 100 equispaced points and the method gives a solution within less than a second.

### (c) Results

We shall compare the results for binary gas mixtures from this method with those in [20] obtained using a direct discretization of the Boltzmann equation. Note that the rarefaction parameter _{1} by *χ*_{α}=*χ*_{β}=1 derived in this paper. Following [20], we calculate the dimensionless total heat flux between the plates

In tables 1–3, we present the dimensionless heat fluxes for Ne–Ar, He–Ar and He–Xe gas mixtures, respectively, for the same values of the rarefaction parameter *δ* used in [20] except *δ*=0. The results are quite close to [20] even for high Knudsen numbers, however, our results slightly overpredict the results of [20]. Furthermore, the decrease in the heat flux with increase in *δ* (or decrease in the Knudsen number) is verified also by our results for all the mixtures and any mass ratio. Possible reasons for the deviations are the use of an interaction potential based on Maxwell molecules which influences the form of the productions terms in the moment system, and the use of a simplified expression for the mean free path (4.12)_{2} implying the viscosity formula (5.3).

Figures 2 and 3 illustrate the variation of the total heat flux with change in mole fraction of the lighter gas in each mixture for *Kn*=0.0707 (*δ*=10) and *Kn*=0.7071 (*δ*=1), respectively. The small circles in figures 2 and 3 denote the data from [20] obtained using the realistic potential. Clearly, our results are in good agreement with those of [20]. The results are more close to those of [20] for small Knudsen number (*Kn*=0.0707). Similar to [20], our results also confirm that the maximum heat flux is observed around mole fraction

## 6. Conclusion and outlook

In this paper, the fully nonlinear *N*×G26 equations for a gas mixture of *N*-gases made up of Maxwell molecules have been derived and the associated nonlinear boundary conditions have been obtained via Maxwell’s accommodation model. The equations have been specialized to 2×G26 equations for a binary mixture and linearized in order to analyse the linear stability. By considering the plane wave solution, it has been shown that the 2×G26 equations for binary gas mixtures of Maxwell molecules are linearly stable in time for spatial perturbation of any wavenumber.

The linearized 2×G26 equations along with the linear boundary conditions from the Maxwell accommodation model have been solved numerically by finite difference method to study the heat flux problem for the binary gas mixture confined between two parallel plates at different temperatures. The results have been compared with the results of [20] obtained with realistic potential, and a good agreement between the two has been found. We conclude that the 2×G26 equations are a promising continuum model for the description of rarefied gas mixtures.

Using the findings of this paper, the method will be used to study more complicated problems, for which the direct solutions of the Boltzmann equation or DSMC solution are computationally too expensive. Furthermore, the complete analysis of the asymptotic order of magnitude will be performed and the regularized moment equations for binary gas mixtures will be developed for studying nonlinear processes in the future.

## Appendix A. Production terms for Maxwell molecules

In this appendix, we present the nonlinear production terms for the *α*-constituent corresponding to *u*^{2(α)}, using the abbreviations

Production term for stress tensor

Production term for heat flux

Production term for third moment

Production term for one trace of the fourth moment:

Production term for full trace of the fourth moment

- Received October 4, 2014.
- Accepted November 14, 2014.

- © 2014 The Author(s) Published by the Royal Society. All rights reserved.