Grad's 13-moment equations describe transport in mildly rarefied gases well, but are not properly embedded into nonequilibrium thermodynamics since they are not accompanied by a formulation of the second law. In this work, the Grad-13 equations are embedded into the framework of GENERIC (general equation for the nonequilibrium reversible–irreversible coupling), which demands additional contributions in the equations to guarantee thermodynamic structure. As GENERIC building blocks, we use a Poisson matrix for the basic convection behavior and antisymmetric friction matrices to correct for additional convective transport terms. The ensuing GENERIC-13 equations completely match the Grad-13 equations up to second-order terms in the Knudsen number and fulfill all thermodynamic requirements.

## I. INTRODUCTION

Closure is an important issue in developing simplified physical theories from more fundamental levels of description. Closure is often perceived as a problem but, more positively, it can also be seen as the crucial first step toward describing and understanding the essence of certain phenomena of interest by identifying appropriate variables.^{1} Closure is expected to be possible whenever a set of variables allows for an autonomous description of a problem of interest.

In the context of fluid dynamics, the densities of the locally conserved mass, momentum, and energy provide an extremely successful level of description (the five fields of hydrodynamics). However, when we are interested in time and length scales comparable to molecular scales, for example in rarefied gas dynamics, then we need more detailed levels of description.^{2} The Boltzmann equation for the single-particle velocity distribution provides an excellent, very detailed level of description for rarefied gas dynamics. However, for many applications, the Boltzmann equation is unnecessarily detailed and too complicated to solve. As standard hydrodynamics can be recovered from the five moments of the single-particle velocity distribution that correspond to the mass, momentum, and energy densities, it is tempting to use larger sets of moments to construct theories of fluid dynamics going beyond standard hydrodynamics.^{3} Then, the issue of moment closure arises from the fact that the Boltzmann equation is equivalent to an infinite hierarchy of coupled moment equations.^{2,3}

How can we recognize suitable sets of moments and valid closure procedures? An obvious criterion, which we refer to as *empirical adequacy*, requires that all relevant experimental results are reproduced by autonomous moment equations. Instead of experimental data, one might use the well-established Boltzmann equation or computer simulations of rarefied gases as a reference for empirical adequacy. However, there are further requirements from fundamental principles, such as the validity of a suitable version of the second law of thermodynamics for nonequilibrium systems. We refer to this type of criteria as *thermodynamic admissibility*, where we here rely on general equation for the nonequilibrium reversible–irreversible coupling (GENERIC)^{4–6} as a thermodynamic framework.

We consider both empirical adequacy and thermodynamic admissibility as necessary conditions for successful closure. Whereas none of these criteria by itself provides a sufficient condition, the combination of these two criteria offers a promising indicator for identifying successful closures. Of course, thermodynamic admissibility as a further constraint restricts the options for fitting experimental data. Yet, empirical adequacy should never be inferred from a perfect fit of a particular set of experimental data; it should be recognized in the overall performance of a model for all kinds of possible experiments and, in particular, in a complete coverage of all relevant qualitative phenomena.

In producing simplified theories of rarefied gas dynamics from the Boltzmann equation, it is natural to supplement the five hydrodynamic fields associated with the local conservation laws by five-second moments characterizing the momentum flux, and three third moments characterizing the energy flux. Whereas a lot is known about the empirical adequacy of the resulting 13-moment equations, starting with the work of Grad,^{2,7,8} very little is known about their thermodynamic admissibility. The construction of 13-moment equations that are consistent with the GENERIC framework of nonequilibrium thermodynamics is the topic of the present paper.

One key issue for thermodynamic admissibility is the existence of an entropy function. In systematic expansion techniques for the Boltzmann equation, such as Grad's moment method,^{7,8} one typically loses the positivity of the probability density for the velocity distribution and hence also Boltzmann's entropy involving the logarithm of that probability density. By using anisotropic Gaussian trial functions, this problem can be solved for 10-moment equations.^{9} Also for the linearized 13-moment equations an entropy has been found.^{10} While the framework of Rational Extended Thermodynamics^{3} aims at a thermodynamic structure for transport equations, approximations are used that are akin to the Grad method, so that the second law is met only approximately. The successful construction of an entropy for the nonlinear 13-moment equations requires not only suitable trial functions but also proper treatment of the general tensor structure of all the moments (see Ref. 11, refined through the discussion in Refs. 12 and 13).

