Chapman-Enskog expansion for Navier-Stokes equations

We detail here the algebraic calculations of Chapman-Enskog expansion of Lattice Boltzmann Equation for simulating Navier-Stokes (NS) equations. Here, the expansions are inspired from [1] (pp 332–336). An introduction of that method can be found in [2] (page 106). An alternative method, but equivalent to C-E expansions, is presented in [3]. Here we assume that the equilibrium distribution function and the lattice are known to focus on derivation and approximation. The intermediate calculations of Chapman-Enskog expansion are rarely detailed in publications except in appendices in few of them which improve the numerical methods. Here we present the classical procedure without force term in NS.

Starting point

The LBM aims to simulate the discrete Boltzmann equation:

(292)\[f_{i}(\mathbf{x}+\mathbf{e}_{i}\Delta t,\,t+\Delta t)=f_{i}(\mathbf{x},\,t)-\frac{1}{\lambda}\left[f_{i}(\mathbf{x},\,t)-f_{i}^{(0)}(\mathbf{x},\,t)\right]\]

where \(\mathbf{x}=(x,y,z)^{T}\) is the position, \(\Delta x\) is the spatial discretization, \(\Delta t\) is the time-step and \(\lambda\) is a relaxation coefficient of distribution function \(f_{i}(\mathbf{x},\,t)\) toward an equilibrium function \(f_{i}^{(0)}(\mathbf{x},\,t)\). Those functions depend on the lattice speeds \(\mathbf{e}_{i}\) with \(i=1,...,N\). The equilibrium distribution function is defined by:

(293)\[f_{i}^{(0)}(\mathbf{x},\,t)=\rho(\mathbf{x},\,t)w_{i}\left[1+3\frac{(\mathbf{e}_{i}\cdot\mathbf{u}(\mathbf{x},\,t))}{c^{2}}+\frac{9}{2}\frac{(\mathbf{e}_{i}\cdot\mathbf{u}(\mathbf{x},\,t))^{2}}{c^{4}}-\frac{3}{2}\frac{\mathbf{u}^{2}(\mathbf{x},\,t)}{c^{2}}\right]\]

where \(\rho(\mathbf{x},\,t)\) is the density, \(\mathbf{u}=(u_{x},\,u_{y},\,u_{z})^{T}\) is the fluid velocity, \(c=\frac{\Delta x}{\Delta t}\) is the lattice speed and \(w_{i}\) are weights depending on the lattice.

In the remainder of this page, in order to alleviate the notations, the dependences with \(\mathbf{x}\) and \(t\) will be canceled in distribution functions \(f_{i}\) et \(f_{i}^{(0)}\) and macroscopic variables \(\rho\) et \(\mathbf{u}\). We also choose the popular lattice D2Q9 defined by ( for \(i=0,...,8\)):

(294)\[\begin{split}\mathbf{e}_{0}=\left(\begin{array}{c} 0\\ 0 \end{array}\right)\end{split}\]

(295)\[\begin{split}\mathbf{e}_{1}=\left(\begin{array}{c} 1\\ 0 \end{array}\right)c,\,\mathbf{e}_{2}=\left(\begin{array}{c} 0\\ 1 \end{array}\right)c,\,\mathbf{e}_{3}=\left(\begin{array}{c} -1\\ 0 \end{array}\right)c,\,\mathbf{e}_{4}=\left(\begin{array}{c} 0\\ -1 \end{array}\right)c\end{split}\]

(296)\[\begin{split}\mathbf{e}_{5}=\left(\begin{array}{c} 1\\ 1 \end{array}\right)c,\,\mathbf{e}_{6}=\left(\begin{array}{c} -1\\ 1 \end{array}\right)c,\,\mathbf{e}_{7}=\left(\begin{array}{c} -1\\ -1 \end{array}\right)c,\,\mathbf{e}_{8}=\left(\begin{array}{c} 1\\ -1 \end{array}\right)c\end{split}\]

The weights \(w_{i}\) are:

(297)\[w_{0}=\frac{4}{9},\qquad w_{1,\,...,\,4}=\frac{1}{9},\qquad w_{5,\,...,\,8}=\frac{1}{36}\]

