**5. Nine Velocity DV Boltzmann Model for** *D***2***Q***9 Lattice**

We consider the DV Boltzmann model on a nine-velocity lattice (Figure 1). This lattice is popular in LB theory [3], since the corresponding LB model recovers hydrodynamics at small Mach numbers limit and, in addition, its numerical implementation is very simple. For this model, we have three types of discrete velocities: zero velocity *c*<sup>0</sup> = (0, 0) with the weight *<sup>w</sup>*<sup>0</sup> <sup>=</sup> 16/36; four velocities, parallel to *<sup>x</sup>*, *<sup>y</sup>* axes, i.e., *<sup>c</sup>*±<sup>1</sup> = (±1, 0)*c*, *<sup>c</sup>*±<sup>2</sup> = (0, <sup>±</sup>1)*<sup>c</sup>* with the weight *<sup>w</sup>*<sup>0</sup> <sup>=</sup> 4/36; four diagonal velocities *<sup>c</sup>*±<sup>3</sup> = (±1, <sup>±</sup>1)*c*, *<sup>c</sup>*±<sup>4</sup> = (±1, <sup>∓</sup>1)*<sup>c</sup>* with the weight *w*<sup>0</sup> = 1/36—here, *c* is the positive constant. The lattice velocity magnitudes for these three groups are 0, *c*, <sup>√</sup>2*c*. Moreover, *<sup>θ</sup>*<sup>0</sup> <sup>=</sup> <sup>∑</sup>*<sup>i</sup> wic*<sup>2</sup> *<sup>i</sup>* = *<sup>c</sup>*2/3.

**Figure 1.** Two-dimensional nine-velocity lattice (*D*2*Q*9). Lattice velocities are labeled by red color.

It is well-known that these lattice velocities and weights satisfy the conditions (3), (20) and (21); therefore, if it is possible to construct the collisions in such a way that the mass and momentum are conserved then the Euler equations are satisfied. We mention that the lattices and collision rules for DV Boltzmann models, which can potentially recover the hydrodynamics, have been considered previously [46,47]—for instance, the model with single-relaxation time describing Navier–Stokes equations has been proposed [46,47]. In here we consider only the collisions for the nine-bit lattice in a more detailed way; the considered model is of the multiple-relaxation-time type:

a. **Broadwell type collision** is the reaction between the particles 1 and −1, which turn into the particles 2 and −2 (Figure 1); schematically, we can denote this reaction as (1, −1) −→ (2, −2). The contribution of this collision to right side of (1) denoted as *J*<sup>0</sup> is as follows

$$\text{Jo} = f \text{--} 2f \text{--} f \text{--}\_{-1} f\_1;\tag{25}$$

b. **the collisions linking all three different energy states**, they define transitions between the particle's states with different kinetic energies, and evidently can not be excluded. We have four different reactions (1, 2) −→ (0, 3), (1, −2) −→ (0, 4), (−1, −2) −→ (0, −3), (−1, 2) −→ (0, −4). The corresponding contributions to the collision kernel are

$$\begin{aligned} f\_1 &= f\_0 f\_3 - f\_1 f\_{2\prime} & f\_2 &= f\_0 f\_4 - f\_1 f\_{-2\prime} \\ f\_3 &= f\_0 f\_{-4} - f\_{-1} f\_{2\prime} & f\_4 &= f\_0 f\_{-3} - f\_{-1} f\_{-2\prime} \end{aligned} \tag{26}$$

c. **Broadwell type collision between the particles with the velocity magnitudes** <sup>√</sup>2*<sup>c</sup>* is defined by the reaction (3, −3) −→ (4, −4), the contributions to the collision kernel are

$$J\_5 = f\_{-4}f\_4 - f\_{-3}f\_3;\tag{27}$$

d. **the collisions between the particles with the velocity magnitudes** <sup>√</sup>2*<sup>c</sup>* **and** *<sup>c</sup>*, we have four different reactions (−4, 1) −→ (−1, 3),(−3, 1) −→ (−1, 4),(−3, 2) −→ (−4, −2),(4, 2) −→ (−2, 3), the contributions to the collision kernel are

$$\begin{aligned} \mathbf{J}\_6 &= f\_3 f\_{-1} - f\_{-4} f\_{1\prime} & \mathbf{J}\_7 &= f\_{-1} f\_4 - f\_{-3} f\_{1\prime} \\ \mathbf{J}\_8 &= f\_{-4} f\_{-2} - f\_{-3} f\_{2\prime} & \mathbf{J}\_9 &= f\_3 f\_{-2} - f\_4 f\_{2\prime} \end{aligned} \tag{28}$$