Another, more subtle key issue is the Hamiltonian structure of the convection mechanism. Assuming the existence of a Poisson bracket for the 10-moment case, its necessary form has been derived on p. 309 of Ref. 6. However, the Jacobi identity for the resulting Poisson bracket had not been tested at the time. Later tests showed that, contrary to what was expected in Ref. 6, one actually does not arrive at a valid Poisson bracket. For the 10-moment case, this problem has recently been solved in Ref. 14 by introducing an irreversible contribution with Casimir symmetry to correct for velocity slip. In the present paper, the same strategy is applied to the 13-moment equations. The structure of the resulting thermodynamically consistent 13-moment equations is then compared to that of previous empirically adequate 13-moment equations lacking thermodynamic admissibility. Specifically, the GENERIC building blocks are designed such that they reproduce most elements of Grad's 13-moment equations within second-order accuracy, but contain additional higher order terms that give thermodynamic structure, which the Grad equations do not possess.

The structure of this paper is as follows: In Sec. II, we review Grad's original 13-moment equations, including a brief discussion of entropy and Knudsen scaling. Section III begins with a summary of the essentials of the GENERIC framework and then provides all the building blocks for a new, thermodynamically admissible version of the 13-moment equations. A detailed discussion of the resulting equations is offered in Sec. IV, including the entropy balance, the proper hydrodynamic limit, and second-order matching in the Knudsen scaling with Grad's original equations, which appears to not be completely possible. A short discussion and outlook concludes the paper. A class of positive-semidefinite matrices that is useful for the modeling of relaxation is introduced in Appendix A. Finally, Appendix B offers a brief discussion of Onsager–Casimir symmetry.

## II. GRAD'S 13-MOMENT EQUATIONS

### A. Variables

We shall not describe the derivation of the Grad 13-moment equations, but give only a summary of their features. Detailed accounts of derivation and properties of the equations are available in the literature.^{2,3,7,8} We will use tensor index notation with the Einstein summation throughout the paper. The 13 variables are moments of the distribution function $f(r,t,c)$, where ** r** and

*t*are the space and time variables, respectively, and

**denotes the microscopic velocities of particles. Specifically, we consider mass density**

*c**ρ*, momentum density $Mi=\rho vi$, where

*v*is center of mass velocity of the gas, and specific internal energy

_{i}*ε*, which are given as

Here, *m* is the mass of a gas particle, $\theta =kBmT$ is the temperature in specific energy units with the Boltzmann constant $kB$ and the thermodynamic temperature *T*, *v _{k}* is center of mass velocity, and

*C*=

_{k}*c*−

_{k}*v*is the so-called peculiar velocity that is the particle velocity in the co-moving frame. The trace-free and symmetric part of the second moment (indicated by angular brackets) defines the kinetic stress tensor as

_{k}and the heat flux is the trace of the third moment

The 13 variables of the original Grad-13 theory are

Within GENERIC, we shall use the equivalent set of variables

where the vector *w _{i}* will be related to the heat flux

*q*, and the temperature tensor Θ

_{i}_{ij}is defined through the pressure tensor $Pij$ as

which is conveniently split into its trace and trace-free parts. Here, $p=\rho \theta $ is the ideal gas pressure; we note that temperature is the trace of the temperature tensor

and that Θ_{ij} and *P _{ij}* are positive definite by definition.

The total energy density is the sum of internal and kinetic contributions,

The closure of the moment equations relies on an approximation of the distribution function that is entirely specified through the chosen variables ($u\xafA$ or *u _{A}*). Grad constructed his approximation such that it reproduces the 13 variables $u\xafA$, as the product of the equilibrium phase density—the Maxwellian—and a polynomial in

*C*that accounts for deviation from equilibrium

_{i}With this, all other moments become constitutive functions of the chosen variables.

By definition, a distribution function is positive. Due to its polynomial structure, the Grad approximation is not strictly positive, since the bracket in (9) becomes negative for large values of *C _{i}*. This behavior is compensated somewhat by the Maxwellian, which suppresses the polynomial for large

*C*. With this, the Grad approximation remains meaningful as long as stress

_{i}*σ*and heat flux

_{ik}*q*remain sufficiently low.

_{i}^{3}

### B. Grad-13 transport equations

The Grad-13 moment equations arise by multiplying the Boltzmann equation with the appropriate velocity polynomials and subsequent integration, using the Grad-13 moment distribution to determine all moments that appear.^{2,3,7,8,14} This yields the conservation laws for mass and momentum, and the balance law for internal energy (which together with momentum conservation yields conservation of total energy), in their well-known forms^{22}

and full balance laws for stress and heat flux, as

The production terms on the right-hand side are for Maxwell molecules, and *μ* is the shear viscosity of the gas.

### C. Entropy and H-theorem

Within kinetic gas theory, the Boltzmann entropy density, the non-convective entropy flux of the gas, and the entropy production rate are given by^{2,15,16}