Objectives

The purpose of this section is to demonstrate with the Chapman-Enskog expansion that Eq. (292) with equilibrium (293) and lattice D2Q9 (relations (294) – (296) and (297)) simulate the low-Mach Navier-Stokes equations:

(298)\[\partial_{t}\rho+\partial_{\alpha}(\rho u_{\alpha})=0\]

(299)\[\partial_{t}(\rho u_{\alpha})+\partial_{\beta}(\rho u_{\alpha}u_{\beta})=-\partial_{\alpha}p+\partial_{\beta}\left[\nu\rho(\partial_{\beta}u_{\alpha}+\partial_{\alpha}u_{\beta})\right]+\mathcal{O}(\rho u^{3})\]

Eq. (298) is the mass balance equation and (299) and the impulsion balance equation. In those equations, the greek indices are \(\alpha=x,y\) and \(\beta=x,y\) (2D) respectively. The Einstein’s convention of summation aver repeated indices is used (for instance \(\partial_{\alpha}u_{\alpha}=\partial_{x}u_{x}+\partial_{y}u_{y})\). The symbol \(\partial_{\alpha}\) is a notation to express the partial derivative \(\frac{\partial}{\partial\alpha}\). For example \(\partial_{x}\equiv\frac{\partial}{\partial x}\), \(\partial_{y}\equiv\frac{\partial}{\partial y}\) and \(\partial_{t}\equiv\frac{\partial}{\partial t}\). In Eq. (299) \(p\) is the (thermodynamic) pressure, \(\nu\) is the kinematic viscosity.

Conservative quantities and moments of equilibrium distribution functions

The macroscopic quantities which have to be conserved are the density and the impulsion:

(300)\[\begin{split}\begin{cases} \sum_{i}f_{i} & =\rho\\ \sum_{i}f_{i}\mathbf{e}_{i} & =\rho\mathbf{u} \end{cases}\end{split}\]

The conservation of those quantities during a collision implies:

(301)\[\begin{split}\begin{cases} \sum_{i}f_{i}^{(0)} & =\rho\\ \sum_{i}f_{i}^{(0)}\mathbf{e}_{i} & =\rho\mathbf{u} \end{cases}\end{split}\]

Without loss of generality, we set \(\Delta x=\Delta t=1\) (lattice units) in order to simplify the notations. In general case, \(\Delta x\) and \(\Delta t\) can take values different of unity. We get in Eq. (293) \(c=1\). We also set:

(302)\[c_{s}^{2}=\frac{1}{3}\]

With greek indices, the equilibrium distribution function writes:

(303)\[f_{i}^{(0)}=\rho w_{i}\left[1+3e_{i\gamma}u_{\gamma}+\frac{9}{2}e_{i\gamma}e_{i\theta}u_{\gamma}u_{\theta}-\frac{3}{2}u_{\gamma}u_{\theta}\delta_{\gamma\theta}\right]\]

On a D2Q9 lattice, we can show that the moments of order 0, 1, 2 et 3 of \(f_{i}^{(0)}\) are:

(304)\[\sum_{i}f_{i}^{(0)}=\rho\]

(305)\[\sum_{i}f_{i}^{(0)}e_{i\alpha}=\rho u_{\alpha}\]

(306)\[\sum_{i}f_{i}^{(0)}e_{i\alpha}e_{i\beta}=c_{s}^{2}\rho\delta_{\alpha\beta}+\rho u_{\alpha}u_{\beta}\]

(307)\[\sum_{i}f_{i}^{(0)}e_{i\alpha}e_{i\beta}e_{i\gamma}=c_{s}^{2}\rho(u_{\alpha}\delta_{\beta\gamma}+u_{\beta}\delta_{\alpha\gamma}+u_{\gamma}\delta_{\alpha\beta})\]

Those results will be used in calculations of following section.

Remark

When we seek to simulate a specific macroscopic PDE, we need to establish the equilibrium distribution function \(f_{i}^{(0)}\) and the number of moving directions, such as the moments of order zero and one \(f_{i}^{(0)}\) are equal to conserved macroscopic quantities.

Details of Chapman-Enskog expansion

