François Dubois

Tony Février
email: tony.février@math.u-psud.fr

Benjamin Graille
email: benjamin.graille@math.u-psud.fr

STABILITY OF A BIDIMENSIONAL RELATIVE VELOCITY LATTICE BOLTZMANN SCHEME

In this contribution, we study the theoretical and numerical stability of a bidimensional relative velocity lattice Boltzmann scheme. These relative velocity schemes introduce a velocity field parameter called "relative velocity" function of space and time. They generalize the d'Humières multiple relaxation times scheme and the cascaded automaton. This contribution studies the stability of a four velocities scheme applied to a single linear advection equation according to the value of this relative velocity. We especially compare when it is equal to 0 (multiple relaxation times scheme) or to the advection velocity ("cascaded like" scheme). The comparison is made in terms of L ∞ and L 2 stability. The L ∞ stability area is fully described in terms of relaxation parameters and advection velocity for the two choices of relative velocity. These results establish that no hierarchy of these two choices exists for the L ∞ notion. Instead, choosing the parameter equal to the advection velocity improves the numerical L 2 stability of the scheme. This choice cancels some dispersive terms and improve the numerical stability on a representative test case. We theoretically strengthen these results with a weighted L 2 notion of stability.

Introduction

Lattice Boltzmann schemes are applicable to many different fields such as hydrodynamics, acoustics, magnetohydrodynamics and multiphase fluids [START_REF] Humières | Generalized Lattice-Boltzmann equations[END_REF][START_REF] Lallemand | Theory of the Lattice Boltzmann method: Dispersion, dissipation, isotropy, Galilean invariance and stability[END_REF][START_REF] Marié | Comparison between Lattice Boltzmann method and Navier-Stokes high order schemes for computational aeroacoustics[END_REF][START_REF] Dellar | Lattice kinetic schemes for magnetohydrodynamics[END_REF][START_REF] He | A lattice Boltzmann scheme for incompressible multiphase flow and its application of Rayleigh-Taylor instability[END_REF], because the associated algorithm is simple, fast and flexible. This algorithm, associated with a cartesian lattice and a finite set of velocities, consists in computing some particle distribution functions at discrete values of time: they are updated through two successive phases of transport (exact) and relaxation (also called collision).

In the d'Humières multiple relaxation times (MRT) framework [START_REF] Humières | Generalized Lattice-Boltzmann equations[END_REF], the relaxation is made thanks to a set of moments, linear combinations of the particle distributions. Each of these moments relaxes towards an equilibrium with a proper time scale. The MRT schemes have been widely studied in terms of consistency [START_REF] Dubois | Equivalent partial differential equations of a lattice Boltzmann scheme[END_REF][START_REF] Dubois | Third order equivalent equation of Lattice Boltzmann scheme[END_REF][START_REF] Junk | Asymptotic analysis of the lattice Boltzmann equation[END_REF] but slightly in terms of stability [START_REF] Lallemand | Theory of the Lattice Boltzmann method: Dispersion, dissipation, isotropy, Galilean invariance and stability[END_REF][START_REF] Ginzburg | Optimal Stability of Advection-Diffusion Lattice Boltzmann Models with Two Relaxation Times for Positive/Negative Equilibrium[END_REF]. They can be particularly used to solve the Navier-Stokes equations [START_REF] Humières | Lattice BGK models for Navier-Stokes equation[END_REF][START_REF] Chen | Lattice Boltzmann Method for fluid flows[END_REF][START_REF] Humières | Generalized Lattice-Boltzmann equations[END_REF][START_REF] Lallemand | Theory of the Lattice Boltzmann method: Dispersion, dissipation, isotropy, Galilean invariance and stability[END_REF][START_REF] Geier | Cascaded digital Lattice Boltzmann automata for high Reynolds number flow[END_REF]. They however undergo some instabilities in the small viscosity limit [START_REF] Lallemand | Theory of the Lattice Boltzmann method: Dispersion, dissipation, isotropy, Galilean invariance and stability[END_REF].

A "cascaded automaton" [START_REF] Geier | Ab initio derivation of the cascaded Lattice Boltzmann Automaton[END_REF], introduced in two and three dimensions for the Navier-Stokes equations, reduces these instabilities [START_REF] Geier | Cascaded digital Lattice Boltzmann automata for high Reynolds number flow[END_REF]. Some studies aim to understand Date: June 8, 2015. the algorithm of this automaton and its improvements in terms of stability. First, it has been written in the MRT framework with a generalized equilibrium depending on the non conserved quantities [START_REF] Asinari | Generalized local equilibrium in the Cascaded Lattice Boltzmann method[END_REF]. It has then been included in a new class of relative velocity schemes [START_REF] Dubois | Lattice Boltzmann schemes with relative velocities[END_REF][START_REF] Février | Extension et analyse des schémas de Boltzmann sur réseau : les schémas à vitesse relative[END_REF] that are the center of interest of the authors. Both MRT and cascaded are particular cases of the relative velocity schemes.

These relative velocity schemes use a set of moments depending on a velocity field parameter called relative velocity for the relaxation. This parameter is a function of space and time. Their consistency has been studied for one and two conservation laws [START_REF] Dubois | Lattice Boltzmann schemes with relative velocities[END_REF][START_REF] Février | Extension et analyse des schémas de Boltzmann sur réseau : les schémas à vitesse relative[END_REF]. Their stability has also been investigated numerically in the case of the D 2 Q 9 scheme for the Navier-Stokes equations [START_REF] Dubois | On the stability of a relative velocity lattice Boltzmann scheme for compressible Navier-Stokes equations[END_REF]. This study, using the L 2 von Neumann stability [START_REF] Lallemand | Theory of the Lattice Boltzmann method: Dispersion, dissipation, isotropy, Galilean invariance and stability[END_REF][START_REF] Sterling | Stability Analysis of Lattice Boltzmann methods[END_REF], shows the importance of the choice of the polynomials defining the moments, that justifies the stability gain obtained by the cascaded automaton.

This contribution studies the theoretical and the numerical stability of a bidimensional scheme with four velocities for the simulation of a single linear advection equation. We expect to improve the stability for a relative velocity equal to the advection velocity ("cascaded like" scheme) instead of 0 (MRT scheme). This study is based on two stability notions. The first is a L ∞ notion that has been used for one dimensional and vectorial schemes [START_REF] Dellacherie | Construction and Analysis of Lattice Boltzmann Methods Applied to 1D Convection-Diffusion Equation[END_REF][START_REF] Graille | Approximation of mono-dimensional hyperbolic systems : A lattice Boltzmann scheme as a relaxation method[END_REF]. The second is a weighted L 2 notion introduced in [START_REF] Banda | A Stability Notion for Lattice Boltzmann Equations[END_REF]. It has been applied to MRT schemes linearized around the zero velocity or in the BGK (Bhatnagar, Gross, Krook) case [START_REF] Bhatnagar | The Lattice Boltzmann equation: theory and applications[END_REF] for advection equations [START_REF] Banda | A Stability Notion for Lattice Boltzmann Equations[END_REF][START_REF] Junk | Weighted L 2 -stability of the lattice Boltzmann method[END_REF][START_REF] Rheinländer | On the stability structure for lattice Boltzmann schemes[END_REF].

In the first part, we briefly recall the framework of the relative velocity schemes. Then we present the four velocities scheme of interest. In the second part, the effect of the relative velocity on the stability is illustrated thanks to the method of the equivalent equations [START_REF] Warming | The modified equation approach to the stability and accuracy analysis of finite difference methods[END_REF][START_REF] Dubois | Equivalent partial differential equations of a lattice Boltzmann scheme[END_REF]. We particularly focus on the structure of the dispersion terms for two choices of relative velocity (0 or the advection velocity). These results are then illustrated by a numerical L 2 von Neumann stability study. The two last parts give theoretical stability results. In the third part, the L ∞ stability area is fully described in terms of relaxation parameters and advection velocity for the two choices of relative velocity. In the fourth one, a weighted L 2 notion [START_REF] Junk | Weighted L 2 -stability of the lattice Boltzmann method[END_REF] is used to validate the numerical results obtained in the second part according to the relative velocity.

Framework

This section presents first the relative velocity schemes in a general framework. It is then particularized to a bidimensional scheme with four velocities for the simulation of a linear advection equation.

Presentation of the relative velocity schemes

For d ∈ N , we consider a cartesian lattice L of R d associated with a space step ∆x. The space time is given by the acoustic scaling ∆t = ∆x/λ for λ ∈ R a velocity scale. A set of q ∈ N velocities denoted by v = {v 0 , . . . , v q-1 } depending on the velocity scale is chosen such that for each node x ∈ L, x + v j ∆t is still a node of L. The lattice Boltzmann method computes the discrete values of several particle distributions denoted by f j , 0 j q -1, the distribution f j being associated with the velocity v j . We denote by f = (f 0 , . . . , f q-1 ) the vector of the particle distributions. An iteration of the algorithm consists in the succession of a phase of relaxation, nonlinear and local in space, and a phase of transport distributing the particles on the neighbouring nodes.

The relaxation operator of the relative velocity schemes originates from the framework of the MRT schemes [START_REF] Humières | Generalized Lattice-Boltzmann equations[END_REF]. Some moments, linear combinations of the particle distributions, relax with a proper relaxation rate. We define a matrix of moments depending on a velocity field u function of space and time. This matrix, supposed to be invertible, is characterized by a set of polynomials P k , 0 k q-1, M ( u) kj = P k (v j -u), 0 k, j q-1.

Let us note that the MRT scheme corresponds to u = 0. The vector of moments is then obtained from the particle distributions through a linear transformation

m( u) = M ( u) f . (1)
The relaxation then carries on the components m k ( u), 0 k q -1, of m( u) = (m 0 ( u), . . . , m q-1 ( u))

m k ( u) = m k ( u)+s k (m eq k ( u)-m k ( u)), 0 k q-1, (2)
where m eq k ( u) are the moments at equilibrium and s k the relaxation parameters. Some relaxation parameters are chosen null to enforce some conservation laws. The equilibrium is supposed to derive from an equilibrium distribution f eq ∈ R q independent of u m eq ( u) = M ( u)f eq .

The equilibrium f eq is an a priori non linear function of these conserved moments. The postcollision distributions are obtained by the linear transformation

f = M -1 ( u)m ( u). (4)
Then, the particles are advected from node to node f j (x, t + ∆t) = f j (x -v j ∆t, t), 0 j q-1.

In the following, the relative velocity u is a real constant. When u is specified, the relative velocity scheme is called scheme relative to u.

1.2. The twisted relative velocity D 2 Q 4 scheme for a linear advection equation

In this section, we present a relative velocity D 2 Q 4 scheme for a linear advection equation. We expect the scheme relative to the advection velocity ( u = V ), also called "cascaded like" scheme in the following, to be more stable than the MRT scheme u = 0 when this velocity V increases.

We want to approximate the following bidimensional advection equation

∂ t ρ(x, t) + V • ∇ρ(x, t) = 0, x ∈ R 2 , t ∈ R + , (5)
where V = (V x , V y ) ∈ R 2 is the advection velocity. With this aim, we consider the four velocity scheme called twisted D 2 Q 4 defined by the set of velocities v = {(λ, λ), (-λ, λ), (-λ, -λ), (λ, -λ)}, for λ = ∆x/∆t the velocity scale (figure 1). The moments are associated with the polynomials 1, X, Y, XY . Another choice with four velocities could be the D 2 Q 4 defined by the set v = {(λ, 0), (0, λ), (-λ, 0), (0, -λ)}, and by the moments associated with 1, X, Y, X 2 -Y 2 . However, we have chosen to work with the twisted scheme because this scheme is exactly a D 1 Q 2 2 in terms of velocities and moments. This feature simplifies the proofs of consistency and stability. Moreover, the stability results of the D 2 Q 4 scheme can be deduced of the results for the twisted D 2 Q 4 scheme: the stability areas correspond thanks to the composition of a rotation and a homothety (Appendix A).

The framework containing only one equation, we choose one conservation law on the density

ρ = 3 j=0 f j .
The other moments are relaxed with a relaxation parameter s q ∈ R for the first order moments and s xy ∈ R for the second order one. To approach the equation ( 5), the equilibrium must read

m eq ( u) = (ρ, (V x -u x )ρ, (V y -u y )ρ, m eq 3 ( u)), (6)
where m eq 3 ( u) is the second order moment equilibrium determining the diffusion terms. Two different equilibria are chosen: the non intrinsic equilibrium

m eq 3 ( u) = ρ(V x -u x )(V y -u y ), (7)
and the intrinsic equilibrium

m eq 3 ( u) = ρ( u x u y -u x V y -u y V x ). (8)
The first equilibrium results from the will to see the twisted D 2 Q 4 as the product of the D 1 Q 2 by himself: it could be seen as the product of two one dimensional equilibria. The second equilibrium is called intrinsic because it leads to a diffusion term that is independent of the basis of writing (proposition 2.2) but is not isotropic.

It is important to note that these equilibria are chosen such that the relation ( 3) is satisfied. Then the second order asymptotics do not depend on the relative velocity u [START_REF] Dubois | Lattice Boltzmann schemes with relative velocities[END_REF][START_REF] Février | Extension et analyse des schémas de Boltzmann sur réseau : les schémas à vitesse relative[END_REF].