with a scaling factor $y=e(mh)3$, where *e* and *h* are Euler's and Planck's constants, respectively, and the Boltzmann collision term $S(f,f)$. The corresponding balance of entropy, also known as H-theorem, reads

Evaluation with the Maxwellian gives the well-known equilibrium entropy density of the ideal gas

with the constant $\eta 0=[52+ln\u2009(2\pi 3m4h3)]$, while flux and production vanish (strictly speaking, the logarithmic terms should be combined because the argument of the ln-function must be dimensionless).

The entropy (13) is defined for arbitrary nonequilibrium states, which have positive distribution functions. With its loss of positivity, the logarithm of the Grad distribution $f|13$ in (9) does not exist for all *C _{i}* and the entropy density can only be determined by approximation for small stress and heat flux through Taylor expansion of $ln\u2009f|13$. To second order, entropy, entropy flux, and entropy production are found as

^{3}

Due to the approximations, insertion of these expressions into the entropy balance (14) does not produce a proper second law for the Grad-13 equations. Indeed, the above quadratic approximations in the nonequilibrium variables (*σ _{kl}*,

*q*) yield a suitable second law for the

_{k}*linearized*Grad-13 moment equations,

^{10,17}but not for the full nonlinear Grad-13 system (10)–(12).

We are not aware of a successful attempt to find a complete entropy inequality to accompany the Grad-13 equations (or any other Grad-type system). The loss of positivity of the distribution function, and the corresponding loss of entropy, removes proper physics when the equations are overstretched. Indeed, while close to equilibrium the Grad-13 equations are symmetric hyperbolic, far from equilibrium, that is, for large stresses and heat fluxes, they lose hyperbolicity and instead have complex eigenvalues,^{3} which leads to unphysical solutions. This implies that the Grad equations do not exclude solutions, which might violate the H-theorem. This does not render the equations meaningless from the perspective of empirical adequacy: The Grad system aims at approximating the Boltzmann equation through moment equations, and with this, also the H-theorem is only approximated. The equations describe gas rarefaction effects such as thermal stresses, or non-gradient heat fluxes, within limits of sufficiently small stress *σ _{kl}* and heat flux $qk$.

For meaningful solutions, the Grad equations should remain hyperbolic. Müller and Ruggeri^{3} demonstrate that hyperbolicity is lost for large (*σ _{kl}*,

*q*) and determine a hyperbolicity radius that estimates limits on (

_{k}*σ*,

_{kl}*q*). Proper thermodynamic structure—a valid H-theorem—and mathematical structure—hyperbolicity—are related; hence, the breakdown of hyperbolicity gives a strong hint that far from equilibrium Grad-13 will violate the second law.

_{k}### D. Knudsen scaling

Scaling arguments are helpful for a deeper understanding of the limitations of Grad-13. Indeed, considering scaling with the Knudsen number $Kn$ (the ratio of mean free path and desired length scale of resolution) as smallness parameter, it can be shown that Grad-13 is the proper model when one is interested in second-order accuracy $O(Kn2)$, while higher order accuracy must account for additional moments and their equations.^{2,18,19} We here present elements of the arguments as required to proceed.

The Grad-13 equations describe transport in mildly rarefied gases. They reduce to the classical Navier–Stokes equations by means of the Chapman–Enskog (CE) expansion, which relies on small mean free path of the gas molecules. Since the mean free path is proportional to the viscosity $\mu $, we can rescale the production terms for stress and heat flux as $\u22121\varpi p\mu \sigma ij,\u2212231\varpi p\mu qi$, with a formal smallness parameter ϖ that plays the role of the Knudsen number; the parameter can be set to unity after the expansion. For the CE expansion, stress and heat flux are written as a series

Next, the expansions are inserted into the respective transport equations, and matching of powers of ϖ yields the leading contributions for stress and heat flux as

When considering stress and heat flux only to first order, we obtain the well-known equations of Navier–Stokes and Fourier with shear viscosity *μ* and heat conductivity $154\mu $. Note that both transport coefficients have the same dimension due to our use of $\theta =kBmT$ as measure for temperature.

For the following, it will be important to note that both, *σ _{ij}* and

*q*, have no zeroth order—that is, equilibrium—contributions, but non-vanishing first-order contributions. Hence, we can say, both are of first order in the Knudsen number, that is, $O(\varpi )$. This understanding will now be employed to discuss the full Grad's 13-moment equations (11) and (12) for stress and heat flux.