First part: Taylor’s expansion and scale separation

The expansions presented in this first part are classical and can be found in references cited above. We can refer to them for physical interpretation of small parameter \(\varepsilon\) in kinetic theory of gases (Knudsen numner).

By setting \(\Delta x=\Delta t=1\) and using the greek indices, the LBE writes:

(308)\[f_{i}(x+e_{i\beta},\,t+1)=f_{i}(x,\,t)-\frac{1}{\lambda}\left[f_{i}(x,\,t)-f_{i}^{(0)}(x,\,t)\right]\]

The Taylor’s expansion of left-hand side of Eq. (308) up to second order in space and times yields:

\[f_{i}(x+e_{i\beta},\,t+1)\simeq f_{i}+e_{i\beta}\partial_{\beta}f_{i}+\frac{1}{2}e_{i\alpha}e_{i\beta}\partial_{\alpha}\partial_{\beta}f_{i}+\partial_{t}f_{i}+\frac{1}{2}\partial_{t}^{2}f_{i}+e_{i\beta}\partial_{\beta}\partial_{t}f_{i}\]

By replacing in Eq. (308), we obtain:

(309)\[e_{i\beta}\partial_{\beta}f_{i}+\frac{1}{2}e_{i\alpha}e_{i\beta}\partial_{\alpha}\partial_{\beta}f_{i}+\partial_{t}f_{i}+\frac{1}{2}\partial_{t}^{2}f_{i}+e_{i\beta}\partial_{\beta}\partial_{t}f_{i}=-\frac{1}{\lambda}\left[f_{i}-f_{i}^{(0)}\right]\]

We perform a scale separation in space (\(x_{1}=\varepsilon x\)) and time. For the latter, we distinguish one characteristic time of convection \(t_{0}=\varepsilon t\) and one characteristic time of diffusion \(t_{1}=\varepsilon^{2}t\) such that:

(310)\[\begin{split}\begin{cases} t & =\frac{1}{\varepsilon}t_{0}+\frac{1}{\varepsilon^{2}}t_{1}\\ x & =\frac{1}{\varepsilon}x_{1} \end{cases}\end{split}\]

That leads to following partial derivatives:

(311)\[\begin{split}\begin{cases} \partial_{t} & =\varepsilon\partial_{t_{0}}+\varepsilon^{2}\partial_{t_{1}}\\ \partial_{\alpha} & =\varepsilon\partial_{\alpha}^{(1)} \end{cases}\end{split}\]

The distribution function is expanded up to the second order:

(312)\[f_{i}\simeq f_{i}^{(0)}+\varepsilon f_{i}^{(1)}+\varepsilon f_{i}^{(2)}\]

That expansion and Eqs. (300) and (301) imply that:

(313)\[\sum_{i}f_{i}^{(1)}=\sum_{i}f_{i}^{(2)}=0\]

(314)\[\sum_{i}f_{i}^{(1)}e_{i\alpha}=\sum_{i}f_{i}^{(2)}e_{i\alpha}=0\]

Remark

Eq. (314) est derived because the impulsion is a conservative local quantity for Navier-Stokes equations. For advection-diffusion equation, this is not any more the case, and it will be necessary to calculate the first order moment of \(f_{i}^{(1)}\).

Moments des termes en \(\varepsilon\)

The partial derivative in space and times are replaced in (309) by using Eq (311) and we also use the expansion of distribution function iwth Eq. (312). We gather the terms of order \(\varepsilon\) and \(\varepsilon^{2}\).

The gathering of first order terms in \(\varepsilon\) yields the following equality:

(315)\[\partial_{\alpha}^{(1)}e_{i\alpha}f_{i}^{(0)}+\partial_{t_{0}}f_{i}^{(0)}=-\frac{1}{\lambda}f_{i}^{(1)}\]

With preliminary results Eqs. (304) and (305), the moment of zeroth-order (sum over \(i\)) of that equation is:

(316)\[\partial_{t_{0}}\rho+\partial_{\alpha}^{(1)}(\rho u_{\alpha})=0\]

Similarily, by using Eqs. (305) and (306), the moment of first order (multiplication of Eq. (315) by \(e_{i\beta}\) and sum over \(i\)) leads to:

(317)\[\partial_{t_{0}}(\rho u_{\alpha})+\partial_{\beta}^{(1)}(c_{s}^{2}\rho\delta_{\alpha\beta}+\rho u_{\alpha}u_{\beta})=0\]

Those two moments of order 0 and 1 of terms in \(\varepsilon\) (Eqs. (316) and (317)) will be used for calculation of the 2nd order moment of \(f_{i}^{(1)}\) in the second part of this expansion.

Moments of terms of order \(\varepsilon^{2}\)

Gathering the \(\varepsilon^{2}\) terms leads to the relation:

(318)\[e_{i\beta}\partial_{\beta}^{(1)}f_{i}^{(1)}+\frac{1}{2}e_{i\alpha}e_{i\beta}\partial_{\alpha}^{(1)}\partial_{\beta}^{(1)}f_{i}^{(0)}+\partial_{t_{0}}f_{i}^{(1)}+\partial_{t_{1}}f_{i}^{(0)}+\frac{1}{2}\partial_{t_{0}}^{2}f_{i}^{(0)}+e_{i\beta}\partial_{\beta}^{(1)}\partial_{t_{0}}f_{i}^{(0)}=-\frac{1}{\lambda}f_{i}^{(2)}\]

Applying separately \(\frac{1}{2}e_{i\beta}\partial_{\beta}\) and \(\frac{1}{2}\partial_{t_{0}}\) on Eq. (315) and summing those two results, that Eq. (318) se simplifies into:

\[\partial_{t_{1}}f_{i}^{(0)}+\left(1-\frac{1}{2\lambda}\right)\left[\partial_{\beta}^{(1)}e_{i\beta}f_{i}^{(1)}+\partial_{t_{0}}f_{i}^{(1)}\right]=-\frac{1}{\lambda}f_{i}^{(2)}\]

With equations (304), (313) and (314), the moment of order 0 of that equation yields:

(319)\[\partial_{t_{1}}\rho=0\]

and its first-order moment:

(320)\[\partial_{t_{1}}(\rho u_{\alpha})+\left(1-\frac{1}{2\lambda}\right)\left[\partial_{\beta}^{(1)}\sum_{i}f_{i}^{(1)}e_{i\alpha}e_{i\beta}\right]=0\]

Gathering moments

We gather the moments of zeroth- and first-order that have been derived separately for terms in \(\varepsilon\) and \(\varepsilon^{2}\). We carry out \(\varepsilon\times\)Eq (316) + \(\varepsilon^{2}\times\)Eq (319), the moment of order zero is:

(321)\[\partial_{t}\rho+\partial_{\alpha}(\rho u_{\alpha})=0\]

which is the first equation of NS system i.e. the mass balance equation. We perform now \(\varepsilon\times\)Eq (317) + \(\varepsilon^{2}\times\)Eq (320), the first-order moment is:

(322)\[\partial_{t}(\rho u_{\alpha})+\partial_{\beta}(\Pi_{\alpha\beta}^{(0)}+\Pi_{\alpha\beta}^{(1)})=0\]

where we set:

(323)\[\Pi_{\alpha\beta}^{(0)}=\sum_{i}f_{i}^{(0)}e_{i\alpha}e_{i\beta}\]

(324)\[\Pi_{\alpha\beta}^{(1)}=\left(1-\frac{1}{2\lambda}\right)\varepsilon\sum_{i}f_{i}^{(1)}e_{i\alpha}e_{i\beta}\]

The notations \(\Pi_{\alpha\beta}^{(0)}\) and \(\Pi_{\alpha\beta}^{(1)}\) are the second-order moments of \(f_{i}^{(0)}\) and \(f_{i}^{(1)}\) respectively. \(\Pi_{\alpha\beta}^{(0)}\) is known with definition of EDF and preliminary results (306). We have to derive now the second-order moment of \(f_{i}^{(1)}\). It is the purpose of the following second part.

Second part: viscosity coefficient and impulsion balance

Viscosity