Exploitation of the equivalent equations for the stability

The purpose of this section is to predict the stability behaviour of the twisted scheme defined in the section 1 according to the choice of u. To do so, we present the third order asymptotics of the scheme: we discuss on the definite positivity of the diffusion tensor and on the dispersive third order terms. The discussion is then illustrated on a numerical study of L 2 stability. We show that the cascaded like scheme contrary to the MRT scheme cancels some dispersive terms and stabilizes the scheme in the same time.

Discussion for the scheme with a non intrinsic diffusion

We study here the third order equivalent equations of the twisted D 2 Q 4 with a non intrinsic diffusion. The following proposition results from the particularization of a general derivation of the equivalent equations for a D d Q q scheme with one conservation law, d, q ∈ N [START_REF] Février | Extension et analyse des schémas de Boltzmann sur réseau : les schémas à vitesse relative[END_REF]. We are interested in the structure of the upper orders according to the velocity field u and their influence on the stability.

Proposition 2.1 (Diffusion and dispersion operators). Given V ∈ R 2 , (s q , s xy ) ∈ R 2 , the twisted D 2 Q 4 scheme relative to a constant velocity u ∈ R 2 associated with the equilibrium (6,7) and with the relaxation parameters (0, s q , s q , s xy ), approaches the third order equation

∂ t ρ(x, t) + V • ∇ρ(x, t) -∆t D 2 : ∇ 2 ρ(x, t) + ∆t 2 D 3 : ∇ 3 ρ(x, t) = O(∆t 3 ), x ∈ R 2 , t ∈ R + , (9)
where the diffusion matrix D 2 is given by

D 2 = σ q (λ 2 -(V x ) 2 ) 0 0 (λ 2 -(V y ) 2 ) , (10)
and the dispersion matrix D 3 by

D 3 = V x 6 (λ 2 -(V x ) 2 )(1 -12σ 2 q ) 2σ q (λ 2 -(V y ) 2 )( u x -V x )(σ q -σ xy ) 2σ q (λ 2 -(V x ) 2 )( u y -V y )(σ q -σ xy ) V y 6 (λ 2 -(V y ) 2 )(1 -12σ 2 q ) . (11)
The Hénon's parameters reas σ q = 1/s q -1/2 and σ xy = 1/s xy -1/2, the second and third order derivative operators are given by

∇ 2 ρ = ∂ 2 x ρ ∂ 2 xy ρ ∂ 2 xy ρ ∂ 2 y ρ , ∇ 3 ρ = ∂ 3 x ρ ∂ 3 xyy ρ ∂ 3 xxy ρ ∂ 3 y ρ , ( 12
)
and : is the scalar product for the matrices viewed as vectors of R 4 .

These equations are obtained thanks to a formal calculus software, that implements an algorithm to obtain the equivalent equations in the linear case [START_REF] Augier | On rotational invariance of lattice Boltzmann schemes[END_REF]. The following lemma exhibits the area of well-posedness of the second order truncation of (9).

Lemma 2.1 (Well-posedness of the second order equation). The diffusion matrix D defined by [START_REF] Dubois | Equivalent partial differential equations of a lattice Boltzmann scheme[END_REF] belongs to the space S ++ 2 (R) of symmetric positive definite matrices if and only if |V | ∞ < λ.

Proof. The spectrum of the matrix

D 2 is {λ 2 -(V x ) 2 , λ 2 -(V y ) 2 } that closes the proof.
Because of the relation (3), the velocity parameter u appears only at the third order of the equations. Let us discuss of the structure of these third order dispersive terms according to the choice of u.

First we choose the velocity field u equal to the advection velocity V ("cascaded like" scheme). The dispersion terms are then small compared to the diffusion terms. Indeed, this choice cancels the two off-diagonal components of D 3 depending on the parameter σ xy . Thus resting dispersion terms are

V x (λ 2 -(V x ) 2 )(1 -12σ 2
q )/6 and its symmetric in y. When σ q tends to 0, those terms are equivalent to V (λ 2 -(V x ) 2 ) multiplied by a constant. Thus when the diffusion (λ 2 -(V x ) 2 ) decreases, the dispersion decreases at least with the same speed.

On the contrary, the dispersion terms of the MRT scheme ( u = 0) create instabilities when the diffusion (represented by σ q ) is weak. Indeed, the off-diagonal terms are conserved for this choice of velocity parameter: they are given by -2σ q V y (λ 2 -(V x ) 2 )(σ q -σ xy ) and its symmetric in y. This term is equivalent to 2σ xy σ q V y (λ 2 -(V x ) 2 ), when σ q gets close to 0. The dispersion terms are thus depending on the size of the parameter σ xy . If σ xy is important and σ q tends to 0, the numerical stability should be deteriorated by dispersion phenomena for V x or V y close to λ: indeed, the dispersion, behaving as 2σ xy σ q λ(λ 2 -(V x ) 2 ), becomes greater than the diffusion given by σ q (λ 2 -(V x ) 2 ). Instead, taking σ xy close to 0 should limit the dispersion effects.

Remark 2.1 (Exact third order scheme). We note that the choices u = V and σ q = 1/ √ 12 cancel the dispersion matrix D 3 given by [START_REF] Dubois | Third order equivalent equation of Lattice Boltzmann scheme[END_REF]. Thus the scheme is consistent with the second order advection diffusion truncation of (9) at the third order. Moreover, there is still a degree of freedom with the parameter σ xy . This is impossible for the MRT scheme ( u = 0) excepted when there is only one relaxation parameter (BGK scheme: σ q = σ xy ).

Discussion for the scheme with an intrinsic diffusion

We discuss on the structure of the third order equivalent equations for the twisted D 2 Q 4 intrinsic scheme obtained thanks to the linear algorithm derived in [START_REF] Augier | On rotational invariance of lattice Boltzmann schemes[END_REF]. The study is completely analogous to the one of the section 2.1.

Proposition 2.2 (Diffusion and dispersion operators

). Given V ∈ R 2 , (s q , s xy ) ∈ R 2
, the twisted D 2 Q 4 scheme relative to a constant velocity u associated with the equilibrium [START_REF] Dellacherie | Construction and Analysis of Lattice Boltzmann Methods Applied to 1D Convection-Diffusion Equation[END_REF][START_REF] Humières | Generalized Lattice-Boltzmann equations[END_REF] and with the relaxation parameters (0, s q , s q , s xy ), approaches the following third order equation

∂ t ρ(x, t) + V • ∇ρ(x, t) -∆t D 2 : ∇ 2 ρ(x, t) + ∆t 2 D 3 : ∇ 3 ρ(x, t) = O(∆t 3 ), x ∈ R 2 , t ∈ R + ,
where the diffusion matrix D 2 is given by

D 2 = σ q (λ 2 -(V x ) 2 ) -σ q V x V y -σ q V x V y σ q (λ 2 -(V y ) 2 ) , (13)
and the dispersion matrix D 3 by

