\documentclass[12pt]{article}
\usepackage{graphicx}
\usepackage{amsmath}
\usepackage{amsfonts}
\usepackage{amssymb}
%\usepackage{refcheck}
%TCIDATA{OutputFilter=latex2.dll}
%TCIDATA{CSTFile=LaTeX article (bright).cst}
%TCIDATA{Created=Tue May 23 10:02:18 2000}
%TCIDATA{LastRevised=Fri May 26 12:44:50 2000}
%TCIDATA{<META NAME="GraphicsSave" CONTENT="32">}
%TCIDATA{<META NAME="DocumentShell" CONTENT="Journal Articles\Standard LaTeX Article">}
%TCIDATA{Language=American English}
\newtheorem{theorem}{Theorem}
\newtheorem{acknowledgement}[theorem]{Acknowledgement}
\newtheorem{algorithm}[theorem]{Algorithm}
\newtheorem{axiom}[theorem]{Axiom}
\newtheorem{case}[theorem]{Case}
\newtheorem{claim}[theorem]{Claim}
\newtheorem{conclusion}[theorem]{Conclusion}
\newtheorem{condition}[theorem]{Condition}
\newtheorem{conjecture}[theorem]{Conjecture}
\newtheorem{corollary}[theorem]{Corollary}
\newtheorem{criterion}[theorem]{Criterion}
\newtheorem{definition}[theorem]{Definition}
\newtheorem{example}[theorem]{Example}
\newtheorem{exercise}[theorem]{Exercise}
\newtheorem{lemma}[theorem]{Lemma}
\newtheorem{notation}[theorem]{Notation}
\newtheorem{problem}[theorem]{Problem}
\newtheorem{proposition}[theorem]{Proposition}
\newtheorem{remark}{Remark}
\newtheorem{solution}[theorem]{Solution}
\newtheorem{summary}[theorem]{Summary}
\newcommand{\ba}{\boldsymbol{a}}
\newcommand{\bA}{\boldsymbol{A}}
\newcommand{\bb}{\boldsymbol{b}}
\newcommand{\bc}{\boldsymbol{c}}
\newcommand{\bd}{\boldsymbol{d}}
\newcommand{\bJ}{\boldsymbol{J}}
\newcommand{\bM}{\boldsymbol{M}}
\newcommand{\bN}{\boldsymbol{N}}
\newcommand{\bQ}{\boldsymbol{Q}}
\newcommand{\bR}{\boldsymbol{R}}
\newcommand{\bS}{\boldsymbol{S}}
\newcommand{\bv}{\boldsymbol{v}}
%\newcommand{\bx}{\boldsymbol{x}}
\newcommand{\ds}{\displaystyle}
\newcommand{\C}{\mathbb{C}} %numeri complessi
%\newcommand{\CA}{\mathcal{A}}
%\newcommand{\CB}{\mathcal{B}}
\newcommand{\CF}{\mathcal{F}}
\newcommand{\CH}{\mathcal{H}}
%\newcommand{\CK}{\mathcal{K}}
\newcommand{\CL}{\mathcal{L}}
\newcommand{\CW}{\mathcal{W}}
\newcommand{\CX}{\mathcal{X}}
\newcommand{\CY}{\mathcal{Y}}
\newcommand{\CZ}{\mathcal{Z}}
\newcommand{\D}{\mathcal{D}} %Il comando \CD e` destinato ai diagrammi
\newcommand{\deff}{\stackrel{\text{\tiny def}}{=}}
\newcommand{\N}{\mathbb{N}} %numeri naturali
%\newcommand{\Q}{\mathbb{Q}} %numeri razionali
\newcommand{\R}{\mathbb{R}} %numeri reali
%\newcommand{\T}{\mathbb{T}} %cerchio unitario nel piano complesso
\newcommand{\Z}{\mathbb{Z}} %numeri interi
\newcommand{\za}{\alpha}
\newcommand{\zb}{\beta}
\newcommand{\zg}{\gamma}
\newcommand{\zG}{\Gamma}
\newcommand{\zd}{\delta}
%\newcommand{\zD}{\Delta}
\newcommand{\ze}{\varepsilon}
%\newcommand{\zz}{\zeta}
\newcommand{\tht}{\theta}
%\newcommand{\zth}{\vartheta}
\newcommand{\zk}{\kappa}
\newcommand{\zl}{\lambda}
%\newcommand{\zL}{\Lambda}
\newcommand{\zp}{\psi}
%\newcommand{\zr}{\rho}
\newcommand{\zs}{\sigma}
\newcommand{\zS}{\Sigma}
%\newcommand{\zt}{\tau}
%\newcommand{\zf}{\varphi}
%\newcommand{\zF}{\Phi}
%\newcommand{\zq}{\chi}
\newcommand{\zo}{\omega}
\newcommand{\zO}{\Omega}
\newcommand{\ti}{\tilde}
\newcommand{\zvh}{\varrho}


\newenvironment{proof}[1][Proof]{\textbf{#1.} }{\ \rule{0.5em}{0.5em}}
\setlength{\textwidth}{16.6cm} \setlength{\textheight}{22.6cm}
\addtolength{\hoffset}{-1.9cm}
\begin{document}
	
	\title{A model for polyatomic gases with hyperboliciy and H-theorem satisfied up to whatever order}
	
		\author{S. Pennisi\\
	\\
	\footnotesize 
	\hspace{-1.8 cm} \footnotesize 	${}^2$	Dipartimento di Matematica ed Informatica, Universit\'{a} di Cagliari, 	Cagliari, Italy; \\
	\footnotesize spennisi@unica.it} 
	
	\date{}
	\maketitle \vspace{0.5 cm}
	
	

	
	
	\abstract{It is well known in literature a relativistic model for polyatomic gases, with an arbitrary but fixed number of moments, whose balance equations have the symmetric hyperbolic form in their left hand sides because the tensors derivated with respect to $x^\alpha$ are gradients of a 4-potential. Here the symmetric form  and a 4-potential is obtained also for their right hand sides, i.e., the production terms. Moreover, this fact will allow to prove the H-theorem up to whatever order, while in other articles present in literature this result was achieved only up to second order with respect to equilibrium. These results are here obtained whether  following the Eckart approach or the Landau-Lifshitz one.}
	
\section{Introduction}
The most diffused classical approaches to model a dissipative relativistic gas are those of Eckart \cite{Eck} and Landau-Lifshitz \cite{Lif}. In \cite{CrisPen} it has been proved that the two approaches are equivalent; they differ only for a different definition of the deviations from equilibrium and consequent expansion performed to obtain a linear model. We will prove here that, if these expansions would be done up to whatever order, the two approaches would give the same result. These expansions are used to do a transformation of the independent variables from apparently anonimous Lagrange Multipliers
to variables which have an immediate physical meaning. However this is just a stylistic issue, as in elementary geometry the algebraic expression of a surface isn't more significant of its parametric expression which uses "anonimous" parameters. So, before to convert the Lagrange Multipliers to the other variables, the two approaches are the same.
The  Landau-Lifshitz approach has the advantage to obtain immediately a zero production of mass and of momentum-energy; here we will see that it satisfies  the H-theorem up to whatever order. We will see that the same result can be obtained also with the Eckart  approach. 
Regarding the hyperboliciy, it has already been proved that it holds up tp whatever order if the Lagrange Multipliers are used as independent variables (\cite{CrisPen2}, \cite{Pen}). Whatever restriction on the so called hyperbolicity region is due to the low degree of approximation which is used for the above unnecessary change of variables. As a confirmation of this fact we see that Profs. Brini and Ruggeri implemented the articles  \cite{Bri}, \cite{R-T} by considering a better approximation (The second order one) in
\cite{Bri2} and obtained a bigger hyperbolicity zone. This zone would cover all the set of possible values of independent variables if no approximation was used; this is not possible for calculations difficulties, but we cannot expect nature to bow to our mathematical difficulties. \\
Other doubts about the validity of Extended Thermodinamics with many moments arose following Struchtrup's articles \cite{Struc1}, \cite{Struc2} and similar, which studied a transition to Ordinary Thermodynamics and found some weak points. But Extended Thermodynamics arose exactly to overcome the problems of Ordinary Thermodynamics, such as the parabolic equations and the propagation waves with infinite speeds; so it makes no sense to test a better theory through the worse theory that she passed. Moreover, Struchtrup too uses approximations around equilibrium and it makes no sense to increase the number of moments without increasing the order of approximations around equilibrium.  Finally, the methods used for the transition to Ordinary Thermodynamic are mathematically well defined, but there is nothing physical which ensures that the result of these methods really lead to Ordinary Thermodynamics; this doubt is reinforced by the fact the results for  bulk viscosity,  heat conductivity and  shear viscosity  are different depending on wheter the Chapman-Enskog Method or the Maxwellian Iteration are used (in this latter case it even depends on the number $N$ of moments being used). This is proved in \cite{DP}. So we don't take into account these non-existing weakness in the sequel. \\
The starting point to obtain the balance equations is the Boltzmann–Chernikov equation
\begin{equation}\label{Bolt}
p^\alpha \, \partial_\alpha \, f =  Q \, ,
\end{equation}
where  $f$ is the distribution function and $Q$ the production term. The expression of the first one of these was found in \cite{Annals2017}, \cite{Pen}, \cite{ACPR} and reads
\begin{equation}\label{f}
f= e^{-1 \, - \, \frac{\chi}{k_B}} \, , \quad \chi= \sum_{n=0}^{N} \frac{1}{m^{n-1}} \, \lambda_{\alpha_1 \alpha_2 \cdots \alpha_n} p^{\alpha_1} p^{\alpha_2} \cdots p^{\alpha_n} \,\left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right)^n \, , 
\end{equation}
where $k_B$ is the Boltzmann constant, $m$ is the particle mass, $N$ denotes the number of moments which is used, $\lambda_{\alpha_1 \alpha_2 \cdots \alpha_n}$ are the Lagrange multipliers, $p^\alpha$ is the 4-momentum of a particle (satisfying the relation $p^\alpha \, p_\alpha= m^2 c^2$), $\mathcal{I}$ is the internal energy of a particle due to its internal modes (rotations and vibrations). \\
Regarding the production term $Q$, if we adopt the Landau-Liftschiz approach, it has the form 
\begin{equation}\label{prod}
Q = - \, \frac{1}{c^2 \tau} \, U_{L \alpha} \, p^\alpha \left( f \, - \, f_{eq.} \right) \, , 
\end{equation}
where $\tau$ is a relaxation time and $U_{L \alpha}$ is the Landau-Lifshitz
4-velocity whose physical meaning can be found in \cite{19,20}. We will see here that it generates an expression valid also for the Eckart  approach. In particular, we will see that, with the Landau-Liftschiz approach, the balance equations take the form
\begin{align}\label{last}
\boxed{\partial_\alpha \, \frac{\partial \, h'^\alpha}{\partial \, \lambda} = 0 \, , \quad \partial_\alpha \, \frac{\partial \, h'^\alpha}{\partial \, \lambda_{\beta}} = 0 \, , \quad \partial_\alpha \, \frac{\partial \, h'^\alpha}{\partial \, \lambda_{\beta_1 \cdots \beta_n}} = \frac{U_\alpha }{c \,   \tau} \, \frac{\partial \, Q^\alpha}{\partial \, \lambda_{\beta_1 \cdots \beta_n}} \quad \mbox{for} \quad n \geq 2} \, ,
\end{align}
where the gradient of a 4-potential is present both in the left hand sides (here the 4-potential is $h'^\alpha$) than in the right hand sides (here the 4-potential is $Q^\alpha$). A less elegant form, but with the same properties, will be obtained by following the Eckart approach and it reads:
\begin{align}\label{last2}
\partial_\alpha \, \frac{\partial \, h'^\alpha}{\partial \, \lambda_{\beta_1 \cdots \beta_n}} = \frac{U_\alpha }{c \,   \tau} \, \left( \frac{\partial \, \tilde{Q}^\alpha}{\partial \, \lambda} \frac{\partial \, g}{\partial \, \lambda_{\beta_1 \cdots \beta_n}} \, + \, 
\frac{\partial \, \tilde{Q}^\alpha}{\partial \, \lambda_\mu} 
\frac{\partial \, g_\mu}{\partial \, \lambda_{\beta_1 \cdots \beta_n}} + 
\frac{\partial \, \tilde{Q}^\alpha}{\partial \, \lambda_{\beta_1 \cdots \beta_n}} \right)\quad \mbox{for} \quad n \geq 2 \, ,
\end{align}
where the 4-potential in the right hand side is  $\tilde{Q}^\alpha$ while $g$ and $g_\mu$ are known functions which will be presented below.
The plan of this article is the following: In the next section we will see the expressions at equilibrium which is the same for both approaches. In sect. 3 we will see the balance equations outside equilibrium by separating the two approaches into 2 subsections; the above expressions \eqref{last} and \eqref{last2} will be found and the H-theorem will be proved for them up to whatever order with respect equilibrium (obviously, refraining from making approximations). \\
In sect. 4 a transformation wil be found which allows to obtain the variables of the Landau-Liftschiz approach in terms of those in Eckart approach. Obviously, the same transformation can be done in the inverse direction but we will refrain to do it for the sake of brevity. Approximations will be used, up to first order with respect to equilibrium, only in its subsection 1 to see how previously result in literature can be recovered from the present one. In its subsection 2, approximations will be used up to second order only to show that the results previously obtained in literature don't hold at any order.
 
 
 

 
 
 
 
 
 
 
 
 



\section{The model for relativistic polyatomic gas at equilibrium}
At equilibrium the Landau-Liftschiz and the Eckart approach give the same expressions. 
In particular we have the balance equations
\begin{align*}
\begin{split}
& \partial_\alpha \, V_E^\alpha = 0 \, , \quad \partial_\alpha \, T_E^{\alpha \beta} = 0  \, , \quad \mbox{where} \quad V_E^\alpha = \rho \, U^\alpha \, , \quad T_E^{\alpha \beta} =
	\frac{e}{c^2} \, U^\alpha U^\beta \, + \, p \, h^{\alpha \beta} \, , \\
& U_{\alpha} \, U^{\alpha} =c^2 \, , \,	h^{\alpha \beta} =- g^{\alpha \beta} +  \frac{U^\alpha U^\beta}{c^2} \, ( \mbox{The projector into the subspace orthogonal to} \, U^\alpha ) \, . 
\end{split}
\end{align*}
If we want to find  an approach which holds for whatever type of gas, we can introduce  the 4-potential
\begin{align*}
h'^\alpha = - \, 4 \, \pi \, m^3 c^5 \, h_0 \left( \lambda^E \, , \, \gamma \right) \frac{\lambda^\alpha_E}{\gamma} \, , \quad \mbox{where} \quad \gamma = \frac{m \, c}{k_B} \sqrt{\lambda^\alpha_E \, \lambda_\alpha^E} \, \rightarrow \, \lambda^\alpha_E \, \lambda_\alpha^E = \left( \frac{k_B \gamma}{m \, c} \right)^2 \, . 
\end{align*}
(The constant coefficients have been introduced for an easier comparison with expressions previously known in literature). It follows
\begin{align*}
V^\alpha_E = \frac{\partial \, h'^\alpha}{\partial \, \lambda_E} = - \, 4 \, \pi \, m^3 c^5 \, \frac{\partial \, h_0}{\partial \, \lambda_E} \, \frac{\lambda_E^\alpha}{\gamma} \, , \, T^{\alpha \beta}_E = \frac{\partial \, h'^\alpha}{\partial \, \lambda^E_\beta} = - \, 4 \, \pi \, m^3 c^5 \, \left( \frac{h_0}{\gamma} \, h^{\alpha \beta} \, + \, \frac{\partial \, h_0}{\partial \, \gamma} \, \frac{\lambda_E^\alpha \, \lambda_E^\beta}{\lambda_E^\mu \, \lambda^E_\mu} \right) \, . 
\end{align*}
As consequence of these results we have 
\begin{align}\label{ro}
\rho = - \, 4 \, k_B \, \pi \, m^2 c^3 \, \frac{\partial \, h_0}{\partial \, \lambda_E} \, , \, U^\alpha = \frac{m \, c^2}{k_B} \,  \frac{\lambda_E^\alpha}{\gamma} \, , \, p= 4 \, \pi \, m^3 c^5 \, \frac{h_0}{\gamma} \, , \, e= - \,  4 \, \pi \, m^3 c^5 \, \frac{\partial \, h_0}{\partial \, \gamma}  \, .
\end{align}
If we know the constitutive function $e=e(\rho \, , \, \gamma)$, then \eqref{ro}$_4$ becomes 
\begin{align*}
e\left( - \, 4 \, k_B \, \pi \, m^2 c^3 \, \frac{\partial \, h_0}{\partial \, \lambda_E} \, , \, \gamma \right) = - \,  4 \, \pi \, m^3 c^5 \, \frac{\partial \, h_0}{\partial \, \gamma} \, ,
\end{align*}
which is a differential equation from which to deduce $h_0$. For example, if $\frac{e}{\rho \, c^2} = \epsilon (\gamma)$ and we define $\eta (\gamma)$  from $\eta ' (\gamma) =\epsilon (\gamma)$, this equation becomes $\frac{\partial \, h_0}{\partial \, \gamma} \, - \, \frac{k_B}{m} \, \eta' \,  \frac{\partial \, h_0}{\partial \, \lambda_E}=0$. By considering $h_0$ a composite function of $H_0(X \, , \, Y)$ and of $X=\lambda^E + \frac{k_B}{m} \, \eta (\gamma)$, $Y= \gamma$, this differential equation becomes $\frac{\partial \, H_0}{\partial \, Y}=0$ so that  the general solution is \\
$h_0= H_0 \left( \lambda^E + \frac{k_B}{m} \, \eta (\gamma) \right)$ for whatever single variable function $H_0$. A further subcase is that when
\begin{align*}
\epsilon = \frac{e}{\rho \, c^2} = \frac{\int_{0}^{+ \, \infty}  J_{22} \left( \gamma* \right) \left( 1 + \frac{\mathcal{I}}{m \, c^2} \right) \varphi (\mathcal{I}) \, d \, \mathcal{I} }{\int_{0}^{+ \, \infty}  J_{21} \left( \gamma* \right)  \varphi (\mathcal{I}) \, d \, \mathcal{I} } \rightarrow \eta(\gamma)= - \, \ln \, \int_{0}^{+ \, \infty}  J_{21} \left( \gamma* \right)  \varphi (\mathcal{I}) \, d \, \mathcal{I} \, . 
\end{align*}
By taking $H_0(X)= e^{-1-\frac{m \, X}{k_B}}$ we find $h_0= e^{-1-\frac{m \, \lambda^E}{k_B}}  \int_{0}^{+ \, \infty}  J_{21} \left( \gamma* \right)  \varphi (\mathcal{I}) \, d \, \mathcal{I}$ and \eqref{ro}$_{1,3}$ become 
\begin{align*}
\rho = 4  \, \pi \, m^3 c^3 \, e^{-1-\frac{m \, \lambda^E}{k_B}}  \int_{0}^{+ \, \infty}  J_{21} \left( \gamma* \right)  \varphi (\mathcal{I}) \, d \, \mathcal{I} \, , \quad \frac{p}{\rho} = \frac{c^2}{\gamma} \, , 
\end{align*}
as in eqs. (26) and (41) of \cite{Annals2017}. For the sake of simplicity, we will use in the sequel this simpler expression.

\section{The dissipative case and the production term}

\subsection{The Eckart approach}
In \cite{ACPR}, by multiplying eq. \eqref{Bolt} with  $\frac{c}{m^{n-1}} \, p^{\alpha_1} p^{\alpha_2} \cdots p^{\alpha_n} \,\left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right)^n$ and integrating in $d \, \mathcal{I} \, d \, \vec{P}$, the balance equations for this case have been found 
\begin{align}\label{bal}
\partial_\alpha \, A^{\alpha \alpha_1 \cdots \alpha_n} = I^{\alpha_1 \cdots \alpha_n} \, , \quad \mbox{for} \quad n=0 \, , \, 1 \, , \, \cdots \, , \, N \, .
\end{align}
where 
\begin{align}\label{dis}
\begin{split}
& A^{\alpha \alpha_1 \cdots \alpha_n} =  \frac{\partial \, h'^\alpha}{\partial \, \lambda_{\alpha_1 \cdots \alpha_n}} \, , \quad h'^\alpha = - k_B \, c \int_{\Re^3} \int_{0}^{+ \infty} \, f  \,  p^\alpha \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} \, , \\ 
& I^{\alpha_1 \cdots \alpha_n}=  \frac{c}{m^{n-1}} \int_{\Re^3} \int_{0}^{+ \infty} \, Q  \,   p^{\alpha_1} p^{\alpha_2} \cdots p^{\alpha_n} \,\left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right)^n \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} \, .
\end{split}
\end{align}
Hence the left hand side of the balance equation \eqref{bal} takes the elegant expression $\partial_\alpha \, \frac{\partial \, h'^\alpha}{\partial \, \lambda_{\alpha_1 \cdots \alpha_n}}$. We will show now that also for the right hand side we can obtain a similar expression $I^{\alpha_1 \cdots \alpha_n}=\frac{U_\alpha}{c \,   \tau} \, \frac{\partial \, Q^\alpha}{\partial \, \lambda_{\beta_1 \cdots \beta_n}} $ with $Q^\alpha$ which will be found later and $\tau$ a relaxation time. To this end we need to know the production term $Q$. Its expression proposed in \cite{Giap} was an approximated one and, in fact, it gave an entropy production $\sigma$ which was non negative only up to second order with respect to equilibrium. To obtain $\sigma \geq 0$ up to whatever order, we note firstly that, for every expressions of the functions $g \left( \lambda_{\alpha_1 \alpha_2 \cdots \alpha_n} \right)$, $g_\mu \left( \lambda_{\alpha_1 \alpha_2 \cdots \alpha_n} \right)$ with $n \geq 2$, we can define $\lambda^{E}$ , $\lambda_\mu^{E}$ from 
\begin{align}\label{def}
\lambda = \lambda^{E} + g  \left( \lambda_{\alpha_1 \alpha_2 \cdots \alpha_n} \right) \, , \quad \lambda_\mu = \lambda_\mu^{E} + g_\mu  \left( \lambda_{\alpha_1 \alpha_2 \cdots \alpha_n} \right) \, .
\end{align}
This can be also considered as a change of independent variables from $\lambda$,  $\lambda_\mu$,  $\lambda_{\alpha_1 \alpha_2 \cdots \alpha_n} $ to $\lambda^{E} $,  $\lambda_\mu^{E}$,  $\lambda_{\alpha_1 \alpha_2 \cdots \alpha_n}$ for $n \geq 2$. We choose here  $g \left( \lambda_{\alpha_1 \alpha_2 \cdots \alpha_n} \right)$, $g_\mu \left( \lambda_{\alpha_1 \alpha_2 \cdots \alpha_n} \right)$ the solution of the conditions $g^E=0$, $g_\mu^E=0$ and of 
\begin{align}\label{2.0}
\begin{split}
& U_\alpha \int_{\Re^3} \int_{0}^{+ \infty} \, f_E  \left( e^{\frac{-1}{k_B} \left[m \, g + p^\mu g_\mu \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) + \Delta \chi \right]} - 1 \right) \,   p^{\alpha} \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} = 0 \, , \\
& U_\alpha \int_{\Re^3} \int_{0}^{+ \infty} \, f_E  \left( e^{\frac{-1}{k_B} \left[m \, g + p^\mu g_\mu \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) + \Delta \chi \right]} - 1 \right) \,   p^{\alpha} \,   p^{\beta} \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} = 0 \, , 
\end{split}
\end{align}
with $\Delta \chi= \sum_{n=2}^{N} \frac{1}{m^{n-1}} \, \lambda_{\alpha_1 \alpha_2 \cdots \alpha_n} p^{\alpha_1} p^{\alpha_2} \cdots p^{\alpha_n} \,\left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right)^n $. 
We will prove in the appendix  that this problem gives one and only one solution. After that, we propose for $Q$ the following expression 
\begin{align}\label{Q}
Q = - \, \frac{U_\alpha p^\alpha}{c^2 \tau} \, f_E \left( e^{\frac{-1}{k_B} \left[m \, g + p^\mu g_\mu \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} e^{\frac{-\, \Delta \chi}{k_B} } - 1 \right)\, .  
\end{align}
We note that, with this expression of $Q$ we have $Q^E=0$ and \eqref{dis}$_3$ gives $I=0$, $I^\alpha=0$ automatically (Thanks to \eqref{2.0}) so that the conservation laws of mass and of momentum-energy are satisfied up to whatever order; moreover, it follows that 
\begin{align}\label{I}
\begin{split}
& I^{\alpha_1 \cdots \alpha_n}=  \frac{U_\alpha}{c \, \tau}  \frac{\partial \, Q^\alpha}{\partial \,  \lambda_{\alpha_1 \alpha_2 \cdots \alpha_n} }  \, , \quad \mbox{with} \quad Q^\alpha = g \, \frac{\rho}{c} \, U^\alpha +  \, \frac{e}{c} \, \left( U^\mu g_\mu \right) U^\alpha + \\
& + k_B \int_{\Re^3} \int_{0}^{+ \infty} \, f_E  \left( e^{\frac{-1}{k_B} \left[m \, g + p^\mu g_\mu \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) + \Delta \chi \right]} + \, \frac{\Delta \chi}{k_B}  \right) \,   p^{\alpha} \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} \, .
\end{split}
\end{align}
(Note that $\rho$, $U^\alpha$ and the energy $e$ depend on $\lambda^{E} $, $\lambda^{E}_\mu$ and not on  $\lambda_{\alpha_1 \alpha_2 \cdots \alpha_n}$ for $n \geq 2$). To prove the previous result, we can calculate 
\begin{align*}
\begin{split}
& \frac{\partial \, Q^\alpha}{\partial \,  \lambda_{\alpha_1 \alpha_2 \cdots \alpha_n} } = \frac{\rho}{c} \, U^\alpha \frac{\partial \, g}{\partial \,  \lambda_{\alpha_1  \cdots \alpha_n} } + \frac{e}{c} \,  U^\mu \, U^\alpha \frac{\partial \, g_\mu}{\partial \,  \lambda_{\alpha_1  \cdots \alpha_n} } \, - \\
& \int_{\Re^3} \int_{0}^{+ \infty} f_E  \, e^{\frac{-1}{k_B} \left[m \, g + p^\mu g_\mu \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) + \Delta \chi \right]}   \left[ m \, \frac{\partial \, g}{\partial \,  \lambda_{\alpha_1 \alpha_2 \cdots \alpha_n} }  + p^\mu \frac{\partial \, g_\mu}{\partial \,  \lambda_{\alpha_1  \cdots \alpha_n} } \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) \right]  p^{\alpha}  \cdot \\
& \cdot \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} - \int_{\Re^3} \int_{0}^{+ \infty} f_E \left( e^{\frac{-1}{k_B} \left[m \, g + p^\mu g_\mu \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) + \Delta \chi \right]}    -  1 \right)   \frac{\partial \, \Delta \chi}{\partial \,  \lambda_{\alpha_1  \cdots \alpha_n} } 
 p^{\alpha} \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} .
\end{split}
\end{align*}
by contracting with $\frac{U_\alpha}{c \, \tau}  $ we obtain
\begin{align*}
\begin{split}
& \frac{U_\alpha}{c \, \tau} \, \frac{\partial \, Q^\alpha}{\partial \,  \lambda_{\alpha_1 \alpha_2 \cdots \alpha_n} } = \\
& = \frac{1}{c^2 \, \tau} \, \frac{\partial \, g}{\partial \,  \lambda_{\alpha_1  \cdots \alpha_n} } 
\left( \rho \, c^2  \, - \,U_\alpha \,  m \, c \int_{\Re^3} \int_{0}^{+ \infty} f_E  \, e^{\frac{-1}{k_B} \left[m \, g + p^\mu g_\mu \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) + \Delta \chi \right]}   p^{\alpha}  \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} \right) + \\
& + \frac{1}{c^2 \, \tau} \, \frac{\partial \, g_\mu}{\partial \,  \lambda_{\alpha_1  \cdots \alpha_n} }
\left( e  \,  U^\mu - U_\alpha c \int_{\Re^3} \int_{0}^{+ \infty} f_E  \, e^{\frac{-1}{k_B} \left[m \, g + p^\mu g_\mu \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) + \Delta \chi \right]}    p^\mu \right. \cdot \\
& \left. \cdot \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right)  p^{\alpha}  \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} \right)  
 +  c \int_{\Re^3} \int_{0}^{+ \infty} Q \,  \frac{\partial \, \Delta \chi}{\partial \,  \lambda_{\alpha_1  \cdots \alpha_n} } 
 \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} =I^{\alpha_1 \cdots \alpha_n} \, , 
\end{split}
\end{align*}
because the coefficients of $\frac{\partial \, g}{\partial \,  \lambda_{\alpha_1  \cdots \alpha_n} } $ and of $\frac{\partial \, g_\mu}{\partial \,  \lambda_{\alpha_1  \cdots \alpha_n} } $ are identically zero for \eqref{2.0}. In this way the proof of \eqref{I} is completed. Hence the balance equations take the elegant form \eqref{last}
where the gradient form with respect to $\lambda_{\beta_1 \cdots \beta_n} $ is present both in the left and in the right hand side. Obviously, there is in \eqref{last} the drawback that the left hand sides uses the variables $\lambda$, $\lambda_\alpha$, $\lambda_{\alpha_1  \cdots \alpha_n}$ while the right hand sides uses the variables $\lambda^E$, $\lambda_\alpha^E$, $\lambda_{\alpha_1  \cdots \alpha_n}$; however we can express $Q^\alpha$ in terms of the old variables by substituting in it the inverse transformation of \eqref{def}. In this way $Q^\alpha$ is the composite function of $\tilde{Q}^\alpha \left( \lambda \, , \, \lambda_\alpha \, , \, \lambda_{\alpha_1  \cdots \alpha_n} \right)$ and of  $\lambda =\lambda^E + g$, $\lambda_\alpha= \lambda_\alpha^E + g_\alpha$, $\lambda_{\alpha_1  \cdots \alpha_n}$ and \eqref{last}$_3$ becomes \eqref{last2}
where the gradient form appears also in the right hand side even if through 3 terms. \\
There remains to prove in this section the \\
{\bf H-Theorem}: " The entropy production $\sigma=\sum_{n=2}^{N} I^{\beta_1 \cdots \beta_n}  \lambda_{\beta_1 \cdots \beta_n}$ is non negative and is zero only at equilibrium". \\
To prove it we calculate
\begin{align}\label{3n}
\begin{split}
& \sigma = - \,  \frac{U_{\alpha} }{c \,   \tau} \,  \int_{\Re^3} \int_{0}^{+ \infty} \, f_E \left( e^{\frac{-1}{k_B} \left[m \, g + p^\mu g_\mu \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} e^{\frac{-\, \Delta \chi}{k_B} } - 1 \right) \,  \Delta \, \chi  \, p^\alpha \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} = \\
& = \frac{k_B U_{\alpha} }{c \,   \tau} \,  \int_{\Re^3} \int_{0}^{+ \infty} \, f_E \left( - \, e^{\frac{-1}{k_B} \left[m \, g + p^\mu g_\mu \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} e^{\frac{-\, \Delta \chi}{k_B} } + 1 \right)  \, p^\alpha \, +
\\
& +   \frac{1}{k_B} \, \left[ \Delta \, \chi +\underline{ m \, g \, + \, g_\mu \, p^\mu  \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right)} \quad \right]
\varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P}  \, . 
\end{split}
\end{align}
Here the underlined terms give a zero contribute thanks to \eqref{2.0}; they have been included for convenience of calculations. In fact, in this way the function to be integrated has the form 
\begin{align*}
F(x)= \left( 1 \, - \, e^{-x} \right) x \quad \mbox{with} \quad x= \frac{1}{k_B} \, \left[ \Delta \, \chi + \, m \, g \, + \, g_\mu \, p^\mu  \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right)  \right]
\end{align*}
and we have 
\begin{align*}
F'(x) = 1 \, +  e^{- \, x} (x-1) \, , \, F''(x) =  e^{- \, x} (-x+2) \, , \,  F'(0) = 0 \, , \,
  \lim_{x \rightarrow+ \infty} F'(x) = 1 \, .
\end{align*}
From these calculations we note that for $x<2$ the function $F'(x)$ is increasing so that it can be have only a root which is $x=0$; for $x \geq 2$ it is a decreasing function and goes from $1 + e^{-2}$ to $1$ so that $F'(x) > 0$ for $x \geq 2$. It follows that the function $F(x)$ is decreasing for $x<0$ and increasing for $x>0$; therefore it has a minimum value in $x=0$. Since $g(0)=0$, it follows that $g(0)>0$ $\forall$ $x \neq 0$, as we wanted to prove. \\
In order to compare the present results with those of \cite{Giap}, we conclude this subsection by seeing its implication on the linear expressions. At first order with respect to equilibrium, eqs. \eqref{2.0} become
\begin{align}\label{ecK3}
\begin{split}
& \rho \, c^2 g^{(1)} + e \, g^{(1)}_\mu \, U^\mu + U_\alpha \sum_{n=2}^{N} A_E^{\alpha \alpha_1  \cdots \alpha_n } \lambda_{\alpha_1  \cdots \alpha_n} = 0 \, , \\
& e \, U^\beta g^{(1)} + U_\alpha A_E^{\mu \alpha \beta} \, g^{(1)}_\mu  + U_\alpha \sum_{n=2}^{N} A_E^{\alpha \beta \alpha_1  \cdots \alpha_n } \lambda_{\alpha_1  \cdots \alpha_n} = 0 \, ,
\end{split}
\end{align}
where $g^{(1)}$ and  $g^{(1)}_\mu$ denote the first order terms of $g$ and  $g_\mu$ respectively. But in the Eckart approach we have also the following conditions \eqref{Eck}$_{1,2}$:
\begin{align}\label{Eck}
\begin{split}
& V^\alpha - \, V^\alpha_E = 0 \, , \quad  U_\alpha U_\beta \left( T^{\alpha \beta} - \, T^{\alpha \beta}_E \right) = 0 \, , \\  
& h_\alpha^\delta U_\beta \left( T^{\alpha \beta} - \, T^{\alpha \beta}_E \right) = - \, q^\delta \, \rightarrow \,  U_\beta \left( T^{\alpha \beta} - \, T^{\alpha \beta}_E \right) = q^\alpha \, , 
\end{split}
\end{align}
where \eqref{Eck}$_{3}$ isn't a condition but only the definition of the heat flux $q^\alpha$; moreover, \eqref{Eck}$_{2,3}$ are consequences of \eqref{Eck}$_{4}$ and of $U_\alpha \,  q^\alpha=0$. The conditions \eqref{Eck}$_{1,4}$ at first order with respect to equilibrium become 
\begin{align}\label{Eck2}
\begin{split}
& \rho \, U^\alpha \, \left( \lambda -  \lambda^E \right)^{(1)} + \left( \frac{e}{c^2} \, U^\alpha U^\mu + p h^{\alpha \mu} \right)\left( \lambda_\mu -  \lambda^E_\mu \right)^{(1)}  + \sum_{n=2}^{N} A_E^{\alpha \alpha_1  \cdots \alpha_n } \lambda_{\alpha_1  \cdots \alpha_n} = 0 \, , \\
& e \, U^\alpha \, \left( \lambda -  \lambda^E \right)^{(1)} + U_\beta \, A_E^{\alpha \beta \mu} \left( \lambda_\mu -  \lambda^E_\mu \right)^{(1)}  + U_\beta \sum_{n=2}^{N} A_E^{\alpha \beta \alpha_1  \cdots \alpha_n } \lambda_{\alpha_1  \cdots \alpha_n} = - \, \frac{k_B}{m} \, q^\alpha   \, .
\end{split}
\end{align}
We can desume from these equations $\sum_{n=2}^{N} A_E^{\alpha \alpha_1  \cdots \alpha_n } \lambda_{\alpha_1  \cdots \alpha_n}$ and $U_\beta \sum_{n=2}^{N} A_E^{\alpha \beta \alpha_1  \cdots \alpha_n } \lambda_{\alpha_1  \cdots \alpha_n} $; by substituting them in 
\eqref{ecK3} this equation becomes 
\begin{align}\label{ecK4}
\begin{split}
& \rho \, c^2 \left[ g^{(1)} - \left( \lambda -  \lambda^E \right)^{(1)} \right] + e \, \left[ g^{(1)}_\mu - \left( \lambda_\mu -  \lambda^E_\mu \right)^{(1)} \right] \, U^\mu= 0 \, , \\
& e \, U^\beta \left[ g^{(1)} - \left( \lambda -  \lambda^E \right)^{(1)} \right]  + U_\alpha A_E^{\mu \alpha \beta} \, \left[ g^{(1)}_\mu - \left( \lambda_\mu -  \lambda^E_\mu \right)^{(1)} \right]     - \, \frac{k_B}{m} \, q^\beta = 0 \, .
\end{split}
\end{align}
Now eq. \eqref{ecK4}$_{1}$, \eqref{ecK4}$_{2}$ contracted with $U_\beta$ and 
\eqref{ecK4}$_{2}$ contracted with $h^\delta_\beta$ give 
\begin{align}\label{ecK5}
\begin{split}
& g^{(1)}=  \left( \lambda -  \lambda^E \right)^{(1)} \, , \quad  g^{(1)}_\mu \, U^\mu = \left( \lambda_\mu -  \lambda^E_\mu \right)^{(1)} \, U^\mu \, , \\ 
& h^{\mu \delta}  \, \left[ g^{(1)}_\mu - \left( \lambda_\mu -  \lambda^E_\mu \right)^{(1)} \right]     = \frac{3}{\rho c^4 \theta_{1,2}}\, \frac{k_B}{m} \, q^\delta \, , 
\end{split}
\end{align}
where eq. (14) of \cite{ACPR} has been used. Now \eqref{Q} can be written as 
\begin{align*}
Q = - \, \frac{U_\alpha p^\alpha}{c^2 \tau} \, f_E \left[ \frac{f}{f_E} \, - \, 1 \, + \frac{f}{f_E} \left( e^{\frac{-1}{k_B} \left[m \, \left( g - \lambda + \lambda^E \right) + p^\mu \left( g_\mu - \lambda_\mu + \lambda_\mu^E \right) \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} - 1 \right) \right] \, .  
\end{align*}
Here we can linearize the term $\frac{f}{f_E} \left( e^{\frac{-1}{k_B} \left[m \, \left( g - \lambda + \lambda^E \right) + p^\mu \left( g_\mu - \lambda_\mu + \lambda_\mu^E \right) \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} - 1 \right)$, i.e., substitute it with 
\begin{align*}
\begin{split}
& \left( \frac{f}{f_E} \right)^{(0)} \left( e^{\frac{-1}{k_B} \left[m \, \left( g - \lambda + \lambda^E \right) + p^\mu \left( g_\mu - \lambda_\mu + \lambda_\mu^E \right) \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} - 1 \right)^{(1)} + \\
& + \left( \frac{f}{f_E} \right)^{(1)} \left( e^{\frac{-1}{k_B} \left[m \, \left( g - \lambda + \lambda^E \right) + p^\mu \left( g_\mu - \lambda_\mu + \lambda_\mu^E \right) \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} - 1 \right)^{(0)} = \\
& = \frac{-1}{k_B} \left[m \, \left( g - \lambda + \lambda^E \right)^{(1)} + p^\mu \left( g_\mu - \lambda_\mu + \lambda_\mu^E \right)^{(1)} \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) \right] \, . 
\end{split}
\end{align*}
Moreover, we use the result \eqref{ecK5} and find 
\begin{align*}
Q = - \, \frac{U_\alpha p^\alpha}{c^2 \tau} \, f_E \left[ \frac{f}{f_E} \, - \, 1 \, +     p^\mu \, q_\mu \, \frac{3}{m \, \rho c^4 \theta_{1,2}}  
 \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) \right] \, ,
\end{align*}
which is the expression found in eq. (43) of \cite{ACPR}.


\subsection{The Landau-Liftschiz approach}
We note that all  the considerations in the previous subsection, up to its equation \eqref{last}, can be repeated also in the present approach and eqs. \eqref{2.0} simply means that 
\begin{align}\label{def2}
U_\alpha \left( V^\alpha - \, V^\alpha_{eq} \right) = 0 \, , \quad U_\alpha \left( T^{\alpha \beta} - \, T^{\alpha \beta}_{eq} \right) = 0 \, ,
\end{align}
because by using \eqref{def}, we have 
\begin{align*}
 f_E  \left( e^{\frac{-1}{k_B} \left[m \, g + p^\mu g_\mu \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) + \Delta \chi \right]}  - 1 \right) = f \, - \, f_E \, .
\end{align*}
Moreover, the expression \eqref{Q} of the production term can be written simply as 
\begin{align*}
Q = \frac{- \, U_{\alpha} \, p^\alpha }{c^2  \tau} \, \left( f -f_E \right) \, .
\end{align*}
From this viewpoint it seems that the definition of deviations from equilibrium, which is present in the Landau-Liftschiz approach, was purposely made to automatically achieve zero mass production and zero momentum-energy production:
\begin{align*}
I = - \, \frac{U_{L \alpha} }{c^2  \tau} \, \left( V^\alpha -V_E^\alpha \right) = 0 \, , \,  I^{\beta} = - \, \frac{U_{L \alpha} }{c^2  \tau} \, \left( T^{\alpha \beta} -T_E^{\alpha \beta} \right) = 0 \, , 
\end{align*}
So also in this approach we obtain the equation \eqref{last}, 
where the gradient form with respect to $\lambda_{\beta_1 \cdots \beta_n} $ is present both in the left and in the right hand side. Moreover, we don't have here the drawback which was present  in \eqref{last} which forced us to transform it in the less elegant form \eqref{last2}. It is true that also in the present approach the left hand sides of \eqref{last} uses the variables $\lambda$, $\lambda_\alpha$, $\lambda_{\alpha_1  \cdots \alpha_n}$ while the right hand sides uses the variables $\lambda^E$, $\lambda_\alpha^E$, $\lambda_{\alpha_1  \cdots \alpha_n}$; but, instead of  expressing $Q^\alpha$ in terms of the old variables, we can express the left hand sides (i.e., $h'^\alpha$) in terms of the new variables. Since an invertible change of independent variables maintains the hyperbolicity requirement, as long as this change of independent variables is not done in an approximated way, this requirement is here preserved. Also the proof of the 
{\bf H-Theorem} which is present in the previous subsection, between eqs. \eqref{last2} and \eqref{ecK3} still holds also in the present case. \\
We see also that the expression of $Q$ proposed in eq. (7) of \cite{CrisPen2} (with $a_1=0$, $a_2= 1/ \tau $) for the Eckart approach is the same of the present one \eqref{Q}, except to identify the functions $\psi$ and $\psi_\mu$ of \cite{CrisPen2} respectively with the  funtions $g$ and $g_\mu$ which are here present. However, they were introduced in \cite{CrisPen2}  ad hoc, as a mathematical tool; instead here we have seen that they come from trying to adapt the Eckart approach to the Landau-Liftschiz approach. 

\section{Determination of $q^\alpha_L$ and $P^{\alpha \beta}$ in terms of the variables in the Eckart approach}

Both the approaches are the same at equilibrium so that they give the same values for $\rho$, $U^\alpha$, $e$ and $p$; in particular they give 
\begin{align*}
\begin{split}
& \rho \, U^\alpha = \frac{\partial \, h'^\alpha_E}{\partial \, \lambda^E} \, , \, 
T^{\alpha \beta}_E = \frac{e}{c^2} \, U^\alpha U^\beta \, + \, p \, h^{\alpha \beta} = \frac{\partial \, h'^\alpha_E}{\partial \, \lambda_\beta^E} \quad \mbox{with} \\
& h'^\alpha_E = - \,  k_B c \int_{\Re^3} \int_{0}^{+ \infty} \, f_E \,  p^\alpha \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} \, , \, f_E = e^{- 1 - \, \frac{1}{k_B} \, \left[ m \, \lambda^E + \lambda^E_\mu p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} \, .
\end{split}
\end{align*}
Outside equilibrium, with the Eckart approach we have $f=f^{Ek}$ with 
\begin{align*}
f^{Ek} = f_E \,  e^{- \, \frac{1}{k_B} \, \left[ m \, \left(\lambda \, - \, \lambda^E \right) + \left(\lambda_\mu \, - \, \lambda^E_\mu\right) p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \, + \, \Delta \chi \right]} \, , \\
\Delta \chi= \sum_{n=2}^{N} \frac{1}{m^{n-1}} \lambda_{\alpha_1  \cdots \alpha_n} p^{\alpha_1}  \cdots p^{\alpha_n} \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right)^n \, ,  
\end{align*}
while  with the Landau-Liftschiz approach we have $f=f^{L}$ with 
\begin{align}\label{A1}
\begin{split}
& f^{L} = f_E \,  e^{- \, \frac{1}{k_B} \, \left[ m \, g + g_\mu p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \, + \, \Delta \chi \right]} \,  \, \rightarrow \\ 
& \frac{f^{L}}{f^{Ek}} =   e^{- \, \frac{1}{k_B} \, \left[ m \, \left( g \, - \, \lambda \, + \, \lambda^E \right) + \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right) p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} \, , \\
& \frac{f^{L} - f_E}{f_E} = \frac{f^{Ek}}{f_E} \, \left( e^{- \, \frac{1}{k_B} \, \left[ m \, \left( g \, - \, \lambda \, + \, \lambda^E \right) + \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right) p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} - 1  \right) \, + \, \frac{f^{Ek}}{f_E} \, - \, 1 \, .  
\end{split}
\end{align}
Now $q^\alpha_L$ and $P^{\alpha \beta}$ are defined by $q^\alpha_L= V^\alpha_L - \rho \, U^\alpha_L$ and $P^{\alpha \beta}=  T^{\alpha \beta}_L \, - \, T^{\alpha \beta}_E$ so that we have the system 
\begin{align}\label{A2}
\begin{split}
& q^\alpha_L = m \, c \int_{\Re^3} \int_{0}^{+ \infty} \, f_E \, \frac{f^{L} - f_E}{f_E} \,  p^\alpha \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} \, , \\
& P^{\alpha \beta} =  c \int_{\Re^3} \int_{0}^{+ \infty} \, f_E \, \frac{f^{L} - f_E}{f_E} \,  p^\alpha \, p^\beta \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} \, , \\
& U_\alpha \, q^\alpha_L=0 \, , \quad U_\alpha \, P^{\alpha \beta} =0 \, , 
\end{split}
\end{align}
for the determination of $g \, - \, \lambda \, + \, \lambda^E$, $g_\mu \, - \, \lambda_\mu
\, + \, \lambda_\mu^E$, $q^\alpha_L$,  $P^{\alpha \beta}$. To this end, from \eqref{A1}$_3$ we have that $\frac{f^{L} - f_E}{f_E}$ at the order zero is zero, while at the orders 1 and 2 is respectively given by
\begin{align}\label{A3}
\begin{split}
& \left( \frac{f^{L} - f_E}{f_E} \right)^{(1)} = \left( \frac{f^{Ek}}{f_E} \right)^{(0)} \, \left( e^{- \, \frac{1}{k_B} \, \left[ m \, \left( g \, - \, \lambda \, + \, \lambda^E \right) + \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right) p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} - 1  \right)^{(1)} \, + \\
& + \left( \frac{f^{Ek}}{f_E} \right)^{(1)} \, \underline{\left( e^{- \, \frac{1}{k_B} \, \left[ m \, \left( g \, - \, \lambda \, + \, \lambda^E \right) + \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right) p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} - 1  \right)^{(0)} } \, + 
\, \left( \frac{f^{Ek}}{f_E} \, - \, 1 \right)^{(1)} = \\
& =  - \, \frac{1}{k_B} \, \left[ m \, \left( g \, - \, \lambda \, + \, \lambda^E \right)^{(1)} + \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(1)} p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right] \, + 
\, \left( \frac{f^{Ek}}{f_E} \, - \, 1 \right)^{(1)} \, , 
\end{split}
\end{align} 
\begin{align*}
\begin{split}
& \left( \frac{f^{L} - f_E}{f_E} \right)^{(2)} = \left( \frac{f^{Ek}}{f_E} \right)^{(0)} \, \left( e^{- \, \frac{1}{k_B} \, \left[ m \, \left( g \, - \, \lambda \, + \, \lambda^E \right) + \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right) p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} - 1  \right)^{(2)} \, + \\
& + \, \left( \frac{f^{Ek}}{f_E} \right)^{(1)} \, \left( e^{- \, \frac{1}{k_B} \, \left[ m \, \left( g \, - \, \lambda \, + \, \lambda^E \right) + \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right) p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} - 1  \right)^{(1)} \, + \\
& + \left( \frac{f^{Ek}}{f_E} \right)^{(2)} \, \underline{\left( e^{- \, \frac{1}{k_B} \, \left[ m \, \left( g \, - \, \lambda \, + \, \lambda^E \right) + \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right) p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} - 1  \right)^{(0)} }\, + 
\, \left( \frac{f^{Ek}}{f_E} \, - \, 1 \right)^{(2)} \, .
\end{split}
\end{align*}
where the underlined terms are zero. 

\subsection{The expressions of  $q^\alpha_L$ and $P^{\alpha \beta}$ at first order with respect to equilibrium}

By using the above passages we see that the system \eqref{A2} at first order becomes
\begin{align}\label{A4}
\begin{split}
& \left( q^\alpha_L \right)^{(1)} = m \, c \int_{\Re^3} \int_{0}^{+ \infty} \, f_E \, \left( \frac{f^{L} - f_E}{f_E} \right)^{(1)} \,  p^\alpha \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} = \\
& = - \, \frac{m}{k_B} \, \left( g \, - \, \lambda \, + \, \lambda^E \right)^{(1)} \rho \, U^\alpha \, - \, \frac{m}{k_B} \,  \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(1)} T^{\alpha \mu}_E \, + \, \underline{ \left( V^\alpha_{Ek} - V^\alpha_E \right)^{(1)} } \, , \\
& \left( P^{\alpha \beta} \right)^{(1)} =  c \int_{\Re^3} \int_{0}^{+ \infty} \, f_E \, \left( \frac{f^{L} - f_E}{f_E} \right)^{(1)} \,  p^\alpha \, p^\beta \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} = \\
& = - \, \frac{m}{k_B} \, \left( g \, - \, \lambda \, + \, \lambda^E \right)^{(1)} T^{\alpha \beta}_E \, - \, \frac{m}{k_B} \,  \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(1)} A^{\alpha \beta \mu}_E \, + \,  \left( T^{\alpha \beta}_{Ek} - T^{\alpha \beta}_E \right)^{(1)}  \, , \\
& U_\alpha \, \left( q^\alpha_L \right)^{(1)}=0 \, , \quad U_\alpha \, \left( P^{\alpha \beta} \right)^{(1)} =0 \, . 
\end{split}
\end{align}
Here too the underlined term is zero, while 
\begin{align*}
\begin{split}
& T^{\alpha \beta}_E = \frac{e}{c^2} \, U^\alpha U^\beta + p \, h^{\alpha \beta} \, , \, A^{\alpha \beta \mu}_E =  \rho \, \theta_{0,2} U^\alpha U^\beta U^\mu \, + \, \rho \, c^2 \theta_{1,2} U^{( \alpha} h^{\beta \mu )} \, , \\
& \left( T^{\alpha \beta}_{Ek} - T^{\alpha \beta}_E \right)^{(1)} = \pi \, h^{\alpha \beta} + \frac{2}{c^2} U^{( \alpha} q^{\beta )} + t^{< \alpha \beta>} \, , 
 \end{split}
\end{align*}
as it can be seen from \cite{ACPR}, in particular from its eq. (14). By contracting eq. \eqref{A4}$_1$ with $U_\alpha$, eq. \eqref{A4}$_2$ with $U_\alpha U_\beta$ and by taking into account \eqref{A4}$_{2,3}$ we obtain
\begin{align}\label{A5}
\left( g \, - \, \lambda \, + \, \lambda^E \right)^{(1)} = 0 \, , \, U^\mu \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(1)} = 0 \, . 
\end{align}
By contracting eq.  \eqref{A4}$_2$ with $U_\alpha h_\beta^\delta$ and by taking into account \eqref{A4}$_{3}$ we obtain
\begin{align}\label{A6}
\left(g_\delta  \, - \, \lambda_\delta  \, + \, \lambda^E_\delta \right)^{(1)} = - \,  \frac{3}{\rho \, c^4 \theta_{1,2}} \, \frac{k_B}{m} \, q_\delta \, . 
\end{align}
There remains to contract eq. \eqref{A4}$_1$ with $h_\alpha^\delta$ and eq. \eqref{A4}$_2$ with $h_\alpha^\delta h_\beta^\psi$; by using eqs. \eqref{A5}, \eqref{A6} the result is 
\begin{align}\label{A7}
\left( q^\delta_L \right)^{(1)} = - \,  \frac{3 \, p}{\rho \, c^4 \theta_{1,2}} \,  q^\delta \, , \quad \left( P^{\delta \psi} \right)^{(1)} = \pi \, h^{\delta \psi} \, + \, t^{< \delta \psi >} \, . 
 \end{align}
Regarding the first term of this equation we note that it can be rewritten by using the recurrence relations (16) and eqs. (12)$_{1,2}$ of \cite{ACPR} as 
\begin{align*}
\left( q^\delta_L \right)^{(1)} = - \,  \frac{\rho}{e+p} \,  q^\delta \, , 
\end{align*}
which is the same value found in (7)$_1$ of \cite{CrisPen} which concerned an approach with a linear deviation from equilibrium. Since we will see in the next subsection that $\left( q^\delta_L \right)^{(2)} \neq 0$, the expression  (7)$_1$ of \cite{CrisPen} cannot be assumed to hold up to whatever order with respect to equilibrium. 
 
 \subsection{The expressions of  $q^\alpha_L$ and $P^{\alpha \beta}$ at second order with respect to equilibrium}
We have to consider eq. \eqref{A3}$_2$; to this end, we need \\  
$\left( e^{- \, \frac{1}{k_B} \, \left[ m \, \left( g \, - \, \lambda \, + \, \lambda^E \right) + \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right) p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} - 1  \right)^{(2)}$ and \\
$\left( e^{- \, \frac{1}{k_B} \, \left[ m \, \left( g \, - \, \lambda \, + \, \lambda^E \right) + \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right) p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} - 1  \right)^{(1)}$. \\
Since the Taylor's expansion of the function $e^x$ around $x=0$ and up to second order is $1+x+ \frac{x^2}{2}$, we have 
\begin{align*}
\begin{split}
& \left( e^{- \, \frac{1}{k_B} \, \left[ m \, \left( g \, - \, \lambda \, + \, \lambda^E \right) + \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right) p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} - 1  \right)^{(2)} = - \, \frac{m}{k_B} \, \left( g \, - \, \lambda \, + \, \lambda^E \right)^{(2)} \, - \\
&  \, \frac{1}{k_B} \, \left[ \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(2)}  p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right] + \frac{1}{2} \left(  \frac{ m \, \left( g \, - \, \lambda \, + \, \lambda^E \right)^{(1)}}{k_B} \right)^2  + \\
& +  \frac{1}{\left( k_B \right)^2} \, \left[  m \, \left( g \, - \, \lambda \, + \, \lambda^E \right)^{(1)} \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(1)} p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right] \\
& + \frac{1}{2 \left( k_B \right)^2} \, \left[ \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(1)} p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right] \left[ \left(g_\nu \, - \, \lambda_\nu \, + \, \lambda^E_\nu\right)^{(1)} p^\nu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right] \, , 
\end{split}
\end{align*}
\begin{align*}
\begin{split}
& \left( e^{- \, \frac{1}{k_B} \, \left[ m \, \left( g \, - \, \lambda \, + \, \lambda^E \right) + \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right) p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} - 1  \right)^{(1)} = - \, \frac{m}{k_B} \, \left( g \, - \, \lambda \, + \, \lambda^E \right)^{(1)} - \\
& \frac{1}{k_B} \, \left[ \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(1)} p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right] \, .
\end{split}
\end{align*}
By using \eqref{A5} and \eqref{A6}, we obtain
\begin{align*}
\begin{split}
& \left( e^{- \, \frac{1}{k_B} \, \left[ m \, \left( g \, - \, \lambda \, + \, \lambda^E \right) + \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right) p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} - 1  \right)^{(2)} = - \, \frac{m}{k_B} \, \left( g \, - \, \lambda \, + \, \lambda^E \right)^{(2)} \, - \\
&  \, \frac{1}{k_B} \, \left[ \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(2)}  p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right]   + \, \frac{1}{2} \, \left(\frac{3}{\rho \, c^4 \theta_{1,2}} \, \frac{1}{m} \right)^2 \, \left( q_\mu  p^\mu \right)^2 \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right)^2  \, , 
\end{split}
\end{align*}
\begin{align*}
\begin{split}
& \left( e^{- \, \frac{1}{k_B} \, \left[ m \, \left( g \, - \, \lambda \, + \, \lambda^E \right) + \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right) p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right]} - 1  \right)^{(1)} =  \frac{3}{\rho \, c^4 \theta_{1,2}} \, \frac{1}{m} \, q_\mu p^\mu \,  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \, . 
\end{split}
\end{align*}
By using these expressions in  \eqref{A3}$_2$ it becomes 
\begin{align*}
\begin{split}
& \left( \frac{f^{L} - f_E}{f_E} \right)^{(2)} =  - \, \frac{m}{k_B} \, \left( g \, - \, \lambda \, + \, \lambda^E \right)^{(2)} \, - \\
&  \, \frac{1}{k_B} \, \left[ \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(2)}  p^\mu  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \right]   + \, \frac{1}{2} \, \left(\frac{3}{\rho \, c^4 \theta_{1,2}} \, \frac{1}{m} \right)^2 \, \left( q_\mu  p^\mu \right)^2 \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right)^2 \, + \\
& + \, \left( \frac{f^{Ek}}{f_E} \right)^{(1)} \, \frac{3}{\rho \, c^4 \theta_{1,2}} \, \frac{1}{m} \, q_\mu p^\mu \,  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \,  + \, \left( \frac{f^{Ek}}{f_E} \, - \, 1 \right)^{(2)} \, .
\end{split}
\end{align*}
We can now substitute this result in \eqref{A2} so that the homogeneous second order part of this system is
\begin{align}\label{A8}
\begin{split}
& \left( q^\alpha_L \right)^{(2)} = m \, c \int_{\Re^3} \int_{0}^{+ \infty} \, f_E \, \left( \frac{f^{L} - f_E}{f_E} \right)^{(2)} \,  p^\alpha \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} = \\
& = - \, \frac{m}{k_B} \, \left( g \, - \, \lambda \, + \, \lambda^E \right)^{(2)} \, \rho \, U^\alpha - \, \frac{m}{k_B} \, \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(2)}  T^{\alpha \mu}_E   + \\
&  + \, \frac{1}{2} \, \left(\frac{3}{\rho \, c^4 \theta_{1,2}} \right)^2 \,   A^{\alpha \mu \nu}_E \, q_\mu \, q_\nu \,  +  \, \frac{3}{\rho \, c^4 \theta_{1,2}}  \, q_\mu  \, \left( T^{\alpha \mu}_{Ek} \right)^{(1)} + \underline{\left( V^\alpha_{Ek} \right)^{(2)}} = \\
& = - \, \frac{m}{k_B} \, \left( g \, - \, \lambda \, + \, \lambda^E \right)^{(2)} \, \rho \, U^\alpha - \, \frac{m}{k_B} \, \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(2)}  T^{\alpha \mu}_E   + \\
& +  \frac{3}{\rho \, c^4 \theta_{1,2}}  \, q_\mu  \, \left( \pi \, h^{\alpha \mu} \, + \, \frac{1}{2 \, c^2} \, U^{\alpha} q^{\mu } + \, t^{< \alpha \mu >} \right)
\, , 
\end{split}
\end{align}
\begin{align*}
\begin{split}
& \left( P^{\alpha \beta} \right)^{(2)} =  c \int_{\Re^3} \int_{0}^{+ \infty} \, f_E \, \frac{f^{L} - f_E}{f_E} \,  p^\alpha \, p^\beta \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) p^\mu \,  \left( 1 + \, \frac{\mathcal{I}}{m \, c^2} \right) \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} = \\
& = - \, \frac{m}{k_B} \, \left( g \, - \, \lambda \, + \, \lambda^E \right)^{(2)} \, T^{\alpha \beta}_E - \, \frac{m}{k_B} \, \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(2)}  A^{\alpha \beta \mu}_E   + \\
&  + \, \frac{1}{2} \, \left(\frac{3}{\rho \, c^4 \theta_{1,2}} \right)^2 \,   A^{\alpha \beta \mu \nu}_E \, q_\mu \, q_\nu \,  +  \, \frac{3}{\rho \, c^4 \theta_{1,2}}  \, q_\mu  \, \left( A^{\alpha \beta \mu}_{Ek} \right)^{(1)} + \underline{\left( T^{\alpha \beta}_{Ek} \right)^{(2)}} \, , \\
& U_\alpha \, \left( q^\alpha_L \right)^{(2)}=0 \, , \quad U_\alpha \, \left( P^{\alpha \beta} \right)^{(2)} =0 \, , 
\end{split}
\end{align*}
where the underlined terms are zero, while 
\begin{align*}
\begin{split}
&  A^{\alpha \beta \mu \nu}_E =  \rho \, \theta_{0,3} \, U^\alpha U^\beta U^\mu U^\nu \, + \, \rho \, c^2 \theta_{1,3} \, U^{( \alpha} U^\beta h^{\mu \nu)} \, + \, \rho \, c^4 \theta_{2,3} \, h^{( \alpha \beta} h^{\mu \nu)} \, , \\
& \left( A^{\alpha \beta \mu}_{Ek} \right)^{(1)} =  \frac{1}{4c^4} \, \Delta  U^{\alpha } U^{\beta} U^\mu - \left(   \frac{3}{4c^2}\, \frac{N^\Delta}{D_4} \, \Delta \, + 3 \, \frac{N^\Pi}{D_4} \, \Pi \right)
\,  h^{(\alpha\beta} U^{\mu)} \\
& + \frac{3}{c^2} \frac{N_3}{D_3}  \, q^{(\alpha} U^\beta U^{\mu)} + \frac{3}{5} \frac{N_{31}}{D_3} h^{(\alpha
	\beta} q^{\mu)} + 3  C_5 t^{(<\alpha \beta>_3} U^{\mu)}
\, ,
\end{split}
\end{align*}
as it can be seen from eqs. (14) and (35) of \cite{ACPR} ($\Delta$ is the 15$^{\mbox{th}}$ variable). By contracting eq. \eqref{A8}$_1$ with $\frac{U_\alpha}{c^2}$, eq. \eqref{A8}$_2$ with $\frac{U_\alpha U_\beta }{c^4}$ and by taking into account \eqref{A8}$_{3,4}$ we obtain
\begin{align*}
\begin{split}
& 0 = - \, \frac{m}{k_B} \, \rho \, \left( g \, - \, \lambda \, + \, \lambda^E \right)^{(2)}  - \, \frac{m}{k_B} \, \frac{e}{c^2} \, U^\mu \, \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(2)}     + \,
  \frac{3}{2 \, \rho \, c^6 \theta_{1,2}}  \, q_\mu  \,  q^{\mu }  \, , \\
& 0 =  - \, \frac{m}{k_B} \, \, \frac{e}{c^2} \, \left( g \, - \, \lambda \, + \, \lambda^E \right)^{(2)} \, - \, \frac{m}{k_B} \, \rho \, \theta_{0,2} \, U^\mu  \, \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(2)}     + \\
&  + \left( - \, \frac{3}{4} \, \frac{\theta_{1,3}}{\rho \, c^6 \, \left( \theta_{1,2} \right)^2}   \,  +  \,  \frac{3}{\rho \, c^6 \theta_{1,2}}   \, \frac{N_3}{D_3} \right)  \, q^{\mu} q_{\mu}   \, . 
\end{split}
\end{align*}
This system fully determine $ \left( g \, - \, \lambda \, + \, \lambda^E \right)^{(2)}$ and $U^\mu  \, \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(2)}$. \\
By contracting eq.  \eqref{A8}$_2$ with $U_\alpha h_\beta^\delta$ and by taking into account \eqref{A8}$_{3,4}$ we obtain
\begin{align}\label{A9}
\begin{split}
&    h^{\mu \delta} \, \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(2)}   =   - \,  \frac{k_B}{m}   \,  \frac{9}{\rho^2 \, c^8 \,  \left( \theta_{1,2} \right)^2 }  \, q_\mu  \, \left( A^{\alpha \beta \mu}_{Ek} \right)^{(1)} U_\alpha h_\beta^\delta = \\
&  = \frac{k_B}{m}   \,  \frac{9}{\rho^2 \, c^8 \,  \left( \theta_{1,2} \right)^2 }    \, \left[  \left(   \frac{1}{4}\, \frac{N^\Delta}{D_4} \, \Delta \, +  c^2 \, \frac{N^\Pi}{D_4} \, \Pi \right)
\,  q^{\delta}  \, +  \,   C_5 \, c^2  \,  t^{<\delta \mu >_3} \, q_\mu \right]  \, . 
\end{split}
\end{align}
There remains to contract eq. \eqref{A8}$_1$ with $h_\alpha^\delta$ and eq. \eqref{A8}$_2$ with $h_\alpha^\delta h_\beta^\psi$;  the result is 

\begin{align}\label{A10}
\begin{split}
& \left( q^\delta_L \right)^{(2)} =  - \,  \frac{m \, p}{k_B} \, h^{\delta \mu} \, \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(2)}      \,
 + \, \frac{3}{\rho \, c^4 \theta_{1,2}}   \, \left( - \, \pi \, q^{\delta} \, + \,  t^{< \delta \mu >} \, q_\mu \right) \, , 
\end{split}
\end{align}
\begin{align*}
\begin{split}
& \left( P^{\delta \psi} \right)^{(2)}  =  
 - \, \frac{m }{k_B} \,  \left[ p \, \left( g \, - \, \lambda \, + \, \lambda^E \right)^{(2)}   \, + \, \frac{1}{3} \, \rho \, c^2 \theta_{1,2}\, U^{ \mu} \, \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(2)} \right] h^{\delta \psi} + \\
&  +  \, \frac{3}{\rho \, c^4 \, \theta_{1,2} }  \,  \left( \frac{1}{2} \,  \frac{\theta_{2,3}}{\theta_{1,2}} \, - \,  \frac{1}{5} \frac{N_{31}}{D_3} \right)\, \left( 2 \, q^\delta \, q^\psi \, -  \, q^\mu\, q_\mu \, h^{\delta \psi} \right)  \, . 
\end{split}
\end{align*}
Obviously, here the above found expressions of $\left( g \, - \, \lambda \, + \, \lambda^E \right)^{(2)} $, $U^{ \mu} \, \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(2)}$ and of $ h^{\delta \mu} \, \left(g_\mu \, - \, \lambda_\mu \, + \, \lambda^E_\mu\right)^{(2)} $ must be used. \\
\\
\\
{\bf \Large Conclusions:} We have obtained the gradient form of the balance equations not only for their left hand sides, but also for their right hand sides; this result has been obtained both with the Eckart approach and the Landau-Lifshitz one. 
In this way the symmetric form  and a 4-potential is also  obtained  for both sides. Moreover, this form allowed to prove the hyperbolicity of the equations and the H-theorem up to whatever order. We have also seen that  the Eckart approach and the Landau-Lifshitz one are equivalent if the Lagrange multipliers $\lambda_{\alpha_1  \cdots \alpha_n}$ are taken as independent variables. The difference comes only after a different definition of the deviations of $\lambda - \lambda^E$ and of $\lambda_\alpha - \lambda_\alpha^E$; in any case, there is an invertible transformation between their independent variables. 


\appendix

\section{Appendix: Proof of the unicity of the solution of eqs. \eqref{2.0} and of $g^E=0$, $g_\mu^E=0$.}
The equations \eqref{2.0} calculated at equilbrium are identically satisfied; so they are equivalent to their derivatives with respect to $\lambda_{\alpha_1  \cdots \alpha_{n_1}}$ (with $n_1 \geq 2$) which now we express, by using a compact notation, as $\lambda_{A_{n_1}}$. These derivatives are
\begin{align}\label{2.1}
\begin{split}
& \frac{U_\alpha}{- \, k_B} \int_{\Re^3} \int_{0}^{+ \infty} \, f_E \,  F \,  \left[m \, \frac{\partial \, g}{\partial \, \lambda_{A_{n_1}}} + p^\mu \frac{\partial \, g_\mu }{\partial \, \lambda_{A_{n_1}}}  \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) + \frac{\partial \, \Delta \chi}{\partial \, \lambda_{A_{n_1}}}  \right] \,  p^{\alpha} \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} = 0 \, , \\
& \frac{U_\alpha}{- \, k_B} \int_{\Re^3} \int_{0}^{+ \infty} \, f_E \,  F \,  \left[m \, \frac{\partial \, g}{\partial \, \lambda_{A_{n_1}}} + p^\mu \frac{\partial \, g_\mu }{\partial \, \lambda_{A_{n_1}}}  \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) + \frac{\partial \, \Delta \chi}{\partial \, \lambda_{A_{n_1}}}  \right] \,    p^{\alpha} \,   p^{\beta} \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) \cdot \\
& \quad \quad \quad \cdot  \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} = 0 \, , \,
\quad  \mbox{with} \quad F=  e^{\frac{-1}{k_B} \left[m \, g + p^\mu g_\mu \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) + \Delta \chi \right]} \, , 
\end{split}
\end{align}
which, calculated at equilibrium, give
\begin{align*}
\begin{split}
& - \, \frac{m}{k_B} \left[ \left( \frac{\partial \, g}{\partial \, \lambda_{A_{n_1}}} \right)_E U_\alpha \, V^\alpha_E  + \,  \left( \frac{\partial \, g_\mu }{\partial \, \lambda_{A_{n_1}}} \right)_E \, U_\alpha \, T^{\alpha \mu}_E   \right]  =  \frac{m}{k_B} \,  U_\alpha A^{\alpha A_{n_1}}_E \, , \\
& - \, \frac{m}{k_B} \left[ \left( \frac{\partial \, g}{\partial \, \lambda_{A_{n_1}}} \right)_E U_\alpha \, T^{\alpha \beta}_E  + \,  \left( \frac{\partial \, g_\mu }{\partial \, \lambda_{A_{n_1}}} \right)_E \, U_\alpha \, A^{\alpha \beta \mu}_E  \right]  =  \frac{m}{k_B} \,  U_\alpha A^{\alpha \beta A_{n_1}}_E \, . 
\end{split}
\end{align*}
By contracting these equations with $\lambda_{A_{n_1}}$ and taking the sum for $n_1 =2, \, \cdots \, , \, N$ they become 
\begin{align}\label{2.2}
\begin{split}
&  g^{(1)} \,  U_\alpha \, V^\alpha_E  + \,  g_\mu^{(1)} \, U_\alpha \, T^{\alpha \mu}_E    = -  \, \sum_{n_1 =2}^{N} U_\alpha A^{\alpha A_{n_1}}_E \, \lambda_{A_{n_1}} \, , \\
& g^{(1)} \,  U_\alpha \, T^{\alpha \beta}_E  + \,  g_\mu^{(1)} \, U_\alpha \, A^{\alpha \beta \mu}_E   = - \, \sum_{n_1 =2}^{N}  U_\alpha A^{\alpha \beta A_{n_1}}_E \, \lambda_{A_{n_1}} \, . 
\end{split}
\end{align}
The first one of these equations, the second one  contracted with $U_\beta$ and the second one contracted with $h_\beta^\delta$ constitute a system by using also the expressions 
\begin{align*}
V^\alpha_E= \rho \, U^\alpha\, , \, T^{\alpha \mu}_E= \frac{e}{c^2} \, U^\alpha U^\mu \, + \, p \, h^{\alpha \mu} \, , \, A^{\alpha \beta \mu}_E = \rho \, \theta_{0,2} \, U^\alpha U^\beta U^\mu \, + \,  \rho \, c^2 \theta_{1,2} \, U^{( \alpha} h^{\beta \mu )} \, . 
\end{align*}
This system has the unique solution 
\begin{align*}
\begin{split}
& g^{(1)}  = \left| \begin{matrix}
\rho && \frac{e}{c^2} \\
&& \\
\frac{e}{c^2} && \rho \, \theta_{0,2}
\end{matrix}\right|^{-1} \sum_{n_1 =2}^{N} \left( \frac{e}{c^6}  \, U_\alpha U_\beta \, A^{\alpha \beta A_{n_1}}_E - \rho \, \, \theta_{0,2}  \, \, U_\alpha \, A^{\alpha  A_{n_1}}_E \right) \, \lambda_{A_{n_1}} \, , \\
&  U^\mu \, g^{(1)}_\mu =  \left| \begin{matrix}
\rho && \frac{e}{c^2} \\
&& \\
\frac{e}{c^2} && \rho \, \theta_{0,2}
\end{matrix}\right|^{-1} \sum_{n_1 =2}^{N} \left( \frac{e}{c^4}  \, U_\alpha  \, A^{\alpha A_{n_1}}_E - \frac{\rho}{c^4}  \,  U_\alpha \, U_\beta A^{\alpha  \beta A_{n_1}}_E \right) \, \lambda_{A_{n_1}} \, , \\
&  h^{\mu \delta} \, g^{(1)}_\mu   = \frac{3}{\rho \, c^4 \theta_{1,2}} U_\alpha \, h_\beta^\delta \, \sum_{n_1 =2}^{N}  A^{\alpha \beta A_{n_1}}_E \, \lambda_{A_{n_1}} \, , 
\end{split}
\end{align*}
which fully determine  $g^{(1)}$ and $g_\mu^{(1)}$. \\
If we want their second order homogeneous parts, we have to take the derivatives of \eqref{2.1} with respect to $\lambda_{A_{n_2}}$. These derivatives are equivalent to
\begin{align*}
\begin{split}
& \frac{U_\alpha}{- \, k_B} \int_{\Re^3} \int_{0}^{+ \infty} \, f_E \,  F \,  \left[m \, \frac{\partial^2 \, g}{\partial \, \lambda_{A_{n_2}} \, \partial \, \lambda_{A_{n_1}}} + p^\mu \frac{\partial^2 \, g_\mu }{\partial \, \lambda_{A_{n_2}} \, \partial \, \lambda_{A_{n_1}}}  \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) + \frac{\partial^2 \, \Delta \chi}{\partial \, \lambda_{A_{n_2}} \, \partial \, \lambda_{A_{n_1}}}  \right] \cdot \\
& \quad \quad \cdot  p^{\alpha} \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} = - \,  U_\alpha \int_{\Re^3} \int_{0}^{+ \infty} \, f_E \, \frac{\partial \, F}{\partial \, \lambda_{A_{n_2}}}   \, \frac{\partial \, F}{\partial \, \lambda_{A_{n_1}}}  \,  p^{\alpha} \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} \, , \\
& \frac{U_\alpha}{- \, k_B} \int_{\Re^3} \int_{0}^{+ \infty} \, f_E \,  F \,  \left[m \, \frac{\partial^2 \, g}{\partial \, \lambda_{A_{n_2}} \, \partial \, \lambda_{A_{n_1}}} + p^\mu \frac{\partial^2 \, g_\mu }{\partial \, \lambda_{A_{n_2}} \, \partial \, \lambda_{A_{n_1}}}  \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) + \frac{\partial^2 \, \Delta \chi}{\partial \, \lambda_{A_{n_2}} \, \partial \, \lambda_{A_{n_1}}}  \right]  \cdot \\
& \quad \quad \quad \cdot  \,    p^{\alpha} \,   p^{\beta} \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P} = \\
& \quad \quad \quad = - \,  U_\alpha \int_{\Re^3} \int_{0}^{+ \infty} \, f_E \, \frac{\partial \, F}{\partial \, \lambda_{A_{n_2}}}   \, \frac{\partial \, F}{\partial \, \lambda_{A_{n_1}}}  \,  p^{\alpha} \,  p^{\beta} \left( 1+ \, \frac{\mathcal{I}}{m \, c^2} \right) \, \varphi (\mathcal{I}) \, d \, \mathcal{I} \, d \, \vec{P}  
 \, . 
\end{split}
\end{align*}
If we calculate these equations at equilibrium, contract the result  with $\frac{1}{2} \, \lambda_{A_{n_2}} \, \lambda_{A_{n_1}}$ and take the sum for $n_1, n_2 =2,  \, \cdots \, , \, N$ they give a system like \eqref{2.2}, but with $g^{(2)}$ and $g_\mu^{(2)}$ instead of $g^{(1)}$ and $g_\mu^{(1)}$ in the left hand sides, while in the right hand sides there are known functions which include polynomials in  $g^{(1)}$ and $g_\mu^{(1)}$ (already known).
Therefore, by proceeding as we have done for \eqref{2.2}, we  fully determine  $g^{(2)}$ and $g_\mu^{(2)}$. \\
It is obvious that, by taking the higher order  derivatives of \eqref{2.1}  and proceeding in the same way, we obtain $g^{(h)}$ and $g_\mu^{(h)}$ $\forall h$. This completes the proof that \eqref{2.0} jointly with $g^E=0$, $g_\mu^E=0$ give one and only one solution.





















\begin{thebibliography}{99}
\bibitem{Eck}  C. Eckart, The thermodynamics of irreversible processes, III. Relativistic theory of the simple fluid, Phys. Rev. 58 (1940) 919–924.
\bibitem{Lif} L.D. Landau, E.M. Lifshitz, Fluid Mechanics, second ed., Butterworth-Heinemann, Oxford, 1997.

\bibitem{CrisPen} M.C. Carrisi, S. Pennisi, Relativistic extended thermodynamics of polyatomic gases in the Landau and Lifshitz description, Int. Journal of Non-Linear Mechanics 135 (2021) 103756, doi; 10.1016/j.ijnonlinmec.2021.103756. 
\bibitem{CrisPen2} M.C. Carrisi, S. Pennisi, Hyperbolicity of a model for polyatomic gases in relativistic extended thermodynamics, Continuum Mech. Thermodyn. (2020) 32:1435–1454
doi: 10.1007/s00161-019-00857-0. 
\bibitem{Pen} Pennisi, S. Consistent Order Approximations in Extended Thermodynamics of Polyatomic Gases. J. Nat. Sci. Tech. 2021, 2, 12–21.
\bibitem{Bri} Brini, F.: Hyperbolicity region in extended thermodynamics with 14 moments. Continuum Mech. Thermodyn. 13, 1–8 (2001). 
\bibitem{R-T} Ruggeri, T., Trovato, M.: Hyperbolicity in extended thermodynamics of fermi and bose gases. Continuum Mech. Thermodyn. 16, 551–576 (2004)
\bibitem{Bri2} Brini, F., Ruggeri, T. Second-order approximation of extended thermodynamics of a monatomic gas and hyperbolicity region. Continuum Mech. Thermodyn. 32, 23–39 (2020). https://doi.org/10.1007/s00161-019-00778-y
\bibitem{Struc1} H. Struchtrup, Heat transfer in the transition regime: Solution of boundary value problems for Grad’s moment equations via kinetic schemes. PHYSICAL REVIEW E, 65, 041204,1-16 (2002)
\bibitem{Struc2} H. Struchtrup, Stable transport equations for rarefied gases at high orders in the Knudsen number. Phisics of Fluids. 16 (11), 3921-3934 (2004).
\bibitem{DP} Demontis F., Pennisi, S., On two possible ways to recover Ordinary Thermodynamics from Extended Thermodynamics of Polyatomic gases, submitted for publication (2023). 
\bibitem{Annals2017} Pennisi, S.; Ruggeri T. Relativistic extended thermodynamics of rarefied polyatomic gas. Ann. Phys. 2017, 377, 415–445.
\bibitem{ACPR} T. Arima, M. C. Carrisi, S. Pennisi, T. Ruggeri, Relativistic Rational Extended Thermodynamics of Polyatomic Gases with a New Hierarchy of Moments, Entropy 2022, 24, 43. Doi: 10.3390/e24010043. 
\bibitem{Giap} Pennisi, S.; Ruggeri T., A  new BGK model for relativistic kinetic theory
of monatomic and polyatomic gases. IOP Conf. Series: Journal of Physics: Conf. Series 1035 (2018) 012005 doi :10.1088/1742-6596/1035/1/012005
\bibitem{19} S.R. De Groot, W.A. van Leeuven, Ch.G. van Weert, Relativistic Kinetic Theory.
Principles and Applications, North Holland, Amsterdam, New York, Oxford, 1980.
\bibitem{20} L. Rezzolla, O. Zanotti, Relativistic Hydrodynamics, Oxford University Press,
Oxford, 2013.



\end{thebibliography}

\end{document}


Struchtrup H.: The BGK-model with velocity-dependent collision frequency. Continuum Mech. Thermodyn. 9, 23–31 (1997)

Struchtrup H.: On the number of moments in radiative transfer problems. Ann. Phys. 266, 1–26 (1998)

Struchtrup H.: Extended moment method for electrons in semiconductors. Physica A 275, 229–255 (2002)

Struchtrup H.: Heat transfer in the transition regime: Solution of boundary value problems for Grad’s moment equations via kinetic schemes. Phys. Rev. E 65, 041204 (2002)

Struchtrup H.: Stable transport equations for rarefied gas at high order in the Knudsen number. Phys. Fluids 16, 3921–3934 (2004)

Struchtrup H.: Derivation of 13 moment equations for rarefied gas flow to second order accuracy for arbitrary interaction potentials. Multiscale Model Simul. 3, 221–243 (2005)

Published: 03 December 2010
Some open problems in non-linear extended thermodynamics and their possible solutions
M. C. Carrisi & S. Pennisi 
Ricerche di Matematica volume 60, pages45–56 (2011)Cite this article