For \(f_{i}^{(1)}\) (\(\sum_{i}f_{i}^{(1)}e_{i\alpha}e_{i\beta}\)), we use the preliminary results (306) and (307) for moments of second- and third-order of EDF \(f_{i}^{(0)}\):

\[\begin{split}\begin{aligned} \sum_{i}f_{i}^{(1)}e_{i\alpha}e_{i\beta} & =-\lambda\left[\partial_{\gamma}^{(1)}(\sum_{i}f_{i}^{(0)}e_{i\alpha}e_{i\beta}e_{i\gamma})+\partial_{t_{0}}(\sum_{i}f_{i}^{(0)}e_{i\alpha}e_{i\beta})\right]\\ & =-\lambda\left[\partial_{\gamma}^{(1)}c_{s}^{2}\rho(u_{\alpha}\delta_{\beta\gamma}+u_{\beta}\delta_{\alpha\gamma}+u_{\gamma}\delta_{\alpha\beta})+\partial_{t_{0}}(c_{s}^{2}\rho\delta_{\alpha\beta}+\rho u_{\alpha}u_{\beta})\right]\\ & =-\lambda c_{s}^{2}\left[\partial_{\beta}^{(1)}(\rho u_{\alpha})+\partial_{\alpha}^{(1)}(\rho u_{\beta})+\partial_{\gamma}^{(1)}(\rho u_{\gamma})\delta_{\alpha\beta}\right]-\lambda\partial_{t_{0}}(c_{s}^{2}\rho\delta_{\alpha\beta})-\lambda\partial_{t_{0}}(\rho u_{\alpha}u_{\beta})\end{aligned}\end{split}\]

For last two terms, the time derivative are transformed into space-derivative by using moments of order zero and one of equation of order \(\varepsilon\) (Eq. (316) and (317) respectively). For the penultimate term we obtain:

\[-\lambda\partial_{t_{0}}(c_{s}^{2}\rho\delta_{\alpha\beta})=+\lambda c_{s}^{2}\partial_{\gamma}^{(1)}(\rho u_{\gamma})\delta_{\alpha\beta}\]

That term cancels with one of opposite sign contained in the bracket. For last term, we use the following result:

(325)\[\partial_{t_{0}}(\rho u_{\alpha}u_{\beta})=-\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\beta}u_{\gamma})-c_{s}^{2}(u_{\alpha}\partial_{\beta}^{(1)}\rho+u_{\beta}\partial_{\alpha}^{(1)}\rho)\]

By expansing the partial derivatives in the bracket, we finally get:

\[\begin{split}\begin{aligned} \sum_{i}f_{i}^{(1)}e_{i\alpha}e_{i\beta} & =-\lambda c_{s}^{2}\left[\partial_{\beta}^{(1)}(\rho u_{\alpha})+\partial_{\alpha}^{(1)}(\rho u_{\beta})\right]+\lambda\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\beta}u_{\gamma})+\lambda c_{s}^{2}(u_{\alpha}\partial_{\beta}^{(1)}\rho+u_{\beta}\partial_{\alpha}^{(1)}\rho)\\ & =-\lambda c_{s}^{2}\left[\rho\partial_{\beta}^{(1)}u_{\alpha}+\rho\partial_{\alpha}^{(1)}u_{\beta}\right]+\lambda\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\beta}u_{\gamma}) \end{aligned}\end{split}\]

The term \(\Pi_{\alpha\beta}^{(1)}\) becomes:

\[\begin{split}\begin{aligned} \Pi_{\alpha\beta}^{(1)} & =\left(1-\frac{1}{2\lambda}\right)\varepsilon\sum_{i}f_{i}^{(1)}e_{i\alpha}e_{i\beta}\\ & =\left(1-\frac{1}{2\lambda}\right)\varepsilon\left\{ -\lambda c_{s}^{2}\left[\rho\partial_{\beta}^{(1)}u_{\alpha}+\rho\partial_{\alpha}^{(1)}u_{\beta}\right]+\lambda\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\beta}u_{\gamma})\right\} \\ & =c_{s}^{2}\left(\lambda-\frac{1}{2}\right)\left\{ -\left[\rho\partial_{\beta}u_{\alpha}+\rho\partial_{\alpha}u_{\beta}\right]+\frac{1}{c_{s}^{2}}\partial_{\gamma}(\rho u_{\alpha}u_{\beta}u_{\gamma})\right\} \end{aligned}\end{split}\]