D 3 = φ 1 (V x ) φ 2 ( u x , u y , V x , V y ) φ 2 ( u y , u x , V y , V x ) φ 1 (V y ) , with φ 1 (V x ) = V x 6 (λ 2 -(V x ) 2 )(1 -12σ 2 q ), φ 2 ( u x , u y , V x , V y ) = 1 2 -V x (V y ) 2 + 4σ 2 q (λ 2 ( u x -V x ) + 3V x (V y ) 2 -u y V x V y -u x (V y ) 2 ) -4σ q σ xy (λ 2 ( u x -V x ) -(V x u y V y + u x (V y ) 2 ).
The Hénon's parameters σ q and σ xy are given by 1/s q -1/2 and 1/s xy -1/2, the operators ∇ 2 and ∇ 3 are given by [START_REF] Dubois | Lattice Boltzmann schemes with relative velocities[END_REF] and : is the scalar product for the matrices viewed as vectors of R 4 .

Remark 2.2. As expected, the diffusion is intrinsic. It is independent of the basis of writing. Indeed it is equal to

D 2 : ∇ 2 = λ 2 (∂ 2 x + ∂ 2 y ) -(V • ∇) 2 .
Lemma 2.2 (Well-posedness of the second order equation). The diffusion matrix D 2 defined by [START_REF] Dubois | On the stability of a relative velocity lattice Boltzmann scheme for compressible Navier-Stokes equations[END_REF] belongs to the space S ++ 2 (R) of symmetric positive definite matrices if and only if |V | 2 < λ.

Proof. The eigenvalues of the matrix D 2 are λ 2 and λ 2 -|V | 2 2 , that closes the proof. Let us compare the well-posedness areas of the schemes with a non intrinsic diffusion (section 2.1) and an intrinsic diffusion. These areas are represented on the figure 2. We note that the scheme with a non intrinsic diffusion has a greater area in V . Consequently, it should authorize larger stable velocities than the scheme with an intrinsic diffusion.

0 λ -λ λ -λ V x V y
We now get back to the intrinsic case and discuss on the dispersion terms. As in the non intrinsic case, the scheme should be stable for larger velocities with the choice u = V than with u = 0. Indeed, when σ q tends to 0, the term φ 2 ( u x , u y , V x , V y ) is of third order in V if u = V . Instead, if u = 0, it is of first order in V . More dispersion is created for u = 0 when V x or V y are close to λ and σ xy is big.

Illustration with numerical stability experiments

In this section, we study the numerical L 2 stability for the twisted D 2 Q 4 scheme according to the choice of the velocity parameter u. The link with the previous consistency study is made: the third order dispersion terms have an influence on the numerical stability. The choice u = V provides larger stability areas in V . We present a first numerical experiment with the advection of a circular spot: the initial conditions are given by

ρ(x, t) = 1 + 1 C (x, t), x ∈ [0, 1] 2 , t ∈ R,
where C is the disc centered in (1/2, 1/2) of radius 0.1 and 1 C is the function equal to 1 on this disc, null elsewhere. This spot moves with an advection velocity V ∈ R 2 in the domain [0, 1] 2 constituted of 128 2 points with periodic boundary conditions. The scheme is considered as stable if it has not broken after 2000 iterations. We present the biggest advection velocities keeping the scheme stable for different choices of relaxation parameters. Two different directions of advection velocities are considered: θ = 0 and θ = π/4. We choose different values for the parameter σ q , the parameter σ xy being set to 1/ √ 3. The results are presented in the tables 1 and 2.

The scheme relative to the advection velocity u = V authorizes larger stable velocities than the MRT scheme ( u = 0). Most of the time, its stability areas do not depend on the relaxation parameters while they are bounded by 0 and 2. Instead, the stability areas of the MRT scheme decreases when σ q tends to 0. This phenomenon occurs for the scheme relative to u = V with an intrinsic diffusion in the direction θ = π/4: however the instabilities appear for larger V than the MRT scheme.

The dispersion terms exhibited in the equivalent equations (section 2) originate these phenomena. Their presence causes instabilities when σ q is too small compared to σ xy and |V |. We have seen that these terms cancel for u = V but are present when u = 0. This explains the constance of the stability area for u = V and its deterioration for u = 0 when σ q decreases.

These dispersive phenomena are illustrated on the figures 3 and 4 for the scheme with an intrinsic diffusion. We choose σ xy = 1/ √ 3, σ q = 1/20 and V = (0.9λ, 0) for these draws. According to the data of the table 1, this configuration is a stable case for the scheme relative to u = V and an unstable case for the MRT scheme. The spot is represented at the times t = 0 and t = 0.4 for 256 2 points. The dispersion for the MRT scheme is clearly visible on the figure 3: some oscillations appear behind the advected spot. The scheme is going to break since the density ρ is increasing. These phenomena are absent for the scheme relative to u = V (figure 4) that remains numerically stable.

Our second numerical experiment uses the L 2 linear stability of von Neumann. We discuss on the spectrum of the amplification matrix of the scheme. This matrix, characterizing an iteration of the scheme, has to be first determined. The equilibrium being linear [START_REF] Dellacherie | Construction and Analysis of Lattice Boltzmann Methods Applied to 1D Convection-Diffusion Equation[END_REF][START_REF] Dellar | Lattice kinetic schemes for magnetohydrodynamics[END_REF][START_REF] Humières | Generalized Lattice-Boltzmann equations[END_REF], there exists a matrix E = E(V , s) for s = (0, s q , s q , s xy ) and V ∈ R 2 so that f eq = Ef .

The relaxation phase of the relative velocity D 2 Q 4 scheme reads

f = (I + M ( u) -1 DM ( u)(E -I))f ,
where D = diag(s) is the diagonal matrix of the relaxation parameters. This expression holds for each node x of the lattice, the relaxation being local in space. The . Velocities V stable in λ scale for the scheme relative to u = 0 (MRT on the left), u = V (on the right), for a non intrinsic diffusion with s q = 1 and s xy = 1.

expression of the distribution after the transport phase is given by

f j (x, t + ∆t) = [(I + M ( u) -1 DM ( u)(E -I))f ] j (x -v j ∆t, t), x ∈ L, t ∈ R. (14)
Taking the Fourier transform of ( 14), the transport operator becomes local in space and is represented by the diagonal matrix A = A(k) for k ∈ R 2 whose diagonal components are given by e i∆tk.v j , 0 j 3. The amplification matrix then reads

L( u) = L( u, V , k, s) = A(I + M ( u) -1 DM ( u)(E -I))
. It characterizes a time iteration of the scheme in the Fourier space

f (k, t + ∆t) = L( u) f (k, t), t ∈ R,
where f is the Fourier transform of f . We want to exhibit the advection velocities V for which the scheme verifies the necessary condition of L 2 stability max

k∈R 2 r(L( u)) 1, ( 15
)
where r is the spectral radius. We are thus interested in the set

{V = (V x , V y ), max k∈R 2 r(L( u)) 1}.
The figures 5 to 10 present this set in the plan (V x , V y ) for the twisted D 2 Q 4 scheme with a non intrinsic diffusion and some relaxation parameters s q and s xy . The left draw is about the MRT scheme ( u = 0) and the right one illustrates the case u = V . We expect the right areas to be larger than the left ones.

Firstly, the draws corresponding to a BGK scheme (s q = s xy = 1) are independent of u (figure 5). This result was expected because the velocity field u does not appear in the scheme for one single relaxation parameter since (3) is verified. Secondly, the scheme relative to u = V verifies the necessary condition of stability on the biggest area |V | ∞ λ (CFL) for all the s chosen (figure 5 to 10). The constance and the optimality of these areas confirm the phenomena observed for the advected spot . Velocities V stable in λ scale for the scheme relative to u = 0 (MRT on the left), u = V (on the right), for a non intrinsic diffusion with s q = 1 and s xy = 1.9.

(tables 1 and 2). Finally, the stability areas of the MRT scheme decrease as the two relaxation parameters are moving away from each other. The stability area is reduced to V = 0 for s q = 2 or s xy = 2. Whatever the choice of s, the areas of the MRT scheme are included in those of the scheme relative to u = V . These results are also consistent with the test case of the advected spot: the deterioration of the areas are due to the third order dispersive terms of the equivalent equations (proposition 2.1).

For the D 2 Q 4 scheme with an intrinsic diffusion (figures 11 to 16), the results are different but the general trend is the same. The scheme relative to u = V verifies [START_REF] Geier | Ab initio derivation of the cascaded Lattice Boltzmann Automaton[END_REF] for larger sets of V than the MRT scheme, excepted when s q = 1, s xy = 1.5 (figure 12). However, the areas associated to u = V are not optimal any more: they decrease when s q and s xy move away from each other. Both velocity fields ( u = 0 and . Velocities V stable in λ scale for the scheme relative to u = 0 (MRT on the left), u = V (on the right), for a non intrinsic diffusion with s q = 1.9 and s xy = 1. u = V ) undergo this phenomenon but the effect of the dispersion is bigger for the MRT scheme. These results confirm the advected spot ones: the same deterioration appears when σ q tends to 0 (table 1 and2) to a lesser extent for u = V than for u = 0.

The main conclusions of this section are the following. Choosing u = V cancels some dispersive terms becoming important as V grows. This cancellation does not occur for u = 0. We must view these results in parallel with the fact that the scheme is more stable (numerical L 2 notion) for u = V than for u = 0. . Velocities V stable in λ scale for the scheme relative to u = 0 (MRT on the left), u = V (on the right), for an intrinsic diffusion with s q = 1 and s xy = 1.

Theoretical L ∞ stability

We study in this section the stability of the relative velocity D 2 Q 4 scheme with respect to a notion of L ∞ stability presented in [START_REF] Dellacherie | Construction and Analysis of Lattice Boltzmann Methods Applied to 1D Convection-Diffusion Equation[END_REF][START_REF] Graille | Approximation of mono-dimensional hyperbolic systems : A lattice Boltzmann scheme as a relaxation method[END_REF]. We use it to demonstrate some maximum principles for u equal to 0 (MRT scheme) and V , the advection velocity ("cascaded like" scheme). The L ∞ stability area is fully described in terms of relaxation parameters and advection velocity for these two choices of relative velocity. We show that this L ∞ notion does not differentiate u = 0 from u = V in terms of stable behaviour. In each case, the parameters (s q , s xy ) corresponding to a non empty area of L ∞ stability in V are included in the square [0, 2] 2 . . Velocities V stable in λ scale for the scheme relative to u = 0 (MRT on the left), u = V (on the right), for an intrinsic diffusion with s q = 1 and s xy = 1.9. Definition 3.1 (L ∞ stability). Let us consider a D d Q q scheme on the cartesian lattice L ∈ R d . The mass at time t is given by

ρ tot (t) = x∈L q-1 j=0 f j (x, t), t ∈ R.
Let C ∈ R + , suppose that the initial distributions of particles are nonnegative and that the initial mass is bounded by C: f j (x, 0) 0, 0 j q -1, x ∈ L, and ρ tot (0) C, the scheme is said to be L ∞ stable if for all time t n ∈ R, n ∈ N, f j (x, t n ) 0, 0 j q -1, x ∈ L, and ρ tot (t n ) C. . Velocities V stable in λ scale for the scheme relative to u = 0 (MRT on the left), u = V (on the right), for an intrinsic diffusion with s q = 1.9 and s xy = 1.

The proofs of this a priori non linear L ∞ stability notion are only depending on the relaxation phase. Indeed, the transport exchanges the distributions between each other: it does not act on their positivity and on their mass. Thus it is sufficient to show that the postcollision distributions remain nonnegative and the mass bounded. Note also that if ρ tot (0) is bounded, ρ tot (t n ) remains bounded for all time t n because this mass is conserved by the scheme. So there is only to study the positivity of the postcollision distributions functions.

To do so, we express the postcollision distributions of particles f as a function of the precollision one f for a general D d Q q scheme with one conservation law. These postcollision distributions are first expressed as functions of the postcollision moments thanks to (4). After applying the relaxation identity (2), the precollision moments are written in the frame of the particle distributions with [START_REF] Asinari | Generalized local equilibrium in the Cascaded Lattice Boltzmann method[END_REF]. These calculations lead to the following expression of the relaxation.

f l = 1 - q-1 k=1 s k M ( u) -1 lk M ( u) kl f l - q-1 j=0 j =l q-1 k=1 s k M ( u) -1 lk M ( u) kj f j + q-1 k=1 s k M ( u) -1 lk m eq k ( u), 0 l q -1. ( 16
)
The framework of study is linear because we are approaching a linear advection equation. As a consequence, the equilibrium term of ( 16) is a linear function of the conserved variable ρ and thus of the distributions f j . Consequently, there is a matricial relation between f and f : sufficient conditions of L ∞ stability are obtained checking when this matrix is nonnegative.

Definition 3.2. For q ∈ N , a square matrix of M q (R) is said to be nonnegative if all its coefficients are nonnegative.

Note that this notion is different from a positive matrix that have all its eigenvalues nonnegative [START_REF] Lascaux | Analyse numérique matricielle appliquée à l'art de l'ingénieur (tome 1 : méthodes directes)[END_REF]. In the following, we prove the results in the twisted case: the analogous ones for the D 2 Q 4 scheme are available in the appendix B and use the proposition A.1 of the appendix A. We first consider the relative velocity D 2 Q 4 scheme with a non intrinsic diffusion (equilibrium [START_REF] Dellar | Lattice kinetic schemes for magnetohydrodynamics[END_REF]). The proofs are quite technical and involves geometric inequalities on V and s.

The non intrinsic diffusion case

Proposition 3.1 (L ∞ stability areas for the MRT scheme). Let V ∈ R 2 , (s q , s xy ) ∈ R 2 , consider the twisted D 2 Q 4 MRT scheme ( u = 0) with the relaxation parameters (0, s q , s q , s xy ), associated with the equilibrium (6,7)

m eq (0) = ρ(1, V x , V y , V x V y ). Note V = (V x , V y ) where V x = V x + V y , V y = V x -V y and γ = s q /s xy :
• if 0 < s xy min(s q , 2 -s q ) (area BCD on the figure 17), the scheme is L ∞ stable for all V so that

(V x ± 2λγ) 2 -(V y ) 2 4λ 2 (γ 2 -1), (17a)
(V y ± 2λγ) 2 -(V x ) 2 4λ 2 (γ 2 -1). (17b)
• if s q s xy 2s q and s q 1 (area ABC on the figure 17), the scheme is L ∞ stable for all V so that

(V x ± 2λγ) 2 -(V y ) 2 4λ 2 (γ -1) 2 , (18a)
(V y ± 2λγ) 2 -(V x ) 2 4λ 2 (γ -1) 2 . (18b)
• if 2 -s q s xy 2(2 -s q ) and s q 1 (area ACD on the figure 17), the scheme is L ∞ stable for all V so that

(V x ± 2λγ) 2 -(V y ) 2 4λ 2 (γ + 1) 2 - 4 s xy , (19a)
(V y ± 2λγ) 2 -(V x ) 2 4λ 2 (γ + 1) 2 - 4 s xy . (19b)
• if s xy = 0 and 0 < s q 2 (ray ]BD] on the figure 17), the scheme is L ∞ stable for V = 0. • if s xy = s q = 0 (point B on the figure 17), the scheme is unconditionally L ∞ stable. • For all other s, there is no V corresponding to a L ∞ stable scheme.

In particular, the parameters (s q , s xy ) corresponding to a non empty area of

L ∞ sta- bility in V are included in the square [0, 2] 2 .
The stability areas in the plan (s q , s xy ) are represented on the figure 17.The limit cases for the parameters s q and s xy lead to consistent inequalities for the velocity conditions.

Proof. The first step consists to express f as a function of f thanks to the relation [START_REF] Geier | Cascaded digital Lattice Boltzmann automata for high Reynolds number flow[END_REF]. If the matrix sending f on f is nonnegative, the scheme is L ∞ stable. A formal calculus software guarantees that it corresponds to ensure the following relations: If s xy = 0, these relations imply V = 0 and 0 s q 2. Otherwise, we factorise by s xy /4λ 2 and express these identities as functions of the coefficient γ = s q /s xy . This leads to

min s xy 4 , 2s q -s xy 4 , 1- 2s q + s xy 4 ± s q 4λ (V x +(-1) i V y )+(-1) i s xy 4λ 2 V x V y 0, i = 0, 1.
λ 2 min 1, 2γ -1, 4 s xy -2γ -1 ±λγ(V x +(-1) i V y )+(-1) i V x V y 0, i = 0, 1, (20)
if s xy > 0 or

λ 2 max 1, 2γ -1, 4 s xy -2γ -1 ± λγ(V x + (-1) i V y ) + (-1) i V x V y 0, i = 0, 1, if s xy < 0.
These inequalities are expressed thanks to quadratic forms in V . To characterize the associated geometry, we write them in an adapted basis. To do so, we use the variable change

V x V y = (V x ) 2 -(V y ) 2 4 , for V x = V x + V y , V y = V x -V y .
The inequalities (20) then read

4λ 2 min 1, 2γ -1, 4 s xy -2γ -1 ± 4λγV x + (V x ) 2 -(V y ) 2 0,
completed by the ones obtained when the roles of V x and V y are exchanged. We then center them in the adapted frame to obtain if s xy > 0

(γ + 1) 2 -4/s xy , area 3 : (γ -1) 2 . (V x ± 2λγ) 2 -(V y ) 2 4λ 2 max γ 2 -1, (γ -1) 2 , (γ + 1) 2 - 4 s xy , (21)
plus the inequalities obtained exchanging the roles of V x and V y . The case s xy < 0 yields to

(V x ± 2λγ) 2 -(V y ) 2 4λ 2 min γ 2 -1, (γ -1) 2 , (γ + 1) 2 - 4 s xy . ( 22
)
We begin by showing that there is no velocity stable in the case s xy < 0. The minimum of ( 22) is equal to (γ -1) 2 (nonnegative) if s q s xy < 0, to γ 2 -1 (negative) if s xy min(s q , 2 -s q ) and (γ + 1) 2 -4/s xy (nonnegative) if 2 -s q s xy < 0. In the first and third cases, [START_REF] Lallemand | Theory of the Lattice Boltzmann method: Dispersion, dissipation, isotropy, Galilean invariance and stability[END_REF] would have a solution if and only if the left pole of the hyperbole centered in 2λ|γ| is negative. We then should have 2λ|γ| -2λ|γ -1| 0, in the first case and 2λ|γ| -2λ((γ + 1) 2 -4/s xy )

1 2 0,
in the third one. Considering that s xy < 0, the first case is equivalent to s xy 2s q that has no intersection with s q s xy < 0 where (γ -1) 2 is the minimum. The third one is equivalent to s xy 2(2 -s q ) whose intersection with 2 -s q s xy < 0 is empty. It remains the case

(V y ) 2 -(V x ± 2λγ) 2 4λ 2 (1 -γ 2 ),
when s xy min(s q , 2 -s q ). These inequations combined with the ones obtained exchanging the roles of V x and V y have no common intersection. This closes the case s xy < 0. We exhibit now the conditions for which [START_REF] Junk | Weighted L 2 -stability of the lattice Boltzmann method[END_REF] has at least a solution V . The eventual stability area is located between the two hyperboles