_{i}We recognize the first-order characteristics formally by writing $\sigma ij=\varpi \sigma \xafij,\u2009qi=\varpi q\xafi,\u2009and\u2009\mu =\varpi \mu \xaf,$ where $\sigma \xafij,\u2009q\xafi,\u2009\mu \xaf$ are considered to be of order $O(1)$, and the factor $\varpi $ explicitly accounts for the actual order. Inserting this into (11) and (12), we obtain scaled transport equations for ($\sigma \xafij,\u2009q\xafk$) as

The Navier–Stokes and Fourier laws result from keeping only the first-order terms with factor $\varpi 1$ on the rhs, that is, from ignoring the terms with factors $\varpi 2$ and $\varpi 3$ on the lhs. The terms with power $\varpi 2$ on the lhs account for the second-order corrections to $(\sigma ij,qk)$. These are complete in the sense that consideration of additional moments will not add any terms at this order.

The underlined term in the heat flux balance is the only third-order term appearing in Grad-13. However, additional third-order terms appear in the equations when higher moments are added to the theory.^{18,19} With this, adding more moments to the list of variables results in higher order accuracy for stress and heat flux, and the predictions. We note that the conservation laws remain unchanged in their structure, and their accuracy is directly related to accuracy of stress and heat flux.

In the present contribution, we focus on the Grad-13 equations, which from the viewpoint of accuracy in terms of Knudsen number powers are limited to second order. Taking this at face value, the (underlined) third-order term should not play a role, so that it can be removed. On the other hand, second-order accuracy is maintained when third-order terms are *added* to the equations. This opens a strategy for reconfiguration of the Grad-13 equations by adding higher order terms such that the revised equations are accompanied by a properly formulated second law. Specifically, we shall aim to add terms that embed the revised Grad-13 equations into GENERIC.

## III. GENERIC FOR GRAD-13

### A. Basic structure of GENERIC

In the framework of GENERIC, irreversible dynamics is generated by the entropy, while reversible dynamics is generated by the energy, as in Hamiltonian dynamics. We give a summary of GENERIC elements that suits our desire to reformulate the Grad-13 equations such that they fit into this general framework.

In GENERIC, when the variables are denoted by *u _{A}*, the balance laws in local form read

where $L$ is the antisymmetric Poisson matrix, and $M$ is the friction matrix, and we assume a strictly local form of the total energy and entropy functionals

where *ϵ* and *η* are the local energy and entropy densities appearing as generators in the evolution equation (23).

Poisson matrix and friction matrix must fulfill a number of requirements:^{6}

The Poisson matrix must obey the degeneracy requirement

to ensure entropy conservation under reversible dynamics.

The antisymmetry of $LAB$ is defined through the Poisson bracket ${A,B}$,

where $A(uA)$ and $B(uB)$ are arbitrary functionals of the variables *u _{A}*. Finally, $LAB$ has to obey the Jacobi identity

The friction matrix $MAB$ appears in the transport equations linearly. Each dissipative mechanism contributes a friction matrix $MABk$ so that the overall matrix $MAB$ is the sum

All friction matrices must obey the degeneracy requirement

which ensures energy conservation in irreversible processes.

The individual $MABk$ must obey Onsager–Casimir relations; that is, they are either symmetric or antisymmetric. The dissipative bracket for $MABk$ is defined as

For a symmetric friction matrix $MABk,sym$, symmetry implies

where non-negative entropy production requires that

For an antisymmetric building block $\u2009MABk,anti$, we have the non-dissipative bracket

with

According to the discussion in Sec. 3.2 of Ref. 6, an antisymmetric contribution $MABk,anti$ corresponds to the irreversible evolution without entropy production associated with the well-known Casimir symmetry in linear irreversible thermodynamics.^{20} Note that a total friction matrix consisting of both symmetric and antisymmetric contributions associated with well-defined irreversible processes has no well-defined symmetry.

### B. GENERIC for 13 variables

While their entropy, entropy flux, and entropy production can be approximated to second order, see Eqs. (16)–(18), the Grad-13 equations do not possess an accompanying second law with positive production. Accordingly, it is not possible to find GENERIC building blocks that just reproduce the full Grad-13 equations (10)–(12).

As already indicated above, instead we aim at finding GENERIC building blocks that reproduce Grad-13 exactly to second order, that is, those with powers $(\varpi ,\varpi 2)$ in Eqs. (21) and (22), while allowing additional higher order terms, with powers $\varpi 3$ or higher. In the spirit of the Knudsen number scaling, these higher terms are expected to contribute very little to the solutions. In other words, we aim at correcting the Grad-13 equations by higher order contributions, such that the resulting equations possess a proper second law, while maintaining their second-order accuracy.