The collisions (25)–(28) conserve mass, momentum and energy; the corresponding *D*2*Q*9 DV Boltzmann model reads as

$$\frac{\partial f\_1}{\partial t} + c \frac{\partial f\_1}{\partial \mathbf{x}} = a f\_0 + \beta (f\_1 + f\_2) + \lambda (f\_6 + f\_7), \tag{29}$$

$$\frac{\partial f\_{-1}}{\partial t} - c \frac{\partial f\_{-1}}{\partial x} = a f\_0 + \beta (f\_3 + f\_4) - \lambda (f\_6 + f\_7), \tag{30}$$

$$\frac{\partial f\_2}{\partial t} + c \frac{\partial f\_2}{\partial y} = -\alpha f\_0 + \beta (f\_1 + f\_3) + \lambda (f\_8 + f\_9), \tag{31}$$

$$\frac{\partial f\_{-2}}{\partial t} - \mathcal{c} \frac{\partial f\_{-2}}{\partial y} = -\kappa f\_0 + \mathcal{J}(f\_2 + f\_4) - \lambda(f\_8 + f\_9),\tag{32}$$

$$c\frac{\partial f\_3}{\partial t} + c\frac{\partial f\_3}{\partial x} + c\frac{\partial f\_3}{\partial y} = \gamma f\_5 - \beta f\_1 - \lambda (f\_6 + f\_9),\tag{33}$$

$$\frac{\partial f\_{-3}}{\partial t} - \mathcal{c} \frac{\partial f\_{-3}}{\partial x} - \mathcal{c} \frac{\partial f\_{-3}}{\partial y} = \gamma f\_5 - \beta f\_4 + \lambda \left(f\_7 + f\_8\right), \tag{34}$$

$$\frac{\partial f\_4}{\partial t} + c \frac{\partial f\_4}{\partial x} - c \frac{\partial f\_4}{\partial y} = -\gamma f\_5 - \beta f\_2 + \lambda (-f\_7 + f\_9), \tag{35}$$

$$\frac{\partial f\_{-4}}{\partial t} - c \frac{\partial f\_{-4}}{\partial x} + c \frac{\partial f\_{-4}}{\partial y} = -\gamma l \mathfrak{s} - \beta l \mathfrak{s} + \lambda (J\mathfrak{s} - J\mathfrak{s}),\tag{36}$$

$$\frac{\partial f\_0}{\partial t} = -\beta (f\_1 + f\_2 + f\_3 + f\_4),\tag{37}$$

where *α*, *β*, *λ*, *γ* in (29)–(37) are positive transition probabilities. Now, we can consider the analogs of the Navier–Stokes equations for the model (29)–(37).

**Proposition 3.** *The Equations (29)–(37) lead to Navier–Stokes equations for nearly incompressible flows with errors of order O*(Δ3) *if*

$$4\alpha = \gamma + 4\beta + 4\lambda\_\prime \tag{38}$$

*the shear viscosity ν equals*

$$\nu = \frac{3}{4\alpha}.\tag{39}$$

**Proof.** From (22), one can deduce that the corrections to the DV distribution function *f* (1) *i* corresponding to the viscous terms can be represented as a linear combination of *d f eq i dt* terms. In the case of nearly incompressible flow, these terms can be represented as (Formula (2.12) in [2])

$$\frac{df\_i^{eq}}{dt} = w\_i \frac{\mathbf{c}\_i \mathbf{c}\_i}{\theta\_0} : \frac{\partial}{\partial r} \Delta \mathbf{u}\_\prime \tag{40}$$

where *<sup>∂</sup> <sup>∂</sup><sup>r</sup>* = ( *<sup>∂</sup> <sup>∂</sup><sup>x</sup>* , *<sup>∂</sup> <sup>∂</sup><sup>y</sup>* ). According (23) we can try add the terms proportional *div*(Δ*u*), but they equal zero for the incompressible limit; then, we seek the solution in the form

$$f\_i^{(1)} = a\_i \mathbf{Q}\_{i\prime} \quad \mathbf{Q}\_i = w\_i \frac{\mathbf{c}\_i \mathbf{c}\_i}{\theta\_0} : \frac{\partial}{\partial r} \Delta \mathbf{u}\_{\prime} \tag{41}$$

where the coefficients *ai* are equal for the indexes *i* corresponding to the discrete velocities *ci* with the same kinetic energy. The substitution of (41) into (29)–(37) leads to three algebraic equations for the coefficients *ai*, (29)–(32) yield the first equation