(V x ± 2λγ) 2 -(V y ) 2 = 4λ 2 max γ 2 -1, (γ -1) 2 , (γ + 1) 2 - 4 s xy ,
plus the ones obtained exchanging the roles of V x and V y . Let us note that this maximum is always nonnegative because (γ -1) 2 is nonnegative. This area is delimited by the points A, B, C, D on the figure 19. The inequalities [START_REF] Junk | Weighted L 2 -stability of the lattice Boltzmann method[END_REF] have a solution V if and only if the left pole P of the hyperbole centered in 2λγ has a nonnegative abscissa (figure 19). When the maximum is

γ 2 -1, this abscissa is 2λ(γ -|γ 2 -1| 1
2 ), that is nonnegative if s q 0. When it is (γ -1) 2 , the abscissa 2λ(2γ -1) is nonnegative if and only if s xy 2s q . Finally when the maximum is (γ + 1) 2 -4/s xy , the abscissa is

equal to 2λ(γ -((γ + 1) 2 -4/s xy ) 1 2 ). The positivity of this coefficient is equivalent to s xy 2(2 -s q ).
It remains to determine the maximum of (21) in the case where there is at least one velocity V stable: three cases, represented on the figure 18, appear. If s xy s q 2 -s xy , the maximum is γ 2 -1, if s xy 2 -s q and s q 1, it is (γ + 1) 2 -4/s xy , if s q min(s xy , 1), it is (γ -1) 2 .

From the proposition 3.1, we deduce the simpler case of the BGK scheme corresponding to one single relaxation parameter.

Corollary 3.1 (The BGK case). Let V ∈ R 2 , u ∈ R 2 , s ∈ R, consider the twisted D 2 Q 4 scheme relative to u, BGK of parameter s, associated with the equilibrium (6,7)

m eq ( u) = ρ 1, V x -u x , V y -u y , (V x -u x )(V y -u y ) . Note V = (V x , V y ) where V x = V x + V y , V y = V x -V y : • if 0 s 1, the scheme is L ∞ stable for |V | ∞ λ. • if 1 s 4/3, the scheme is L ∞ stable for all V so that (V x ± 2λ) 2 -(V y ) 2 16λ 2 1 - 1 s , (23a)
(V y ± 2λ) 2 -(V x ) 2 16λ 2 1 - 1 s . (23b)
• If s > 4/3, no V corresponds to a L ∞ stable scheme.

In particular, the parameters s corresponding to a non empty area of L ∞ stability in V are included in [0, 2].

Remark 3.1. The BGK scheme does not depend on the velocity field u since the relation ( 3) is verified.

Proof. In this context, the parameter γ is equal to 1. Thus the inequalities [START_REF] Junk | Weighted L 2 -stability of the lattice Boltzmann method[END_REF] summarize in

(V x ± 2λ) 2 -(V y ) 2 max 0, 16λ 2 1 - 1 s ,
plus the inequality obtained exchanging the roles of V x and V y . If s 1, it is equivalent to

(V x ± 2λ) 2 -(V y ) 2 0,
and so to |V | ∞ λ, else we obtain (23a,23b).

The inequalities (23a,23b) admit a non empty set of stable velocities V if the abscissa of the left pole of the hyperbole

(V x -2λ) 2 -(V y ) 2 = 16λ 2 1 - 1 s , is nonnegative. This abscissa, equal to 2λ(1 -2(1 -1/s) 1 2
) for s 1, is nonnegative if and only if s 4/3.

For s 1, the BGK scheme is optimal in terms of L ∞ stability (CFL reached). When s q and s xy go far from each other (TRT scheme), the stability area decreases until it reaches V = 0 at the boundaries of the triangle represented on the figure 17. There is no stable velocities outside of this triangle.

We now focus on the scheme relative to u = V with a non intrinsic diffusion and two relaxation parameters. We show that the optimal area of stability spreads on some TRT schemes contrary to the MRT scheme. Proposition 3.2 (L ∞ stability areas for a relative velocity scheme). Let V ∈ R 2 , (s q , s xy ) ∈ R 2 , consider the twisted D 2 Q 4 scheme relative to u = V associated with the relaxation parameters (0, s q , s q , s xy ), and the equilibrium (6,7) m eq (V ) = ρ(1, 0, 0, 0).

Let us note V = (V x , V y ) where V x = V x +V y , V y = V x -V y and γ = s xy /(2s q -s xy ), then:

• if s q s xy min(1, 2s q ) (area ABC on the figure 20), the scheme is

L ∞ stable on |V | ∞ λ. (24)
• if s xy < 2s q , max(s q , 1) s xy 2(2 -s q ) (area BCED on the figure 20), the scheme is L ∞ stable on the intersection of ( 24) with

(V x ± 2λγ) 2 -(V y ) 2 16λ 2 (2s q -s xy ) 2 (s xy -s q (2 -s q )), (25a)
(V y ± 2λγ) 2 -(V x ) 2 16λ 2 (2s q -s xy ) 2 (s xy -s q (2 -s q )). (25b)
• if 0 s xy min(s q , s q (2 -s q )) (area ACF on the figure 20), the scheme is

L ∞ stable on |V | ∞ λγ. (26)
• if s q (2 -s q ) s xy min(s q , 2(2 -s q ))(area CEF on the figure 20), the scheme is L ∞ stable on (25a,25b). • For all other s, no V corresponds to a L ∞ stable scheme.

• if s xy = 2s q and 0 < s xy 2 (ray ]AD] on the figure 20), the scheme is L ∞ stable on the intersection of [START_REF] Marié | Comparison between Lattice Boltzmann method and Navier-Stokes high order schemes for computational aeroacoustics[END_REF] and

|V | 1 2λ 2 -s xy s xy . (27)
• if s xy = s q = 0 (A on the figure 20), the scheme is unconditionally L ∞ stable.

In particular, the parameters (s q , s xy ) corresponding to a non empty area of L ∞ stability in V are included in the square [0, 2] 2 .

Proof. We determine the sufficient conditions for the positivity of the matrix sending f on f . A formal calculus software leads to the following inequalities

(2s q -s xy )(λ ± V x )(λ ± V y ) 0, (28)
(λ ± V x )(±(2s q -s xy )V y + λs xy ) 0, ( 29
)
4λ 2 -((2s q +s xy )λ 2 ±s xy λ(V x +(-1) i V y )+(-1) i (s xy -2s q )V x V y ) 0, i = 0, 1, (30)
plus the inequality coming from (29) when the roles of V x and V y are exchanged.

We begin by the study of (28). Case 1 : 2s q < s xy . These inequalities are equivalent to

(λ ± V x )(λ ± V y ) 0.
No V verifies this inequality: such a set of V would be the intersection of fourth areas that do not intersect. We must assume that s xy 2s q to eventually have a non empty stability area in V .

Case 2 : s xy = 2s q .

If s xy = 0, all the inequalities are obviously verified and the scheme is unconditionally stable in V . If s xy < 0, λ ± V x 0 that is impossible. If s xy > 0, (29) becomes (λ ± V x )λ 0, plus the analogous inequality in V y : this is equivalent to [START_REF] Marié | Comparison between Lattice Boltzmann method and Navier-Stokes high order schemes for computational aeroacoustics[END_REF]. The inequalities (30) read

1 - s xy 4λ (2λ ± (V x ± V y )) 0,
that is equivalent to [START_REF] Warming | The modified equation approach to the stability and accuracy analysis of finite difference methods[END_REF]. The final stability area is the intersection of ( 24) and [START_REF] Warming | The modified equation approach to the stability and accuracy analysis of finite difference methods[END_REF]. Three cases, represented on the figures 21, 22 and 23, are then possible.

Let

s xy 1, then 2λ (2 -s xy ) s xy 2λ, λ -λ 2λ -2λ λ -λ 2λ -2λ 0 P V x V y Figure 21
. L ∞ stability area for the scheme relative to u = V with a non intrinsic diffusion for s xy = 2s q . The point P is of abscissa 2λ(2 -s xy )/s xy 2λ. and the stability area is given by (24) (figure 21), let 1 s xy 4/3, then