We construct the GENERIC version of the 13-moment equations not for the original variables $(\rho ,vi,\theta ,\sigma ij,qi)$ of the Grad model, but instead consider the related set (5) with

Compared to the third moments *q _{i}* defined in Eq. (3), the new variables

*w*have the major advantage that they possess vector transformation behavior under the general space deformations occurring in fluid dynamics, whereas

_{i}*q*defines a Cartesian vector only. This behavior is crucial for finding a Poisson bracket with a degenerate entropy, which, as a scalar, can only depend on a scalar combination of $\Theta jk$ and

_{i}*w*under general space transformations. Evaluation of (35) with the Grad-13 distribution gives

_{j}which to leading order reduces to

The energy density of the monatomic ideal gas is the sum of internal and kinetic energy (8), which for these variables assumes the form

with the gradient

The convex entropy density is postulated to be of the form^{11}

where $S\xaf(\phi )$ as a dimensionless scalar under general space transformations is assumed to be a dimensionless function of the scalar quantity

with negative slope

The desired consistency with the Grad-13 equations is achieved by comparing the second-order expansion of the entropy (40) in *σ _{ij}* and

*q*with Eq. (16) by assuming a linear dependence of $S\xaf(\phi )$ on $\phi $. Note that by its definition (6), the temperature tensor is positive definite, and that by definition the flux variable

_{i}*w*has unit $[\theta ]$, while the unit of heat flux

_{i}*q*is $[p\theta ]$.

_{i}The entropy gradient with respect to the variables *u _{A}* is

To complete the model, we require a Poisson matrix $LAB$ and a set of friction matrices $MABk$ that reproduce the various terms in the Grad-13 equations. These must be inserted in to the balance laws (23) for the variables *u _{A}*. We note that the balance of total energy and entropy balance results as

### C. $L$ matrix

A proper Poisson matrix that fulfills the degeneracy requirement (25) and the Jacobi identity (27) is given by (see Table 2.1 on p. 66 of Ref. 6 for lower convected tensors and vectors)

where indices in round brackets indicate symmetric tensors. It has been argued in Sec. 5(a) of Ref. 10 that lower convected derivatives appear naturally for moments of momentum. The contribution to the transport equations for the variables $uA=(\rho ,Mi,\Theta ij,wi)A$ is found as

where for compact notation, we have introduced the velocity gradient as

The contributions to energy and entropy balance are of divergence form

The degeneracy of the entropy with density (40) for the Poisson matrix (46) can be verified by a lengthy but straightforward calculation. The validity of the entropy degeneracy is the motivation for introducing variables with the proper general transformation behavior for scalars, vectors, and tensors.

### D. $Mcc$ matrix

The convection of moments obtained from the Boltzmann equation does not coincide with the lower convection behavior implemented by the Poisson matrix (46). In our previous contribution on the 10-moment equations,^{14} we hence proposed to introduce a “convective correction” term. Inspired by the thermodynamically consistent realization of Schowalter derivatives describing slip effects in rheology (see pp. 119f of Ref. 6), we proposed to implement the convective correction via an antisymmetric or Casimir-type contribution to the friction matrix, which is irreversible but non-dissipative (no entropy production; see Appendix B for a brief discussion of Onsager–Casimir symmetry). Incorporation of the convective correction as a reversible term is impossible because of a violation of the Jacobi identity, which is widely accepted as a criterion for reversible behavior (in the sense of being under “mechanistic control”). As an extension of the friction matrix $Mcc$ of our previous contribution^{14} from 10 to 13 moments by introducing suitable additional contributions, we propose

with the vorticity tensor

The contribution to the transport equations for the *u _{A}* reads

This matrix fulfills the degeneracy condition (29) that it is has no contribution to the energy balance, and also, it has no contribution to entropy balance

### E. $MABcross$ matrix

The next friction contribution is constructed to give heat flux in the energy balance as well as the cross-coupling between heat flux and stress tensor in their respective moment equations, which we know from the Grad-13 equations. It is of the form

where a natural ansatz for the nonzero matrix elements reads

Here, the dimensionless parameters *Z*_{1}, *Z*_{2,} and *c* can be scalar functions of the variable $\phi $. In these matrix elements, the expression in square brackets has been chosen such that it provides the most general closure approximation for the third-moment tensor $\u222bCiCjCkfdc$ [cf. Eq. (15) of Ref. 11]. The straightforward implementation of this general closure approximation through $Mcross$ completes the discussion of 13-moment equations with entropy initiated in Refs. 11 and 13, which are both based on simplified closures. The corresponding general form of the entropy current (3) of Ref. 13 implied by $Mcross$ is given by