$$3w\_1 \frac{\partial}{\partial \mathfrak{x}} \Delta u\_{\mathfrak{x}} = 2aw\_1 a\_1 \left(\frac{\partial}{\partial y} \Delta u\_{\mathfrak{y}} - \frac{\partial}{\partial \mathfrak{x}} \Delta u\_{\mathfrak{x}}\right) = $$

$$= 2aw\_1 a\_1 \left(\frac{\partial}{\partial y} \Delta u\_{\mathfrak{y}} - \frac{\partial}{\partial \mathfrak{x}} \Delta u\_{\mathfrak{x}} - \operatorname{div}(\Delta \mathfrak{u})\right) = -4aw\_1 a\_1 \frac{\partial}{\partial \mathfrak{x}} \Delta u\_{\mathfrak{x}}.$$

from which we obtain

$$a\_1 = -\frac{3}{4\alpha'} $$

(33)–(36) yield the second equation

$$2\Im w\_2 \left(\frac{\partial}{\partial \mathbf{x}} \Delta u\_y + \frac{\partial}{\partial y} \Delta u\_x\right) = -(4\gamma w\_2 + \beta w\_0 + 4\lambda w\_1)a\_2 \left(\frac{\partial}{\partial \mathbf{x}} \Delta u\_y + \frac{\partial}{\partial y} \Delta u\_x\right),$$

then

$$a\_2 = -\frac{3w\_2}{4\gamma w\_2 + \beta w\_0 + 4\lambda w\_1}.$$

The third equation, which can be obtained from (37) is satisfied automatically. Now, with the exact expressions for *a*1, *a*2, we can evaluate *f* (1) *<sup>i</sup>* and the viscous corrections to the pressure tensor *<sup>P</sup>*(1) = <sup>∑</sup>*<sup>i</sup> <sup>f</sup>* (1) *<sup>i</sup> cici*. Then, the Navier–Stokes viscous terms can be evaluated as

$$-\sum\_{\sigma} \frac{\partial}{\partial r\_{\sigma}} P^{(1)}\_{\eta \sigma} = -\sum\_{\sigma} \frac{\partial}{\partial r\_{\sigma}} \left( \sum\_{i} f^{(1)}\_{i} c\_{i,\eta} c\_{i,\sigma} \right) \,,$$

where *σ*, *η* equal *x* or *y*. For instance,

$$-\frac{\partial}{\partial x}P\_{xx}^{(1)} - \frac{\partial}{\partial y}P\_{xy}^{(1)} = \frac{3}{2a}\frac{\partial^2}{\partial x^2}\Delta u\_x + \frac{12w\_2}{4\gamma w\_2 + \beta w\_0 + 4\lambda w\_1} \left(\frac{\partial^2}{\partial x \partial y}\Delta u\_y + \frac{\partial^2}{\partial y^2}\Delta u\_x\right)$$

we require 4*α* = *γ* + 4*β* + 4*λ*, then by applying *div*(Δ*u*) = 0 we finally obtain

$$-\frac{\partial}{\partial \boldsymbol{x}} P\_{\boldsymbol{\chi}\boldsymbol{x}}^{(1)} - \frac{\partial}{\partial \boldsymbol{y}} P\_{\boldsymbol{\chi}\boldsymbol{y}}^{(1)} = \frac{3}{4a} \left( \frac{\partial^2}{\partial \boldsymbol{x}^2} + \frac{\partial^2}{\partial \boldsymbol{y}^2} \right) \Delta \boldsymbol{u}\_{\boldsymbol{\chi}\boldsymbol{\cdot}}$$

therefore *ν* = 3/4*α*.

For the model (29)–(37), there are ten collisions. If we consider all reaction vectors and the corresponding collision matrix, one can convince that *rank*(*C*) = 5, the number of the discrete velocities *N* = 9. This means that we do not have any collision invariants except mass, momentum, energy (Proposition 2). We can exclude up to five reactions from the model; for instance, we can keep only the Broadwell collisions (type a.) and the collisions of type b., i.e., we set *γ* = *λ* = 0. On the other side, the numerical simulations show that the addition of the collisions from the group c. or d. enhances the stability properties.

Finally, we emphasize that, for the model (29)–(37), all the collisions conserve energy (elastic). Generally speaking, this is not necessary because we are focused on the correct reproduction of the mass and momentum equations. For instance, it is possible to construct the model of DV Boltzmann type in one spatial dimension with inelastic collisions [26] (quasi-chemical model with three discrete velocities) which leads to the correct Navier– Stokes equation at small Mach limit.