By setting:

(326)\[\nu=c_{s}^{2}\left(\lambda-\frac{1}{2}\right)\]

we obtain:

(327)\[\Pi_{\alpha\beta}^{(1)}=-\nu\left[\rho\partial_{\beta}u_{\alpha}+\rho\partial_{\alpha}u_{\beta}\right]+\frac{1}{c_{s}^{2}}\nu\partial_{\gamma}(\rho u_{\alpha}u_{\beta}u_{\gamma})\]

L’équation (326) donne la relation entre le paramètre de relaxation \(\lambda\) de l’équation cinétique de Boltzmann et la viscosité cinématique \(\nu\) du fluide.

Remarque

: en toute rigueur, lors du développement de Taylor, il est nécessaire de considérer \(\Delta x\neq1\) et \(\Delta t\neq1\) et on voit apparaître le rapport \(\Delta x^{2}/\Delta t\) en facteur devant le terme entre parenthèses :

\[\nu=c_{s}^{2}\left(\lambda-\frac{1}{2}\right)\frac{\Delta x^{2}}{\Delta t}\]

Impulsion balance equation

En remplaçant (306) et (327) dans (322), on obtient :

\[\partial_{t}(\rho u_{\alpha})+\partial_{\beta}\left[c_{s}^{2}\rho\delta_{\alpha\beta}+\rho u_{\alpha}u_{\beta}-\nu(\rho\partial_{\beta}u_{\alpha}+\rho\partial_{\alpha}u_{\beta})\right]+\frac{1}{c_{s}^{2}}\partial_{\beta}\nu\partial_{\gamma}(\rho u_{\alpha}u_{\beta}u_{\gamma})=0\]

\[\partial_{t}(\rho u_{\alpha})+\partial_{\beta}(\rho u_{\alpha}u_{\beta})=-\partial_{\alpha}(c_{s}^{2}\rho)+\partial_{\beta}\left[\nu(\rho\partial_{\beta}u_{\alpha}+\rho\partial_{\alpha}u_{\beta})\right]+\mathcal{O}(\rho u^{3})\]

où on a considéré \(-\frac{1}{c_{s}^{2}}\partial_{\beta}\nu\partial_{\gamma}(\rho u_{\alpha}u_{\beta}u_{\gamma})=\mathcal{O}(\rho u^{3})\). Finalement, en posant \(p=c_{s}^{2}\rho\), on obtient l’équation de conservation de la quantité de mouvement :

\[\partial_{t}(\rho u_{\alpha})+\partial_{\beta}(\rho u_{\alpha}u_{\beta})=-\partial_{\alpha}p+\partial_{\beta}\left[\nu\rho(\partial_{\beta}u_{\alpha}+\partial_{\alpha}u_{\beta})\right]+\mathcal{O}(\rho u^{3})\]

qui correspond à la seconde équation recherchée du système d’équations de Navier-Stokes.

Récapitulatif

En effectuant le développement de Chapman-Enskog (ordre deux en espace, ordre deux en temps et en \(\varepsilon^{2}\)) et en calculant les moments d’ordre 0 et 1 des équations en \(\varepsilon\) et \(\varepsilon^{2}\), on a montré que les équations (292) et (293) permettent de résoudre sur un réseau D2Q9 les équations de Navier-Stokes :

(328)\[\partial_{t}\rho+\partial_{\alpha}(\rho u_{\alpha})=0\]

(329)\[\partial_{t}(\rho u_{\alpha})+\partial_{\beta}(\rho u_{\alpha}u_{\beta})=-\partial_{\alpha}p+\partial_{\beta}\left[\nu\rho(\partial_{\beta}u_{\alpha}+\partial_{\alpha}u_{\beta})\right]+\mathcal{O}(u^{3})\]

avec la pression reliée à la densité par la relation \(p=c_{s}^{2}\rho\) et une viscosité cinématique \(\nu\) donnée par :