or, more explicitly,

with

We recognize that the contributions with *Z*_{2} in (56) and (57) are of higher order; that is, they will not play a role in matching Grad-13 to second order in $Kn$. For simplicity of the following line of arguments, we set $Z2=0$, and we moreover assume that the parameters *Z*_{1} and *c* are constants.

By construction, $MABcross$ is energy degenerate, $MABcross\u2202\u03f5\u2202uB=0$, while the contribution to the GENERIC balance laws becomes

The contribution to total energy balance gives the negative divergence of the heat flux as

with the heat flux

To identify the coefficients *c*, *Z*_{1}, and $dS\xaf(\phi )/d\phi $, we study the leading order contributions and match these two terms in the Grad equations. To leading order the heat flux reduces to

and by comparison with (37), we identify

The leading contribution to stress balance is the trace-free part times $\rho $,

and comparison with (11) shows that *c *=* *1, which will be used from now on.

The leading contribution to the *w _{i}*-balance can be written as

where we used $\u2202\Theta kl\u22121\u2202xs=\u2212\Theta kn\u22121\Theta ml\u22121\u2202\Theta nm\u2202xs$. With the temperature tensor given by (6), this can be further reduced to give the leading contribution to the heat flux equation (multiplication with $p/2$)

Comparison with the Grad-13 equation (12) shows that setting $Z1=\u22124$ results in recovery of the first three terms, while the last term $(\u2212\sigma ik\u2202\theta \u2202xk)$ is extra, as will be discussed further below.

Thus, the matching at leading order yields

where the value for *c* implies symmetry of the third-moment tensor in its three indices and the values for $dS\xaf(\phi )/d\phi $ and $Q\xaf1$ coincide with previous findings.^{11} We note that $S\xaf$ itself must vanish in equilibrium, $S\xaf(\phi =0)=0$, so that the proper equilibrium entropy results from (40). We will continue with the simplest choice for $S\xaf$ that also gives the proper equilibrium value and set

While $MABcross$ is antisymmetric in the Casimir sense, it contributes to entropy flux (but not to entropy production). We find

with the non-convective entropy flux

which is consistent with (59).

### F. Dissipation terms

For the relaxation terms in the equations for stress and heat flux, we use

with $Aij\u200akl=\delta ik\u2009\Theta jl$. The matrix $MABrelax$ gives energy degeneracy, $MABrelax\u2202\u03f5\u2202uB=0$, is Onsager symmetric, and its positive semidefiniteness of the friction matrix follows from the results of Appendix A [see, in particular, the examples given in (A9) and (A10)].

The production terms are found as

where we already used (70). Indeed, the matrix was constructed such, that to second order the production terms reduce to the relaxation terms of the Grad equations

The non-linear second-order contribution $[23p\mu (1\theta \Theta \u27e8ij\u27e9)wj]$ to the *w _{i}*-balance is required for full matching, as will become clear further below.

This matrix contributes to entropy production as

where for compact notation, we introduced the abbreviation

The entropy production is non-negative due to the positive definiteness of the matrix. This can also be recognized by considering that $\Theta ik$ is positive definite with trace $\Theta kk=3\theta $, which implies for the trace of the inverse that $\Theta ss\u22121\u22653/\theta $.

## IV. GENERIC-13 EQUATIONS

### A. Collected GENERIC-13 equations

We use the above results to assemble the transport equations for the GENERIC variables $uA=(\rho ,Mi,\Theta ij,wi)A$ as

where we use the abbreviations

#### 1. Mass balance

The mass balance has only a contribution from $LAB$ and assumes the usual form

#### 2. Momentum balance

Also the momentum balance assumes its usual form, with the pressure tensor $\rho \Theta ij$,

#### 3. Temperature tensor balance

The balance for the temperature tensor Θ_{ij} includes the energy balance as its trace, and the balance for anisotropic stress tensor $\sigma ij=\rho \Theta \u27e8ij\u27e9$ as its trace-free part; it reads

We note that this equation is equivalent to a conservation law for the full second-moment $Fij=\rho (\Theta ij+vivj)$ with the same production term

#### 4. w-balance

The balance for the vector *w _{i}* reads

This equation is equivalent to the balance for the heat flux, where both vectors are related in a nonlinear fashion as

#### 5. Entropy balance

Finally, we take a look at the entropy balance, which has the general form

where entropy, entropy flux, and entropy production are given by

### B. Hydrodynamic limit