λ -λ 2λ -2λ λ -λ 2λ -2λ 0 P V x V y
λ 2λ (2 -s xy ) s xy 2λ, λ -λ 2λ -2λ λ -λ 2λ -2λ 0 P V x V y Figure 23
. L ∞ stability area for the scheme relative to u = V with a non intrinsic diffusion for s xy = 2s q . The point P is of abscissa 2λ(2 -s xy )/s xy λ. and the stability area is the intersection of ( 24) and ( 27 and the stability area is given by ( 27) (figure 23). This area is reduced to V = 0 when s xy = 2.

Case 3 : s xy < 2s q .

The inequalities (28) read

(λ ± V x )(λ ± V y ) 0,
that is equivalent to [START_REF] Marié | Comparison between Lattice Boltzmann method and Navier-Stokes high order schemes for computational aeroacoustics[END_REF]. The inequalities (29) write

(λ ± V x )(λγ ± V y ) 0, (31)
for γ = s xy /(2s q -s xy ). The inequalities (31) have solutions if γ 0, that is equivalent to s xy 0. If s xy s q , (31) contains [START_REF] Marié | Comparison between Lattice Boltzmann method and Navier-Stokes high order schemes for computational aeroacoustics[END_REF], that correspond to the same constraints that those imposed by (28). If s xy s q , then the stability area reduces to [START_REF] Sterling | Stability Analysis of Lattice Boltzmann methods[END_REF].

There is still to study the influence of (30) on the stability area. These identities read

λ 2 2s q -s xy (4 -(2s q + s xy )) ± γλ(V x + (-1) i V y ) + (-1) i V x V y 0, i = 0, 1, 0 2λ -2λ 2λ -2λ V x V y Figure 24.
Layout of the areas given by ( 26) and (32). after division by 2s q -s xy . We change of frame,

V x = V x + V y , V y = V x -V y , and obtain 4λ 2 2s q -s xy (4 -(2s q + s xy )) ± 4γλV x + (V x ) 2 -(V y ) 2 0,
and the analogous inequality exchanging the roles of V x and V y . We center these equations to obtain (25a) and (25b).

Sub-case 1 : s xy s q (2 -s q ).

The inequality (25a) reads

(V y ) 2 -(V x ± 2λγ) 2 16λ 2 (2s q -s xy ) 2 (s q (2 -s q ) -s xy ). (32)
We represent the areas corresponding to the inequations ( 26) and (32) on the figure 24. The relation [START_REF] Sterling | Stability Analysis of Lattice Boltzmann methods[END_REF] is equivalent in the new frame to |V| 1 2λγ. It is contained in the area delimited by (32) (figure 24). The final stability area is then given by ( 24) if γ 1 and by ( 26) otherwise.

Sub-case 2 : s xy s q (2 -s q ).

The inequality (25a) has a solution V if the abscissa of the left pole of the hyperbole centered in 2λγ is nonnegative. This abscissa, equal to

2λγ - 4λ (2s q -s xy ) (s xy -s q (2 -s q )) 1 2 , is nonnegative if |2 -s xy | 2|s q -1|. (33)
Assuming that s xy 2, then s xy 2|s q -1| + 2 2s q , that gives an empty L ∞ stability area in V (refer to the beginning of the proof). As a consequence s xy 2. If s q 1, the inequality (33) is equivalent to s xy 2s q , otherwise to

s xy 2(2 -s q ). 2λγ -2λγ 2λγ -2λγ 2λ -2λ 2λ -2λ 0 P V x V y Figure 25.
Layout of the bounds of the L ∞ stability areas of the scheme relative to u = V with a non intrinsic diffusion for γ 1 and s q (2 -s q ) s xy 1.

These conditions carrying on s ensure the existence of a non empty stability area in V .

To obtain the final stability area in V given by the set of inequations from (28) to (30), we compare those associated with (25a) to (24) if γ 1 or (26) if γ 1. The two areas [START_REF] Marié | Comparison between Lattice Boltzmann method and Navier-Stokes high order schemes for computational aeroacoustics[END_REF] and ( 26) read respectively |V| 1 2λ, and

|V| 1 2λγ. ( 34
)
in the frame (0, V x , V y ). If γ 1, the stability area is included in the one given by (34): it is represented by the points A, B, C, D on the figure 26. If γ 1 and s xy 1, the area [START_REF] Marié | Comparison between Lattice Boltzmann method and Navier-Stokes high order schemes for computational aeroacoustics[END_REF] is included in the one given by (25a): the abscissa of the left pole P of the hyperbole centered 2λγ is greater than 2λ (figure 25). If s xy 1, it is not the case any more and the area decreases when s xy increases.

The domain of L ∞ stability in s is identical for u = 0 and u = V . It is given by the triangle of the figures 17 and 20. However, the choice u = V increases the set of the s verifying the optimal stability area |V | ∞ λ: indeed the area (24) spreads to TRT schemes consistent with s q s xy min(1, 2s q ).

We now compare the areas in V associated with u = 0 and u = V for three s that don't verify s q s xy min(1, 2s q ). We want to know if u = V provides a better L ∞ stability behaviour than u = 0. The associated inequalities are presented in

2λγ -2λγ 2λγ -2λγ 0 V x V y P A B D C Figure 26.
Layout of the bounds of the L ∞ stability areas of the scheme relative to u = V with a non intrinsic diffusion for γ 1 and s q (2 -s q ) s xy .

Choice of s (0, 1, 1, 1/2) (0, 1, 1, 3/2) (0, 3/2, 3/2, 3/4) Area for u=0

(V x ±4) 2 -(V y ) 2 12 (V y ±4) 2 -(V x ) 2 12 (V x ±4/3) 2 -(V y ) 2 4/9 (V y ±4/3) 2 -(V x ) 2 4/9 (V x ±4) 2 -(V y ) 2 44/3 (V y ±4) 2 -(V x ) 2 44/3 Area for u=V |V | ∞ 1/3 (V x ± 6) 2 -(V y ) 2 32 (V y ± 6) 2 -(V x ) 2 32 |V | ∞ 1/3 |V | ∞ 1
Table 3. L ∞ stability areas of the schemes with a non intrinsic diffusion for λ = 1.

the table 3. The solutions V are represented on the figures from 27 to 29 in the frame (O, V x , V y ): the areas delimited by the points A, B, C, D are associated with the scheme relative to u = V , and those bounded by E, F, G, H correspond to the MRT scheme ( u = 0).

We can not say that a scheme is better than the other in terms of L ∞ stability: generally the areas just intersect. In the section 2.3, we saw that the velocity field u = V improve the stability associated with the blow up of the scheme that rather corresponds to a L 2 notion than to a L ∞ notion. We confirm this fact theoretically in the section 4.

O A B C D E F G H V x V y Figure 27. L ∞ stability areas for s = (0, 1, 1, 1/2). Scheme relative to u = V (A, B, C, D), to u = 0 (E, F, G, H), with a non intrinsic diffusion. O A B C D E F G H V x V y

The intrinsic diffusion case

As for the non intrinsic case, we first deal with the MRT scheme corresponding to u = 0.

Proposition 3.3 (L ∞ stability areas for the MRT scheme). Let V ∈ R 2 , (s q , s xy ) ∈ R 2 , consider the twisted D 2 Q 4 MRT scheme ( u = 0) associated with the relaxation parameters (0, s q , s q , s xy ), and the equilibrium (6,8)

m eq (0) = ρ(1, V x , V y , 0). Note γ = s xy /s q ,
• if 0 s xy min(s q , 2 -s q ) (area BCD on the figure 30), the scheme is L ∞ stable for all V such that |V | 1 λγ.

(35)

• if s q min(1, s xy ) and s xy 2s q (area ABC on the figure 30), the scheme is • if s q max(1, 2 -s xy ) and s xy 2(2 -s q ) (area ACD on the figure 30), the scheme is L ∞ stable for all V such that

L ∞ stable for all V such that |V | 1 λ(2 -γ). ( 36
)
|V | 1 λ 4 s q -2 -γ . ( 37
)
• if s xy = s q = 0 (point B on the figure 30), the scheme is unconditionally L ∞ stable. • For all other s, no V corresponds to a L ∞ stable scheme.

In particular, the parameters (s q , s xy ) corresponding to a non empty area of L ∞ stability in V are included in the square [0, 2] 2 .

Proof. The positivity of the matrix sending f on f reads min(

s xy 4 , 2s q -s xy 4 , 1 - 2s q + s xy 4 ) ± s q 4λ (V x ± V y ) 0. ( 38
)
If s q = 0, then the inequalities (38) impose s xy = 0 and there is unconditional stability in V . Otherwise these inequalities are equivalent to

|V | 1 min λγ, λ(2 -γ), λ( 4 s q -2 -γ) . ( 39
)
The conditions to have a non empty set of V are to be determined: for that we must study when the minimum in (39) is nonnegative.

γ 0 ⇐⇒ s xy 0, 2 -γ 0 ⇐⇒ s xy 2s q , 4 s q -2 -γ 0 ⇐⇒ s xy 2(2 -s q ).
We now suppose that these conditions are verified and we determine the minimum according to the choice of s: if s xy min(s q , 2 -s q ), the area is given by (35), if s q min(1, s xy ), it is given by (36), if s q max(1, 2 -s xy ), we obtain (37). This closes the proof.

The areas are represented on the figure 30. We can deduce the L ∞ stability areas for the BGK scheme.

Proposition 3.4 (The BGK case). Let V ∈ R 2 , s ∈ R, consider the twisted relative velocity D 2 Q 4 scheme, BGK of relaxation parameters s, and of equilibrium given by (6,8)

• if 0 s 1, the scheme is L ∞ stable for all V such that |V | 1 λ. • if 1 s 4/3, the scheme is L ∞ stable for all V such that |V | 1 λ 4 s -3 .
• For all other s, no V corresponds to a L ∞ stable scheme.

In particular, the parameter s corresponding to a non empty area of L ∞ stability in

V are included in [0, 2].
We obtain a constant area of stability for s 1. When s becomes larger than one (overrelaxation) the area decreases as s increases. For s greater than 4/3, there is no velocity L ∞ stable. Now we do the same job for the TRT scheme with u = V . Proposition 3.5 (L ∞ stability areas for a relative velocity scheme). Let V ∈ R 2 , (s q , s xy ) ∈ R 2 , consider the twisted D 2 Q 4 scheme relative to u = V associated with the relaxation parameters (0, s q , s q , s xy ), and the equilibrium (6,8)

m eq (V ) = ρ(1, 0, 0, -V x V y ).
Note γ 1 = (2s q -s xy )/(s q -s xy ) and γ 2 = s xy /(s q -s xy ),

• if s q < s xy min(2s q , 2(2 -s q )) (area ACD on the figure 31), the scheme is stable for all V such that

(V x ± λγ 1 ) 2 -(V y ) 2 λ 2 γ 1 γ 2 , ( 40a
) (V x ± λ γ 1 + γ 2 2 ) 2 -(V y ± λ) 2 λ 2 γ 2 2 , ( 40b
) (V x ± λ γ 1 + γ 2 2 ) 2 -(V y ∓ λ) 2 λ 2 γ 2 2 , ( 40c
) (V x ± λγ 2 ) 2 -(V y ) 2 λ 2 (s q -s xy ) 2 (4s 2 q -2s q s xy -s 2 xy + 8(s xy -s q )), (40d)
plus the analogous inequalities exchanging V x and V y . • if s xy min(s q , 2(2 -s q )) (area ABD on the figure 31), the scheme is stable for all V such that

(V x ± λγ 1 ) 2 -(V y ) 2 λ 2 γ 1 γ 2 , ( 41a
) (V x ± λ γ 1 + γ 2 2 ) 2 -(V y ± λ) 2 λ 2 γ 2 2 , (41b)
(V x ± λ γ 1 + γ 2 2 ) 2 -(V y ∓ λ) 2 λ 2 γ 2 2 , (41c)
(V x ± λγ 2 ) 2 -(V y ) 2 λ 2 (s q -s xy ) 2 (4s 2 q -2s q s xy -s 2 xy + 8(s xy -s q )), (41d)
plus the analogous inequalities exchanging V x and V y . • if s xy = s q (segment [AD] on the figure 31), the L ∞ stability area is given by the proposition 3.4. • For all other s, no V corresponds to a L ∞ stable scheme.

In particular, the parameters (s q , s xy ) corresponding to a non empty area of L ∞ stability in V are included in the square [0, 2] 2 .

We now compare the L ∞ stability areas of the schemes relative to the velocities u = 0 or u = V with an intrinsic diffusion for different s. We draw the areas given by the propositions 3.3 and 3.5 on the figures 32 and 33.

The results are similar to the ones obtained in the non intrinsic case. We can't distinguish a scheme from an other in terms of L ∞ stability. The MRT scheme is better on the figure 32 but on the figure 33, the areas just intersect.

Weighted L 2 stability

In this section, we present some theoretical weighted L 2 stability results for the relative velocity D 2 Q 4 schemes. These results, based on the notion of stability proposed by Yong and al in [START_REF] Banda | A Stability Notion for Lattice Boltzmann Equations[END_REF], confirm the phenomena numerically observed in the section 2.3. In particular, taking u equal to the advection velocity V ("cascaded like" scheme) provides good stability features. These results also extend the L ∞ stability ones (section 3): the L ∞ stability areas are included in the L 2 ones.

The weighted L 2 stability notion

The general framework of the relative velocity D d Q q schemes, d, q ∈ N * , is chosen to present the stability notion introduced in [START_REF] Banda | A Stability Notion for Lattice Boltzmann Equations[END_REF] and studied in [START_REF] Junk | Weighted L 2 -stability of the lattice Boltzmann method[END_REF][START_REF] Rheinländer | On the stability structure for lattice Boltzmann schemes[END_REF].

An iteration of a lattice Boltzmann scheme splits into a transport step plus a relaxation step that is here linear since the equilibrium is linear [START_REF] Dellacherie | Construction and Analysis of Lattice Boltzmann Methods Applied to 1D Convection-Diffusion Equation[END_REF][START_REF] Dellar | Lattice kinetic schemes for magnetohydrodynamics[END_REF][START_REF] Humières | Generalized Lattice-Boltzmann equations[END_REF]. Given f (., t) the matrix composed of the distribution vectors taken on all the lattice points, we have

f (., t + ∆t) = T R( u)f (., t), t ∈ R,
where R( u) is the relaxation matrix and T the linear non local transport operator. The size of the matrix f (., t) is equal to the number of velocities q multiplied by the number of lattice points. In the case of the relative velocity D d Q q scheme, the collision matrix reads R( u) = I q + J ( u) where J ( u) is given by

J ( u) = M ( u) -1 DM ( u)(B -I q ),
where D = diag(s) is the diagonal matrix of the relaxation parameters and B defines the linear equilibrium thanks to f eq = Bf .

The calculus of the spectrum of the amplification operator T R( u) being difficult because of the matrix size, the idea is to study separately the transport and the collision: we want T and R( u) to be bounded by one in a well-chosen norm. The amplification operator T R( u) is then bounded by one and the scheme is stable for this norm. In the following, we say that an operator is stable in a given norm if the norm of this operator is bounded by one.

The idea proposed in [START_REF] Banda | A Stability Notion for Lattice Boltzmann Equations[END_REF] consists in weighting the L 2 norm defined by

|y| 2 = q-1 j=0 |y j | 2 1
2 , y = (y 0 , . . . , y q-1 ) ∈ R q , so that the transport is an isometry and the collision operator is stable in the weighted norm. The classical L 2 norm keeps the transport as isometric but the collision norm is difficult to evaluate. Introducing a weight allows to overcome this difficulty. First, we define a norm on R q depending on an invertible matrix P ∈ M q (R),

|y| P = |P y| 2 , y ∈ R q .
We can then define a norm for a matrix g of the size q multiplied by the number of lattice points with

|g| P ,L = x∈L |g(x)| 2 P 1 2 ,
where g(x) is the column of g associated with the node x. The necessity to define such a non local norm is due to the non locality of the transport. The matrix P is chosen such that t P P is diagonal. Thus the transport is an isometry for the norm | • | P ,L for the periodic and bounceback boundary conditions [START_REF] Junk | Weighted L 2 -stability of the lattice Boltzmann method[END_REF]. Indeed, the hypothesis of "quasi orthogonality" ( t P P diagonal) carrying on P cancels all the cross terms in the calculation of |f (., t)| P ,L . Then the isometry of the transport is obtained thanks to a simple change of variable: it uses the bijection between the nodes of the mesh and the nodes after the transport at a given velocity.

Contrary to the transport, the collision is a local operator so that its study reduces to the vectorial norm | • | P . If the matrix R( u) is stable in a particular node x of the lattice L for | • | P , it is stable for | • | P ,L . That's why, we use the local operator norm

||R( u)|| P = sup |x| P =1 |R( u)x| P ,
for the collision.

The matrix P is requested to diagonalize the collision: ||R( u)|| P is then easy to evaluate. These requirements define a notion of stability structure first evocated in [START_REF] Banda | A Stability Notion for Lattice Boltzmann Equations[END_REF].

Definition 4.1 ([25]

). A matrix N ∈ M q (R) has a pre-structure of stability if exists an invertible matrix P ∈ M q (R) and some vectors s, p ∈ R q so that P N P -1 = -diag(s),

t P P = diag(p),
where diag(p) is the diagonal matrix whose diagonal is constituted of the coefficients of p. The pre-structure of stability becomes a structure of stability if

s k ∈ [0, 2], 0 k q -1. (42)
Suppose that J ( u) has a pre-structure of stability and consider the collision in norm | • | P . Knowing that the eigenvalues of R( u) are 1 -s k for 0 k q -1,

||R( u)|| P = sup |P x| 2 =1 |P R( u)P -1 P x| 2 = ||P R( u)P -1 || 2 = max(|1 -s|),
where || • || 2 is the operator norm associated with | • | 2 . Thus the collision is stable for || • || P the condition (42) is verified. Under these conditions, the scheme is stable in the norm | • | P ,L since the transport is isometric.

We present a theorem from [START_REF] Rheinländer | On the stability structure for lattice Boltzmann schemes[END_REF] giving a necessary and sufficient condition of existence of a pre-structure of stability. In the following, this theorem is the tool used to obtain some stability results for our relative velocity schemes.

Théorème 4.1 ([25]

). A matrix N ∈ M q (R) has a pre-structure of stability if and only if there exists a diagonal positive definite matrix Λ ∈ M q (R) so that

N Λ = Λ t N . (43)
This criteria gives a practical criteria of existence of a pre-structure of stability through the resolution of a linear system in the coefficients of Λ. The size of this linear system is q(q -1)/2 since the matrix N Λ -Λ t N is antisymmetric.

The matrix P defining the weighted norm is explicitely derived in the proof of this theorem [START_REF] Rheinländer | On the stability structure for lattice Boltzmann schemes[END_REF]. We exhibit P through the proof of the sufficient condition of prestructure. The identity (43) is equivalent to the fact that Λ -1 N Λ is symmetric with respect to the scalar product

(x, y) Λ = q-1 i,j=0 Λ ij x i y j , x, y ∈ R q .
This implies the existence of an orthonormal matrix Q in the sense of (•, •) Λ diagonalizing Λ -1 N Λ. Then, the matrix P = (ΛQ) -1 diagonalizes the collision. The matrix t P P is diagonal because of the orthonormality of Q for the scalar product (•, •) Λ .

Stability results for the twisted relative velocity

D 2 Q 4 scheme
In this section, we apply the stability criteria given by the theorem 4.1 to obtain weighted L 2 results for the twisted relative velocity D 2 Q 4 scheme. Our motivation is to compare the cases u = 0 and u = V . To do so, the scheme must have at least two different relaxation parameters because the BGK scheme does not depend on the relative velocity since (3) is verified. We begin by studying the MRT scheme ( u = 0) whose results are limited. We then obtain more results for the cascaded like scheme ( u = V ).

The MRT scheme

We focus on the MRT scheme corresponding to u = 0. We limit to V = 0 because when u = 0 and V = 0, there is no pre-structure of stability. The framework does not distinguish the non intrinsic and intrinsic case because for V = 0, the equilibria ( 7) and ( 8) are identical.

Proposition 4.1 (L 2 stability for the MRT scheme). Consider the twisted D 2 Q 4 scheme relative to u = 0, of equilibrium m eq (0) = (ρ, 0, 0, 0) associated with the relaxation parameters s = (0, s q , s q , s xy ) ∈ R 4 . The associated matrix J (0) has a pre-structure of stability. Moreover, if 0 s q , s xy 2, then J (0) has a structure of stability. The collision matrix R(0) is then stable for || • || P and the scheme is stable in norm | • | P ,L .

Proof. Since u = V = 0, the matrix J (0) reads

J (0) = M (0) -1 D(E -I 4 )M (0),
with the diagonal matrix E = diag(1, 0, 0, 0). The matrix D(E -I 4 ) is then diagonal as a product of two diagonal matrices. Since M (0) -1 = t M (0)/4, the matrix J (0) = t M (0)D(E -I 4 )M (0)/4 is symmetric. So the identity is solution of the equation (43) and J (0) has a pre-structure of stability. The spectrum of J (0) being composed by 0, -s q , -s q , -s xy , the result on the structure of stability is obvious. 4.2.2. The scheme relative to the advection velocity ("cascaded like" scheme)

We now focus on the choice u = V beginning by the non intrinsic case. Proposition 4.2 (L 2 stability for a non intrinsic relative velocity scheme). Consider the twisted D 2 Q 4 scheme relative to u = V = (V x , V y ) ∈ R 2 , of equilibrium m eq (V ) = ρ(1, 0, 0, 0), associated with the relaxation parameters s = (0, s q , s q , s xy ) ∈ R 4 . The matrix J (V ) has a pre-structure of stability if and only if |V | ∞ < λ. Moreover, if 0 s q , s xy 2, then J (V ) has a structure of stability. The collision matrix R(V ) is then stable for || • || P and the scheme is stable in norm | • | P ,L .

Proof. According to the theorem 4.1, the existence of a pre-structure of stability for N = J (V ) is equivalent to the existence of a diagonal positive definite matrix Λ such that J (V )Λ = Λ t J (V ). A solution of this linear system is given by

Λ =     (λ + V x )(λ + V y ) 0 0 0 0 (λ -V x )(λ + V y ) 0 0 0 0 (λ -V x )(λ -V y ) 0 0 0 0 (λ + V x )(λ -V y )     .
This matrix is positive definite if and only if |V | ∞ < λ, that closes the first part of the proof. The eigenvalues of J (V ) are 0, -s q , -s q , -s xy because the matrix M (0)J (V )M (0) -1 is given by

    0 0 0 0 s q V x -s q 0 0 s q V y 0 -s q 0 (2s q -s xy )V x V y V y (s xy -s q ) V x (s xy -s q ) -s xy     .
We deduce that P J (V )P -1 = diag(0, -s q , -s q , -s xy ), and that J (V ) has a structure of stability if 0 s q , s xy 2. The norm || • || P of the collision operator R(V ) is equal to 1 and the scheme is stable in norm | • | P ,L .

These results extend the L ∞ stability ones obtained in the proposition 3.2. The proposition 3.2 gives the L ∞ stability area |V | ∞ λ for s q s xy min(1, 2s q ). The weighted L 2 notion generalizes it to 0 s q , s xy 2. This result was expected because the numerical experiments of the section 2.3 put in evidence the independence of the area of L 2 stability towards the relaxation parameters. The fact that we can not obtain the same type of result for the MRT scheme as for u = V is an other evidence of the good behaviour of this latter.

Let's also note that these stability conditions are defined by opened sets. The numerical experiments seem to show that the scheme is still L 2 stable on the closure of these sets. However, it is not possible to proove it with this notion because the matrix Λ is null on this closure. This proposition leads to a natural corollary for a single relaxation parameter (BGK) that has already been proved in [START_REF] Rheinländer | On the stability structure for lattice Boltzmann schemes[END_REF]. In the BGK case, the scheme does not depend on u. In particular, MRT and cascaded like scheme are identical. That's why the following result is valid whatever the relative velocity u.

Corollary 4.1 (The BGK non intrinsic case [START_REF] Rheinländer | On the stability structure for lattice Boltzmann schemes[END_REF]). For V = (V x , V y ) ∈ R 2 , consider the twisted relative velocity D 2 Q 4 BGK scheme of equilibrium

m eq ( u) = ρ 1, V x -u x , V y -u y , (V x -u x )(V y -u y ) ,
= V = (V x , V y ) ∈ R 2 , of equilibrium m eq (V ) = ρ(1, 0, 0, -V x V y ),
associated with the relaxation parameters s = (0, s q , s q , s xy ) ∈ R 4 . Suppose that V x = 0 (resp. V y = 0) then the associated matrix J (V ) has a pre-structure of stability if and only if |V y | < λ (resp. |V x | < λ). Moreover, if 0 s q , s xy 2, then J (V ) has a structure of stability. The collision matrix R(V ) is then stable for || • || P and the scheme is stable in norm | • | P ,L .

Proof. A non null solution of the equation ( 43) exists if and only if one of the two components of V is null. Suppose that it is V y , then a solution of ( 43) is given by

Λ =     λ + V x 0 0 0 0 λ -V x 0 0 0 0 λ -V x 0 0 0 0 λ + V x     .
This matrix is definite nonnegative if and only if |V x | < λ. The reasoning to get a structure of stability is identical to the one of the proposition 4.1.

This proposition extends the L ∞ stability results of the proposition 3.5 for the directions V x = 0 and V y = 0. The L ∞ notion provides areas decreasing as much as s q and s xy go far from each other. The | • | P ,L notion gives a constant optimal area in V for these directions while these parameters are bounded by 0 and 2. This phenomenon was observed numerically on the figures 11 to 16.

We close the section by a proposition for the BGK case that is not a corollary of the proposition for the intrinsic cascaded like scheme.

Proposition 4.4 (The BGK intrinsic case).

For V = (V x , V y ) ∈ R 2 , consider the twisted D 2 Q 4 scheme, of equilibrium m eq ( u) = ρ 1, V x -u x , V y -u y , u x u y -u x V y -u y V x ),
BGK of relaxation parameter s ∈ R. We note that this proposition generalizes the L ∞ stability results in the BGK case (proposition 3.4). This weighted L 2 notion extends the set of s corresponding to the stability area |V | 1 < λ to 0 s 2.