(330)\[\nu=c_{s}^{2}\left(\lambda-\frac{1}{2}\right)\frac{\Delta x^{2}}{\Delta t}\]

Dans l’équation (329), le terme en \(\mathcal{O}(u^{3})\) est de la forme \(-\frac{1}{c_{s}^{2}}\partial_{\beta}\nu\partial_{\gamma}(\rho u_{\alpha}u_{\beta}u_{\gamma})\) qui est négligeable pour les faibles nombres de Mach. Lorsque ce dernier devient important, il est nécessaire de corriger la fonction de distribution à l’équilibre pour éliminer/diminuer l’influence de ce terme. On pourra se référer à [4] qui analyse plusieurs \(f_{i}^{(0)}\) ou termes forces en ce sens.

Démonstration de l’égalité (325)

On démontre dans cette annexe la relation :

\[\partial_{t_{0}}(\rho u_{\alpha}u_{\beta})=-\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\beta}u_{\gamma})-c_{s}^{2}(u_{\alpha}\partial_{\beta}^{(1)}\rho+u_{\beta}\partial_{\alpha}^{(1)}\rho)\]

Le développement du membre de gauche donne :

\[\partial_{t_{0}}(\rho u_{\alpha}u_{\beta})=u_{\alpha}u_{\beta}\partial_{t_{0}}\rho+\rho u_{\alpha}\partial_{t_{0}}u_{\beta}+\rho u_{\beta}\partial_{t_{0}}u_{\alpha}\]

En remarquant que \(\rho u_{\alpha}\partial_{t_{0}}u_{\beta}=u_{\alpha}\partial_{t_{0}}(\rho u_{\beta})-u_{\alpha}u_{\beta}\partial_{t_{0}}\rho\) et que \(\rho u_{\beta}\partial_{t_{0}}u_{\alpha}=u_{\beta}\partial_{t_{0}}(\rho u_{\alpha})-u_{\alpha}u_{\beta}\partial_{t_{0}}\rho\), on obtient :

\[\partial_{t_{0}}(\rho u_{\alpha}u_{\beta})=-u_{\alpha}u_{\beta}\partial_{t_{0}}\rho+u_{\alpha}\partial_{t_{0}}(\rho u_{\beta})+u_{\beta}\partial_{t_{0}}(\rho u_{\alpha})\]

Les dérivées partielles en temps sont converties en dérivées partielles en espace à l’aide des équations (316) pour le premier terme et (317) pour les deux derniers :

(331)\[\begin{split}\partial_{t_{0}}(\rho u_{\alpha}u_{\beta}) &=u_{\alpha}u_{\beta}\partial_{\gamma}^{(1)}(\rho u_{\gamma})-u_{\alpha}\left[\partial_{\gamma}^{(1)}(c_{s}^{2}\rho\delta_{\beta\gamma})+\partial_{\gamma}^{(1)}(\rho u_{\beta}u_{\gamma})\right]-u_{\beta}\left[\partial_{\gamma}^{(1)}(c_{s}^{2}\rho\delta_{\alpha\gamma})+\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\gamma})\right]\\ &=u_{\alpha}u_{\beta}\partial_{\gamma}^{(1)}(\rho u_{\gamma})-u_{\alpha}\left[\partial_{\beta}^{(1)}(c_{s}^{2}\rho)+\partial_{\gamma}^{(1)}(\rho u_{\beta}u_{\gamma})\right]-u_{\beta}\left[\partial_{\alpha}^{(1)}(c_{s}^{2}\rho)+\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\gamma})\right]\end{split}\]

Les trois premiers termes contenus dans le premier crochet peuvent être exprimés sous la forme condensée \(\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\beta}u_{\gamma})\). Il existe plusieurs façons de la démontrer, on présente dans cette annexe une façon de procéder.

Un premier développement donne \(\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\beta}u_{\gamma})=u_{\alpha}\partial_{\gamma}^{(1)}(\rho u_{\beta}u_{\gamma})+\rho u_{\beta}u_{\gamma}(\partial_{\gamma}^{(1)}u_{\alpha})\) c’est-à-dire :