The first-order CE expansion of the GENERIC-13 transport equations (83), (85), and (86) results from inserting equilibrium values for trace-free temperature tensor $\Theta \u27e8ij\u27e9$ and vector *w _{i}* on the left hand side of their transport equations. With the equilibrium values $\Theta ij=\theta \delta ij,\u2009\Theta ij\u22121=\delta ij/\theta ,\u2009\Theta kk\u22121=3/\theta $, this yields the Navier–Stokes and Fourier laws as

In this limit, we can also replace $qk=12pwk$, hence stress and heat flux agree with (20),

Then, entropy, entropy flux, and entropy production reduce to the well-known relations

### C. Second-order matching

Due to the nonlinear relation between heat flux *q _{i}* and vector

*w*, it is quite a cumbersome task to rewrite the GENERIC equations for the Grad variables

_{i}Indeed, for further understanding of agreement and differences between the Grad-13 equations and the GENERIC-13 system, it suffices to consider the GENERIC equations within the appropriate order of magnitude of the various terms. Using scaling based on the Knudsen number, we identified both stress $\sigma ij$ and heat flux *q _{k}* as first-order quantities in terms of the Knudsen number, or rather the formal parameter ϖ. Hence, the task is to use the same kind of order expansion on the GENERIC equations, while introducing the Grad-13 variables. For the evaluation, one must use that temperature tensor is related to the stress tensor as

so that within second order, that is, considering up to quadratic terms in the stress, the inverse becomes

Similarly, the vector *w _{i}* can be approximated to second order as

Careful evaluation yields that the conservation laws for mass, momentum, and energy for both sets of equations are identical (this is already guaranteed by the full nonlinear relation between *q _{i}* and

*w*).

_{i}Differences between Grad-13 and GENERIC in the equations for stress $\sigma ij=\rho \Theta \u27e8ij\u27e9$ turn out to be of higher order; that is, to the required second order, the GENERIC balance (83) agrees with the Grad-13 balance.

However, differences at the leading order occur at first in the heat flux equation. After replacing $(\Theta ij,wi)$ by $(\sigma ij,qi)$, multiplication with the inverse of $15\rho (32\theta \delta kl+\Theta kl)$, and removing all terms that are of third or higher orders, we find

This equation is written such that all terms of the Grad-13 heat flux balance (12) appear explicitly, and additional terms are underlined, with the difference to Grad-13 given by

Using the NSF relations (20), we can rewrite this difference within the required second-order accuracy as an algebraic expression

which vanishes within the required order. Indeed, we have constructed the relaxation matrix $MABrelax$ such that the non-linear relaxation term compensates the two gradient terms. With this, also the Grad-13 heat flux equation is included in the GENERIC-13 equations to second order.

## V. DISCUSSION AND OUTLOOK

In this paper, we have developed a set of 13-moment equations that is thermodynamically admissible in the sense of possessing the full GENERIC structure. On the one hand, these new 13-moment equations generalize previous work^{11–13} that developed the ideas for establishing an entropy, which is an important necessary but not sufficient condition for thermodynamic admissibility. On the other hand, we have generalized the full thermodynamic consistency of recently established equations for 10 moments^{14} to the case of 13 moments.

We have completely aligned the well-known 13-moment equations of Grad and the GENERIC formalism of nonequilibrium thermodynamics. Owing to unavoidable approximations in their derivation, the original Grad-13 equations are not accompanied by a proper entropy inequality and exhibit entropic behavior only within the limits of approximation, that is, for sufficiently small deviations from the equilibrium state. Outside the limits of approximation, Grad's equations become problematic (e.g., breakdown of hyperbolicity). However, within the limits of approximation, the Grad-13 equations describe the physics of a mildly rarefied gas quite well.

GENERIC, on the other hand, guarantees thermodynamic structure for all processes. The GENERIC framework requires Poisson and friction matrices as building blocks, which describe the reversible and irreversible transport behavior, respectively. Moreover, the different physical contributions to friction matrices are either of symmetric Onsager type or antisymmetric Casimir type.

In our effort to align both approaches, we aim at constructing a set of 13 equations within GENERIC, with full entropic structure, that agrees to Grad-13 within second-order accuracy, based on Knudsen number scaling. The resulting equations have additional terms of higher order that provide the desired structure. The matching is fully successful.

An important role in the construction of the GENERIC-13 equations is played by friction matrices of Casimir type, on which we rely heavily. Indeed, we use a well-established Poisson matrix from the theory of complex fluids for the basic convection behavior and introduce antisymmetric friction matrices to correct for additional transport terms that occur in the Grad-13 equations. This somewhat nonstandard and pragmatic approach to GENERIC allowed us to almost completely match GENERIC to Grad-13 (to second order), while guaranteeing the proper structure of GENERIC, including the validity of the Jacobi identity for reversible transport.