Interpretation of the results

When we compare these propositions to the numerical results obtained in the section 2.3, we notice that this stability notion is only usable when the parameters s and V are decoupled. For the scheme relative to u = V with an intrinsic diffusion, we have only theoretical results for V x = 0 ou V y = 0, corresponding to a case where the areas in V are independent of s. Instead, when the area in V is a function of s, it is not possible to build a pre-structure of stability any more. For example, there is no pre-structure of stability when V x and V y are different from 0 for the D 2 Q 4 scheme relative to u = V with an intrinsic diffusion: for this scheme, the numerical results of the section 2.3 exhibit a link between V and s. Similarly, the MRT scheme present stability areas in V depending on s and it is impossible to build a pre-structure of stability for V = 0.

The origin of this limitation seems to be the hypothesis of diagonalization of the collision. Indeed, the spectrum of this operator being equal to 0, -s q , -s q , -s xy , the existence of a pre-structure of stability with such a diagonalization implies automatically the stability of the scheme for 0 s q , s xy 2. Then V and s can't be linked because the former hypothesis is too requiring. The linear system of six equations and four unknowns for the D 2 Q 4 scheme is an evidence of the constraints imposed by the notion. The existence of a pre-structure of stability requires its rank to be three that imposes V and s to be decoupled here. For example, when u = 0 and V = 0, there is no solution positive definite to the equation ( 43), because this rank is greater than four. That's why it is not possible to show similar results than in the case u = V which requires the resolution of a rank three linear system.

However, this notion gives some promising results as those exposed here and in [START_REF] Junk | Weighted L 2 -stability of the lattice Boltzmann method[END_REF][START_REF] Rheinländer | On the stability structure for lattice Boltzmann schemes[END_REF]. It seems to describe well the limit of blow-up of a scheme. The main purpose now is to obtain theoretical stability results in more complex cases. Particularly, we want to complete this study for the D 2 Q 4 especially when u = 0. To do so, we must be able to relax the constraint of diagonalization of the collision without penalizing the isometry of the transport. Maybe, it would be possible to find some matrix P so that the norm | • | P ,L of the collision operator is computable without diagonalizing it.

Conclusion

We have studied the stability of a four velocities relative velocity scheme with two different equilibria for a linear advection equation. The discussion is based on two notions of L ∞ and weighted L 2 stability. It carries on the choice of the relative velocity equal to 0 (MRT scheme) or to the advection velocity ("cascaded like" scheme). The main conclusions are the following: comparing MRT and relative velocity schemes, no scheme is "better than the other" in terms of L ∞ stability; the stability areas generally just intersect. Instead, in terms of L 2 norm, the scheme relative to the advection velocity is more stable than the MRT scheme. This improvement is correlated to the cancellation of some dispersive terms of the equivalent equations for the scheme relative to the advection velocity. These terms stay for u = 0 and are linked to instabilities at low diffusion. This result has been justified theoretically for the scheme relative to the advection velocity thanks to the weighted L 2 notion of stability. Further work needs to be done to obtain more precise theoretical results for the MRT scheme and for systems of conservation laws. The study also needs to be generalized to an arbitrary relative velocity u, the present work being adapted to the two choices u = 0 and u = V . and

R(XY )(V) = 3 j=0 (Rv j ) x (Rv j ) y f j = 3 j=0 ((v x j ) 2 -(v y j ) 2 )f j = (X 2 -Y 2 )(V). (46)
Providing these images, we can write the transformation of the relaxation of the D 2 Q 4 relative velocity scheme by R. We here choose to keep the polynomial notations to express the moments. The relaxation of the D 2 Q 4 reads

X * (V) = X(V) + s q ((V x -u x )ρ -X(V)), Y * (V) = Y (V) + s q ((V y -u y )ρ -Y (V)), (X 2 -Y 2 ) * (V) = (X 2 -Y 2 )(V) + s xy ((V x -u x ) 2 -(V y -u y ) 2 )ρ -(X 2 -Y 2 )(V) .
Introducing the relations (45,46), we obtain

R(X) * (V) = R(X)(V) + s q (((RV ) x -(R u) x )ρ -R(X)(V)), R(Y ) * (V) = R(Y )(V) + s q (((RV ) y -(R u) y )ρ -R(Y )(V)), R(XY ) * (V) = R(XY )(V) + s xy ((RV ) x -(R u) x )((RV ) y -(R u) y )ρ-R(XY )(V) ,
that is the relaxation of the twisted D 2 Q 4 relative velocity scheme. Note that this last step uses the fact that X and Y have the same relaxation parameter s q .

Proposition A.1 (Relation between the stability areas). Let (s q , s xy ) ∈ R 2 , note S nt ⊂ R 2 (resp. S t ⊂ R 2 ), the set of the velocities V ∈ R 2 such that the relative velocity D 2 Q 4 scheme (resp. twisted relative velocity D 2 Q 4 scheme) of equilibrium given by (44) (resp. (6,7)) and of relaxation parameters s = (0, s q , s q , s xy ) is L ∞ or L 2 stable. We have

S t = RS nt . (47)
Proof. Note R t (V , u, s) ∈ M 4 (R), the relaxation operator of the twisted relative velocity D 2 Q 4 scheme of equilibrium (6,7): it verifies

f = R t (V , u, s)f . Note R nt (V , u, s) ∈ M 4 (
R) its analogue for the D 2 Q 4 scheme associated with the equilibrium (44). The lemma A.1 means that