(332)\[u_{\alpha}\partial_{\gamma}^{(1)}(\rho u_{\beta}u_{\gamma})=\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\beta}u_{\gamma})-\rho u_{\beta}u_{\gamma}(\partial_{\gamma}^{(1)}u_{\alpha})\]

Un second développement conduit à \(\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\beta}u_{\gamma})=u_{\beta}\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\gamma})+\rho u_{\alpha}u_{\gamma}(\partial_{\gamma}^{(1)}u_{\beta})\) c’est-à-dire :

(333)\[u_{\beta}\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\gamma})=\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\beta}u_{\gamma})-\rho u_{\alpha}u_{\gamma}(\partial_{\gamma}^{(1)}u_{\beta})\]

Enfin, un dernier développement donne \(\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\beta}u_{\gamma})=u_{\alpha}u_{\beta}\partial_{\gamma}^{(1)}(\rho u_{\gamma})+\rho u_{\gamma}\partial_{\gamma}^{(1)}(u_{\alpha}u_{\beta})\) c’est-à-dire :

(334)\[-u_{\alpha}u_{\beta}\partial_{\gamma}^{(1)}(\rho u_{\gamma})=-\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\beta}u_{\gamma})+\rho u_{\gamma}\partial_{\gamma}^{(1)}(u_{\alpha}u_{\beta})\]

En effectuant la somme (332), (333) et (334), on obtient :

\[u_{\alpha}\partial_{\gamma}^{(1)}(\rho u_{\beta}u_{\gamma})+u_{\beta}\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\gamma})-u_{\alpha}u_{\beta}\partial_{\gamma}^{(1)}(\rho u_{\gamma})=\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\beta}u_{\gamma})-\rho u_{\beta}u_{\gamma}(\partial_{\gamma}^{(1)}u_{\alpha})-\rho u_{\alpha}u_{\gamma}(\partial_{\gamma}^{(1)}u_{\beta})+\rho u_{\gamma}\partial_{\gamma}^{(1)}(u_{\alpha}u_{\beta})\]

Dans l’équation ci-dessus, les trois derniers termes du membre de droite s’annulent entre eux, on a donc bien démontré que :

\[\left[u_{\alpha}\partial_{\gamma}^{(1)}(\rho u_{\beta}u_{\gamma})+u_{\beta}\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\gamma})-u_{\alpha}u_{\beta}\partial_{\gamma}^{(1)}(\rho u_{\gamma})\right]=\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\beta}u_{\gamma})\]

Finalement en remplaçant ce résultat dans (331), on obtient la relation recherchée :

\[\partial_{t_{0}}(\rho u_{\alpha}u_{\beta})=-\partial_{\gamma}^{(1)}(\rho u_{\alpha}u_{\beta}u_{\gamma})-c_{s}^{2}\left[u_{\alpha}\partial_{\beta}^{(1)}\rho+u_{\beta}\partial_{\alpha}^{(1)}\rho\right]\]

References

[1]

Shiyi Chen and Gary D. Doolen. Lattice boltzmann method for fluid flows. Annual Review of Fluid Mechanics, 30(Volume 30, 1998):329–364, 1998. URL: https://www.annualreviews.org/content/journals/10.1146/annurev.fluid.30.1.329, doi:https://doi.org/10.1146/annurev.fluid.30.1.329.

[2]

Timm Krüger, Halim Kusumaatmaja, Alexandr Kuzmin, Orest Shardt, Goncalo Silva, Viggen, and Erlend Magnus. The Lattice Boltzmann Method: Principles and Practice. Graduate Texts in Physics. Springer, 2016. ISBN 978-3-319-44647-9.

[3]

François Dubois. Equivalent partial differential equations of a lattice boltzmann scheme. Computers & Mathematics with Applications, 55(7):1441–1449, 2008. Mesoscopic Methods in Engineering and Science. URL: https://www.sciencedirect.com/science/article/pii/S0898122107006268, doi:https://doi.org/10.1016/j.camwa.2007.08.003.

[4]

X.-Y. Lu Huang H., M. Sukop. Multiphase Lattice Boltzmann Methods. Theory and Application. Wiley & Sons, 2015. ISBN 978-1-118-97145-1.

Section author: Alain Cartalade