Here, it must be noted that the Jacobi identity is the biggest obstacle in the development of equations that comply with GENERIC since it is tedious to prove and even harder to achieve in the construction of a Poisson matrix. Friction matrices are simpler—but not trivial—to construct and have fewer restrictions. Hence, their use provides a welcome degree of flexibility.

In this paper, we focus entirely on the thermodynamic admissibility of the 13-moment equations. The empirical adequacy of the resulting admissible equations remains to be investigated by studying their solutions in a variety of flow problems, such as 1D shock structures, and the interesting question of hyperbolicity for the GENERIC-13 equations.

## ACKNOWLEDGMENTS

The topic of the present work lies in the intersection of the molecular theory of gases and transport phenomena, two fields in which R. Byron Bird has broken new ground in his first book^{21} (with J. O. Hirschfelder and C. F. Curtiss, 1954) and in his most popular book^{22} (with W. E. Stewart and E. N. Lightfoot, 1960), also known as BSL or the red bible. Moreover, moment equations are at the origin of modeling viscoelastic fluids (the relaxation time in the famous Maxwell model of rheology is nothing but the collision timescale in a gas); Bob Bird has also shaped this field as the leading author of a two-volume book on the dynamics of polymeric liquids.^{23,24} The general phase-space kinetic theory in Part VII of Ref. 24 is a powerful framework for deriving properly structured theories (including a fluctuation-dissipation theorem incorporated through the assumption of “equilibration in momentum space”), and the most carefully motivated “major approximations” in Figs. 18.0–1 and 19.0–1 of that book simply provides closure. In our own work, we feel the strong direct and indirect influence of Bob Bird and his inspiring books. We will always remember him with admiration and gratitude.

H.S. gratefully acknowledges support from the National Sciences and Engineering Research Council (NSERC), Grant No. RGPIN-2016–03679. H.C.O. is grateful for the freedom to do essential parts of this work during a sabbatical at the University of Chicago.

## AUTHOR DECLARATIONS

### Conflict of Interest

The authors have no conflicts to disclose.

## DATA AVAILABILITY

Data sharing is not applicable to this article as no new data were created or analyzed in this study.

### APPENDIX A: A POSITIVE-SEMIDEFINITE MATRIX

Let Θ_{ij} be a positive-semidefinite, regular matrix, and *w _{i}* a vector. We then define the auxiliary quantities

where *x*, *y*, and $z1\u2026z4$ are real parameters, and we consider the matrix

where the prefactor $N$ is positive. The trace-free symmetrization $\u27e8kl\u27e9$ is required to obtain energy conservation. As a sum of dyadics, this matrix *M* must be positive semidefinite. By inserting the definitions (A1), we obtain the alternative representation

with

If we are interested only in the first-order Knudsen expansion and use $x=\theta \u2212(z1+z2)$ and $y=y\xaf\u2009\theta \u2212(1+z3+z4)$, then we get the simplified representation

with

If we consider the corresponding relaxation matrix

we find

##### 1. Example

We choose $z1=z3=0,\u2009z2=1/2,\u2009z4=\u22121/2,\u2009Aij\u200akl=\delta ik\u2009\Theta jl$, *x *=* *1, *y *=* *1/2 and obtain

which leads to

### APPENDIX B: BARE AND DRESSED ONSAGER–CASIMIR SYMMETRY

The brief discussion of Onsager–Casimir symmetry properties of friction matrices in this appendix is based on Sec. 3.2.1 of Ref. 6. As the variables $uA=(\rho ,Mi,\Theta ij,wi)A$ are alternatingly even and odd moments of velocity, the corresponding bare symmetry properties of the friction matrix under time reversal are given by the checkerboard pattern

Since the bare symmetry depends only on the time-reversal behavior of the independent variables, the checkerboard pattern (B1) actually characterizes the bare symmetry for all contributions to the friction matrix.

As the matrix elements of the friction matrix are not just constants (as they would be in linear irreversible thermodynamics), but they themselves can change sign under time reversal if they contain the third-moment vector $wi$ or the vorticity tensor *ω _{ij}*, the structure of $Mcc$ proposed in (51) leads to the dressed symmetry properties

As zeros are consistent with both symmetry and antisymmetry, the entire matrix $Mcc$ is consistent with antisymmetry, that is, of the Casimir type.

For the cross-coupling between heat flux and stress tensor through friction matrix $Mcross$ in (55), the bare checkerboard pattern is dressed to become

which is also consistent with Casimir antisymmetry. On the other hand, for the relaxation term $Mrelax$ in (73), the dressed symmetry is given by

which is consistent with Onsager symmetry.