R nt (V , u, s) = R t (RV , R u, s). (48)
The transport does not influence the L ∞ stability because it only exchanges the particle distributions. Thus S nt is the set of the velocities V ∈ R 2 so that the matrix R nt (V , u, s) is nonnegative. According to the relation (48), it corresponds to the velocities V ∈ R 2 so that RV belongs to the L ∞ stability area of the twisted scheme: that is equivalent to the relation (47) and the lemma is proven in the L ∞ case.

• if 0 s xy min(s q , 2 -s q ), the scheme is L ∞ stable for all V such that

|V | ∞ λγ 2 .
• if s q min(1, s xy ) and s xy 2s q the scheme is L ∞ stable for all V such that

|V | ∞ λ(2 -γ) 2 . • if s q max(1, 2 -s xy ) and s xy 2(2 -s q ) the scheme is L ∞ stable for all V such that |V | ∞ λ 4 s q -2 -γ 2 .
• if s xy = s q = 0, the scheme is unconditionally L ∞ stable.

• For all other s, no V corresponds to a L ∞ stable scheme.

Proposition B.4. Let V ∈ R 2 , (s q , s xy ) ∈ R 2 , consider the D 2 Q 4 scheme relative to u = V associated with the relaxation parameters (0, s q , s q , s xy ), and the equilibrium

m eq (V ) = ρ(1, 0, 0, (V y ) 2 -(V x ) 2 ).
Note γ 1 = (2s q -s xy )/(s q -s xy ) and γ 2 = s xy /(s q -s xy ),

• if s q < s xy min(2s q , 2(2 -s q )), the scheme is stable for all V such that

(V x ± λγ 1 2 ) 2 -(V y ) 2 λ 2 γ 1 γ 2 4 , (V x ± λ γ 1 + γ 2 4 ) 2 -(V y ∓ λ 2 ) 2 λ 2 γ 2 2 4 , (V x ± λ γ 1 + γ 2 4 ) 2 -(V y ± λ 2 ) 2 λ 2 γ 2 2 4 , (V x ± λγ 2 2 ) 2 -(V y ) 2 λ 2 4(s q -s xy ) 2 (4s 2 q -2s q s xy -s 2 xy + 8(s xy -s q )),
plus the analogous inequalities exchanging V x and V y . • if s xy min(s q , 2(2 -s q )), the scheme is stable for all V such that

(V x ± λγ 1 2 ) 2 -(V y ) 2 λ 2 γ 1 γ 2 4 , (V x ± λ γ 1 + γ 2 4 ) 2 -(V y ∓ λ 2 ) 2 λ 2 γ 2 2 4 , (V x ± λ γ 1 + γ 2 4 ) 2 -(V y ± λ 2 ) 2 λ 2 γ 2 2 4 , (V x ± λγ 2 2 ) 2 -(V y ) 2 λ 2 4(s q -s xy ) 2 (4s 2
q -2s q s xy -s 2 xy + 8(s xy -s q )), plus the analogous inequalities exchanging V x and V y .

• if s xy = s q 1, the scheme is L ∞ stable for all V such that

|V | ∞ λ 2 . • if 1 s xy = s q 4/3, the scheme is L ∞ stable for all V such that |V | ∞ λ 2 4 s - 3 .
• For all other s, no V corresponds to a L ∞ stable scheme.

Appendix C. Proof of the proposition 3.5

The positivity of the matrix shifting from f to f is equivalent to

λ(2s q -s xy )(λ ± (V x + (-1) i V y )) + 2(-1) i (s q -s xy )V x V y 0, (52a)
λ 2 s xy ± λ(s xy V x + (-1) i (2s q -s xy )V y ) + 2(-1) i (s q -s xy )V x V y 0, (52b)
λ 2 s xy ± λ(s xy V y + (-1) i (2s q -s xy )V x ) + 2(-1) i (s q -s xy )V x V y 0, (52c)
4λ 2 -((2s q + s xy )λ 2 ± s xy λ(V x + (-1) i V y ) -2(-1) i (s q -s xy )V x V y ) 0. (52d)
for i = 0, 1. We begin by the case s q = s xy of a BGK scheme corresponding to the proposition 3.4. The four inequations (52a) to (52d) become

λs q (λ ± (V x + (-1) i V y )) 0, (53)
s q (λ ± (V x + (-1) i V y )) 0, (54)
s q (λ ± (V y + (-1) i V x )) 0, (55)
4λ 2 -s q (3λ 2 ± λ(V x + (-1) i V y ) 0. (56)
The case s q 0 is impossible because no velocity V would verify these four inequalities requiring V x ± V y to be simultaneously greater than λ and inferior to -λ. Then s q 0 that leads to

|V | 1 min(λ, λ( 4 s q -3)),
that gives exactly the proposition 3.4 according to the choice of s q . If s xy = s q , we can write the inequalities (52a) to (52d) thanks to γ 1 and γ 2 . If s xy < s q , we obtain

λ γ 1 2 (λ ± (V x + (-1) i V y )) + (-1) i V x V y 0, i = 0, 1, (57a)
λ 2 γ 2 2 ± λ 2 (γ 2 V x + (-1) i γ 1 V y ) + (-1) i V x V y 0, i = 0, 1, (57b)
λ 2 γ 2 2 ± λ 2 (γ 2 V y + (-1) i γ 1 V x ) + (-1) i V x V y 0, i = 0, 1, (57c) 4
-(2s q + s xy ) 2(s q -s xy ) λ 2 ± γ 2 2 λ(V x + (-1) i V y ) + (-1) i V x V y 0, i = 0, 1. (57d)
We do the following change of variable V

x = V x + V y , V y = V x -V y . Noticing that γ 1 -γ 2 = 2, the inequalities (57a) to (57d) are equivalent to 2λγ 1 (λ ± V x ) + (V x ) 2 -(V y ) 2 0, (58a)
2λ 2 γ 2 ± 2λ( γ 1 + γ 2 2 V x -V y ) + (V x ) 2 -(V y ) 2 0, (58b)
2λ 2 γ 2 ± 2λ( γ 1 + γ 2 2 V x + V y ) + (V x ) 2 -(V y ) 2 0, (58c) 8
-(4s q + 2s xy ) s q -s xy λ 2 ± 2γ 2 λV x + (V x ) 2 -(V y ) 2 0, (58d)
for i = 0 in addition to the ones obtained exchanging V x and V y corresponding to i = 1. We center these inequalities to get the identities from (41a) to (41d). The case s q s xy is obtained by changing the sense of the inequalities. It is characterized by the inequalities from (40a) to (40d).

We begin by treating the case s q < s xy corresponding to the inequalities from (40a) to (40d). The case s q > s xy will be considered further. The reasoning first eliminates the couples (s q , s xy ) leading to no stable velocity and then concentrates on stable areas. First, we assume that s xy > max(0, 2s q ) and show that there is no stable velocity V . Since s xy > s q , we have γ 1 γ 2 0 and the equation (40a) reads

(V y ) 2 -(V x ± λγ 1 ) 2 -λ 2 γ 1 γ 2 ,
plus the two ones obtained exchanging V x and V y . The intersection of these four hyperbolic areas is empty.

Secondly we eliminate the case s q < s xy < 0. In this situation, γ 1 γ 2 is nonnegative (γ 1 and γ 2 are both nonnegative) and (40a) has a solution if and only if the abscissa of the right pole of the hyperbole centered in -λγ 1 is nonnegative. This reads -λγ 1 + λ(γ 1 γ 2 )

1 2
0 ⇐⇒ γ 1 γ 2 ⇐⇒ s xy s q .

The last inequality contradicts our hypothesis and there is no stable velocity when s q < s xy < 0.

It remains to treat the case s q < s xy 2s q for which γ 1 γ 2 0. The equation (40a) has some solutions V if and only if the abscissa of the right pole P of the hyperbole centered in λγ 1 is nonnegative

λγ 1 + λ(γ 1 γ 2 ) 1 2

0.

(59)

The area satisfying these inequations corresponds to the points A, B, C, D on the figure 35. It increases as the abscissa of P increases. The inequality (59) is equivalent to γ 1 0, that is verified since s q < s xy 2s q .

For the equation (40d), there are two cases. Hence 4s 2 q -2s q s xy -s 2 xy + 8(s xy -s q ) 0, q -2s q s xy -s 2 xy +8(s xy -s q ) = 0, and straight lines s xy = 2 -s q and s xy = s q .

-λγ 1 λγ 1 -λγ 1 then the reasoning is the same as the one made in the case s xy > max(0, s q ) to show that there is no stable velocity. Hence 4s 2 q -2s q s xy -s 2 xy + 8(s xy -s q ) 0, then the abscissa of the right pole of the hyperbole centered in λγ 2 must be nonnegative: λγ 2 + λ |s q -s xy | (4s 2 q -2s q s xy -s 2 xy + 8(s xy -s q )) 1 2

0.

The parameter γ 2 being negative, this is equivalent to 1 (s q -s xy ) 2 (4s 2 q -2s q s xy -s 2 xy + 8(s xy -s q )) γ 2 2 , and so to (s q -s xy )(4s q + 2s xy -8) 0 ⇐⇒ s xy 2(2 -s q ), since s q s xy .

We now summarize the results we have obtained for s q < s xy . We have proven that if s q < s xy , the inequalities (40a) and (40d) have solutions for if s xy min(2s q , 2(2-s q )). It remains to prove that (40c) and (40b) have also solutions in this case.

The conditions (40c) and (40b) lead to a non empty area in V since the intersection I (figure 36) between the hyperbole given by (40c) centered in (-λ(γ 1 + γ 2 )/2, λ) and the straight line V y = 0 has a negative abscissa. The same condition on the intersection of the hyperbole centered in (-λ(γ 1 +γ 2 )/2, -λ) and V y = 0 is requested for the equation (40b). We focus on the hyperbole centered in (-λ(γ 1 + γ 2 )/2, λ) associated with (40c), because the one associated with (40b) leads to the same result. The intersections of the hyperbole of equation

(V x + λ γ 1 + γ 2 2 ) 2 -(V y -λ) 2 = λ 2 γ 2 2 ,
with V y = 0 verify

(V x + λ γ 1 + γ 2 2 ) 2 -λ 2 = λ 2 γ 2 2 .
Their abscissas are then given by -λ γ 1 + γ 2 2 ± λ(γ 2 2 + 1)

1 2 .
We check when the abscissa of the point I of the figure 36 is negative

-λ γ 1 + γ 2 2 -λ(γ 2 2 + 1) 1 2 0 ⇐⇒ γ 1 + γ 2 2 2 γ 2 2 + 1 ⇐⇒ s 2 q s 2 xy + (s q -s xy ) 2 .
This is equivalent to s xy (s xy -s q ) 0, that is valid because s xy s q .

We now focus on the case s xy < s q corresponding to the inequalities from (41a) to (41d). Let's begin by eliminating the case 2s q s xy < s q . In this case (41a) has solutions if and only if both poles of the hyperbole centered in -λγ 1 (γ 1 0) have nonnegative abscissas that reads -λ(γ 1 + (γ 1 γ 2 ) 1 2 ) 0. Because of the negativity of γ 1 , this is equivalent to s q s xy that contradicts our hypothesis. Second, we eliminate the area 2(2 -s q ) s xy < s q . This is done showing that there is no solution to (41d). In this area, we have 4s 2 q -2s q s xy -s 2 xy + 8(s xy -s q ) 0, (figure 34) so that (41d) has a solution if and only if the abscissas of both poles of the hyperbole centered in λ|γ 2 | are nonnegative. This reads λ|γ 2 | -λ |s q -s xy | (4s 2 q -2s q s xy -s 2 xy + 8(s xy -s q )) 1 2

0,

that is equivalent to s xy 2(2 -s q ). This contradicts our hypothesis.

It remains to eliminate the area s xy < min(0, s q , 2(2 -s q )) corresponding to γ 1 0 and γ 2 0. If s q 0, then γ 1 + γ 2 2 0 and the inequality (41c) is non empty if the abscissas of the intersections of the hyperbole centered in (-λ(γ 1 + γ 2 )/2, λ) with V y = 0 are both nonnegative. This is equivalent to

-λ γ 1 + γ 2 2 -λ(γ 2 2 + 1) 1 2 0 ⇐⇒ γ 1 + γ 2 2 2 γ 2 2 + 1 ⇐⇒ s 2 q s 2 xy + (s q -s xy ) 2 .
This is equivalent to s xy (s xy -s q ) 0, that is verified for s xy s q since s xy < 0. This enters in contradiction with our hypothesis. If s q > 0, the reasoning is symmetric with respect to the ordinates axis (V x = 0) and leads to the same contradiction.

To close the proof, it remains to guarantee the existence of at least one velocity stable in the triangle defined by 0 s xy min(s q , 2(2 -s q )). The equation (41a) has a non empty stability area if and only if the two poles of the hyperbole centered in λγ 1 have a nonnegative abscissa: this case, illustrated by the figure 19, is equivalent to

λ(γ 1 -(γ 1 γ 2 ) 1 2 ) 0.
This inequality is equivalent to γ 1 0, that is true for s xy s q .

For the equation (41d), we distinguish two cases: the first is 4s 2 q -2s q s xy -s 2 xy + 8(s xy -s q ) 0. There is a non empty stability area if and only if the two poles of the hyperbole centered in λγ 2 have nonnegative abscissas (represented on the figure 19): this reads λγ 2 -λ |s q -s xy | (4s 2 q -2s q s xy -s 2 xy + 8(s xy -s q ))

1 2

0.

This inequality is equivalent to s xy 2(2 -s q ). Otherwise, 4s 2 q -2s q s xy -s 2 xy + 8(s xy -s q ) 0, and the equation (41d) reads

(V y ) 2 -(V x ± λγ 2 ) 2 - λ 2 (s q -s xy ) 2 (4s 2
q -2s q s xy -s 2 xy + 8(s xy -s q )).

This inequation has always solutions: the area of stability in V is analogous to the one represented on the figure 24.

The inequations (41c) and (41b), illustrated on the figure 37, have solutions since the two intersections of the hyperbole centered in (λ(γ 1 + γ 2 )/2, -λ) with V y = 0 have nonnegative abscissas λ γ 1 + γ 2 2 -λ(γ 2 2 + 1)

1 2 0 ⇐⇒ γ 1 + γ 2 2 2 γ 2 2 + 1 ⇐⇒ s 2 q s 2 xy + (s q -s xy ) 2 .
This is equivalent to s xy (s xy -s q ) 0, that is verified for s xy s q .

λ γ 1 + γ 2 2 -λ γ 1 + γ 2 2 λ -λ 0 A B C D V x V y
Figure 37. Layout of the bounds of the L ∞ stability areas corresponding to (41c) and (41b) for the twisted scheme relative to u = V with an intrinsic diffusion.

Figure 1 .

 1 Figure 1. Velocities of the twisted D 2 Q 4 scheme.

Figure 2 .

 2 Figure 2. Well posedness areas in V for the twisted D 2 Q 4 schemes: with intrinsic diffusion (circle of radius λ), with a non intrinsic diffusion (square [-λ, λ] 2 ).

Figure 4 .

 4 Figure 4. Advected spot of velocity V = (0.9λ, 0) for the twisted D 2 Q 4 scheme relative to u = V with an intrinsic diffusion at t = 0 (left) and t = 0.4 (right) for σ xy = 1/ √ 3, σ q = 1/20.

  Figure 5. Velocities V stable in λ scale for the scheme relative to u = 0 (MRT on the left), u = V (on the right), for a non intrinsic diffusion with s q = 1 and s xy = 1.

Figure 6 .

 6 Figure 6. Velocities V stable in λ scale for the scheme relative to u = 0 (MRT on the left), u = V (on the right), for a non intrinsic diffusion with s q = 1 and s xy = 1.5.

  Figure 7. Velocities V stable in λ scale for the scheme relative to u = 0 (MRT on the left), u = V (on the right), for a non intrinsic diffusion with s q = 1 and s xy = 1.9.

Figure 8 .

 8 Figure 8. Velocities V stable in λ scale for the scheme relative to u = 0 (MRT on the left), u = V (on the right), for a non intrinsic diffusion with s q = 1 and s xy = 2.

Figure 10 .

 10 Figure 10. Velocities V stable in λ scale for the scheme relative to u = 0 (MRT on the left), u = V (on the right), for a non intrinsic diffusion with s q = 2 and s xy = 1.

Figure 12 .

 12 Figure 12. Velocities V stable in λ scale for the scheme relative to u = 0 (MRT on the left), u = V (on the right), for an intrinsic diffusion with s q = 1 and s xy = 1.5.

  Figure 13. Velocities V stable in λ scale for the scheme relative to u = 0 (MRT on the left), u = V (on the right), for an intrinsic diffusion with s q = 1 and s xy = 1.9.

Figure 14 .

 14 Figure 14. Velocities V stable in λ scale for the scheme relative to u = 0 (MRT on the left), u = V (on the right), for an intrinsic diffusion with s q = 1 and s xy = 2.

Figure 16 .

 16 Figure[START_REF] Geier | Cascaded digital Lattice Boltzmann automata for high Reynolds number flow[END_REF]. Velocities V stable in λ scale for the scheme relative to u = 0 (MRT on the left), u = V (on the right), for an intrinsic diffusion with s q = 2 and s xy = 1.

Figure 17 .

 17 Figure 17. L ∞ stability areas in s of the twisted D 2 Q 4 MRT scheme with anistropic diffusion. BCD : (17a,17b), ABC : (18a,18b), ACD : (19a,19b).

Figure 18 .

 18 Figure 18. Value of the maximum of (21). Area 1 : γ 2 -1, area 2 : (γ + 1) 2 -4/s xy , area 3 : (γ -1) 2 .

Figure 19 .

 19 Figure 19. Layout of the bounds of the stability areas in V for the twisted D 2 Q 4 MRT scheme.

Figure 20 .

 20 Figure 20. L ∞ stability area in s of the twisted D 2 Q 4 scheme relative to u = V with a non intrinsic diffusion. ACF :[START_REF] Sterling | Stability Analysis of Lattice Boltzmann methods[END_REF], ABC :[START_REF] Marié | Comparison between Lattice Boltzmann method and Navier-Stokes high order schemes for computational aeroacoustics[END_REF], BCED : (25a,25b), CEF : (25a).

Figure 22 .

 22 Figure 22. L ∞ stability area for the scheme relative to u = V with a non intrinsic diffusion for s xy = 2s q . The point P is of abscissa λ 2λ(2 -s xy )/s xy 2λ.

  ) represented on the figure 22, let s xy 4/3, then 2λ (2 -s xy ) s xy λ,

Figure 28 .

 28 Figure 28. L ∞ stability areas for s = (0, 1, 1, 3/2). Scheme relative to u = V (A, B, C, D), to u = 0 (E, F, G, H), with a non intrinsic diffusion.

Figure 29 .

 29 Figure 29. L ∞ stability areas for s = (0, 3/2, 3/2, 3/4). Scheme relative to u = V (A, B, C, D), to u = 0 (E, F, G, H), with a non intrinsic diffusion.

Figure 30 .

 30 Figure 30. L ∞ stability area in s of the twisted D 2 Q 4 MRT scheme ( u = 0) with an intrinsic diffusion. BCD : (35), ABC : (36), ACD : (37).

Figure 31 .

 31 Figure 31. L ∞ stability area in s of the twisted D 2 Q 4 scheme relative to u = V with an intrinsic diffusion. ACD : (40a) to (40d), ABD : (41a) to (41d).

Figure 32 .Figure 33 .

 3233 Figure 32. L ∞ stability areas for s = (0, 1, 1, 3/2). Scheme relative to u = V (A, B, C, D), u = 0 (E, F, G, H), with an intrinsic diffusion.

Figure 34 .

 34 Figure 34. Hyperbole 4s 2q -2s q s xy -s 2 xy +8(s xy -s q ) = 0, and straight lines s xy = 2 -s q and s xy = s q .

Figure 35 .

 35 Figure35. Layout of the bounds of the L ∞ stability areas corresponding to (40a) and (40d) for the twisted D 2 Q 4 scheme relative to u = V with an intrinsic diffusion.

Figure 36 .

 36 Figure36. Layout of the bounds of the L ∞ stability areas corresponding to (40c) and (40b) for the twisted D 2 Q 4 scheme relative to u = V with an intrinsic diffusion.

Table 1 .

 1 Maximal |V | stable in λ scale for θ = 0, σ xy = 1/

	σ q	1/10	1/20	1/50	1/100	1/200
	D 2 Q 4 non intrinsic u = 0	1.00	0.80	0.49	0.34	0.23
	D 2 Q 4 non intrinsic u = V	1.00	1.00	1.00	1.00	1.00
	D 2 Q 4 intrinsic u = 0	1.00	0.79	0.48	0.33	0.23
	D 2 Q 4 intrinsic u = V	1.00	1.00	1.00	1.00	1.00
					√	3 and
	128 2 points.
	σ q	1/10	1/20	1/50	1/100	1/200
	D 2 Q 4 non intrinsic u = 0	1.41	0.80	0.42	0.28	0.20
	D 2 Q 4 non intrinsic u = V	1.41	1.41	1.41	1.41	1.41
	D 2 Q 4 intrinsic u = 0	0.75	0.56	0.36	0.26	0.18
	D 2 Q 4 intrinsic u = V	0.86	0.76	0.65	0.59	0.53

Table 2 .

 2 Maximal |V | stable in λ scale for θ = π/4, σ xy = 1/

	√	3 and

  associated with the relaxation parameter s ∈ R. The matrix J ( u) has a pre-structure of stability if and only if |V | ∞ < λ. Moreover, if 0 s 2, then J ( u) has a structure of stability. The collision matrix R( u) is then stable for || • || P and the scheme is stable in norm | • | P ,L . Proposition 4.3 (L 2 stability for an intrinsic relative velocity scheme). Consider the twisted D 2 Q 4 scheme relative to u

  The associated matrix J ( u) has a pre-structure of stability if and only if |V | 1 < λ. Moreover, if 0 s 2, then J ( u) has a structure of stability. The collision matrix R( u) is then stable for || • || P and the scheme is stable in norm | • | P ,L .

Appendix A. Link between the twisted D 2 Q 4 and the D 2 Q 4 schemes

We proove that the L ∞ and L 2 stability areas of the twisted D 2 Q 4 and D 2 Q 4 relative velocity schemes with a non intrinsic diffusion correspond through the composition of a rotation and a homothety. The same result is still true in the intrinsic case.

If we note v the velocity set of the D 2 Q 4 scheme, the twisted set of velocities is given by Rv where

is the composition of a rotation and a homothety of scale factor √ 2. This transformation allows to link the relaxation of both schemes.

Lemma A.1 (Relation between the relaxation operators). Let V ∈ R 2 , u ∈ R 2 , (s q , s xy ) ∈ R 2 , the relative velocity D 2 Q 4 scheme associated with the equilibrium

for the relaxation parameters (0, s q , s q , s xy ) and the relative velocity twisted D 2 Q 4 scheme associated with the equilibrium

and with the relaxation parameters (0, s q , s q , s xy ) have the same relaxation phase.

Proof. We have defined a transformation R sending the velocities of the D 2 Q 4 scheme on the velocities of the twisted D 2 Q 4 scheme. This naturally leads to the following application also denoted by R.

where R(P ) :

Instead, the transport plays a role for the L 2 stability. It becomes local in the Fourier space: the transport matrices associated with the D 2 Q 4 scheme (A nt ) and the twisted one (A t ) are given by

is the composition of a rotation and a homothety of scale factor √ 2. The amplification matrices L nt and L t of both schemes are then defined by

We then have the following identity

The L 2 stability needs to evaluate the maximum of the spectral radius r of this matrix for all the wavenumbers k ∈ R 2 . According to (49),

the last equality coming from a variable change in R 2 . As for the case of the L ∞ stability, only the velocities V finally matter: the relation between the L 2 stability areas is then obtained in the same way.

Appendix B. Theoretical L ∞ stability results for the D 2 Q 4 scheme B.1. For a non intrinsic diffusion Proposition B.1 (L ∞ stability areas for the MRT scheme). Let V ∈ R 2 , (s q , s xy ) ∈ R 2 , consider the D 2 Q 4 of MRT scheme with the relaxation parameters (0, s q , s q , s xy ), associated with the equilibrium

• if 0 < s xy min(s q , 2 -s q ), the scheme is L ∞ stable for all V so that

• if s q s xy 2s q and s q 1, the scheme is L ∞ stable for all V so that

• if 2 -s q s xy 2(2 -s q ) and s q 1, the scheme is L ∞ stable for all V so that

• if s xy = 0 and 0 < s q 2, the scheme is L ∞ stable for V = 0.

• if s xy = s q = 0, the scheme is unconditionally L ∞ stable.

• For all other s, there is no V corresponding to a L ∞ stable scheme.

Proposition B.2. Let V ∈ R 2 , (s q , s xy ) ∈ R 2 , consider the D 2 Q 4 scheme relative to u = V associated with the relaxation parameters (0, s q , s q , s xy ), and the equilibrium m eq (V ) = ρ(1, 0, 0, 0).

Let us note γ = s xy /(2s q -s xy ),

• if s q s xy min(1, 2s q ), the scheme is L ∞ stable on

• if s xy < 2s q , max(s q , 1) s xy 2(2 -s q ), the scheme is L ∞ stable on the intersection of (50) with

• if 0 s xy min(s q , s q (2 -s q )), the scheme is L ∞ stable on |V | 1 λγ.

• if s q (2 -s q ) s xy min(s q , 2(2 -s q )), the scheme is L ∞ stable on (51a,51b).

• if s xy = 2s q and 1 < s xy 2, the scheme is L ∞ stable on the intersection of (50) and

• if s xy = s q = 0, the scheme is unconditionally L ∞ stable.

• For all other s, no V corresponds to a L ∞ stable scheme.

B.2. For an intrinsic diffusion

Proposition B.3. Let V ∈ R 2 , (s q , s xy ) ∈ R 2 , consider the D 2 Q 4 MRT scheme associated with the relaxation parameters (0, s q , s q , s xy ), and the equilibrium m eq (0) = ρ(1, V x , V y , 0). Note γ = s xy /s q ,