A higher-order moment method of the lattice Boltzmann model for the conservation law equation

Applied Mathematical Modelling 34 (2010) 481–494

Contents lists available at ScienceDirect

Applied Mathematical Modelling journal homepage: www.elsevier.com/locate/apm

A higher-order moment method of the lattice Boltzmann model for the conservation law equation q Yinfeng Dong, Jianying Zhang, Guangwu Yan * College of Mathematics, Jilin University, Changchun 130012, PR China

a r t i c l e

i n f o

Article history: Received 21 October 2007 Received in revised form 6 May 2009 Accepted 1 June 2009 Available online 16 June 2009

Keywords: Lattice Boltzmann model Higher-order moment method Conservation law equation

a b s t r a c t In this paper, we proposed a higher-order moment method in the lattice Boltzmann model for the conservation law equation. In contrast to the lattice Bhatnagar–Gross–Krook (BGK) model, the higher-order moment method has a wide flexibility to select equilibrium distribution function. This method is based on so-called a series of partial differential equations obtained by using multi-scale technique and Chapman–Enskog expansion. According to Hirt’s heuristic stability theory, the stability of the scheme can be controlled by modulating some special moments to design the third-order dispersion term and the fourth-order dissipation term. As results, the conservation law equation is recovered with higher-order truncation error. The numerical examples show the higher-order moment method can be used to raise the accuracy of the truncation error of the lattice Boltzmann scheme for the conservation law equation. Ó 2009 Published by Elsevier Inc.

1. Introduction In recent years, computational methods based on the lattice Boltzmann method (LBM) have attracted much attention. They have been developed as an alternative method for computational fluid dynamics (CFD). These lattice Boltzmann methods originated from a Boolean fluid model known as the lattice gas automata (LGA) [1] originally developed to overcome certain drawbacks such as the presence of statistical noise and lack of Galilean invariance of LGA for modeling fluid flow based upon kinetic theory [2–5]. The lattice Boltzmann models have been used to simulate formidable problems such as multiphase flows [6–11], multi-component flows [7,12–14], porous media flows [15,16], flows of suspensions [17,18], and compressible flows [19–28]. Additionally, the lattice Boltzmann models have been developed to simulate linear and nonlinear partial differential equations such as wave motion equation [29], Burgers equation [30], KdV equation [31], Lorenz equations [32], the shallow equation [33], the Richards equation [34], Poisson equation [35], and the nonlinear Schrödinger equation [36–40]. Unlike conventional methods based on macroscopic continuum equation, these LBMs start from mesoscopic kinetic equation, say, the lattice Boltzmann equation, to determine macroscopic fluid flows. Their kinetic nature brings certain advantages over conventional numerical methods, such as their algorithmic simplicity, parallel computation, easy handing of complex boundary conditions, and efficient numerical simulations. Especially their algorithmic simplicity and the flexibility to select equilibrium distribution functions are outstanding advantages. There are two basic problems in the lattice Boltzmann method: (1) how to improve the accuracy of lattice Boltzmann model, (2) how to obtain higher-order accuracy boundary conditions for some complex flows. In this paper, we focus on

q This work is Project 20092005 supported by Gradate Innovation Fund of Jilin University, National Nature Science Foundation of China (Grant No. 90305013), and the Chuangxin Foundation of Jilin University (No. 2004CX041). * Corresponding author. E-mail address: [email protected] (G. Yan).

0307-904X/$ - see front matter Ó 2009 Published by Elsevier Inc. doi:10.1016/j.apm.2009.06.024

482

Y. Dong et al. / Applied Mathematical Modelling 34 (2010) 481–494

the first problem. Furthermore, we choose one of the conversation law equations the mass conservation equation to construct a high accuracy LBM model. The one-dimensional mass conservation equation is

@u @FðuÞ þ ¼ 0; @t @x

ð1Þ

where u expresses mass, F(u) is the mass flux. In Eq. (1), the mass flux F(u) is a nonlinear function of mass u. The variable u can also express momentum and energy, and then flux F(u) represents impulse flux and energy flux respectively. In order to propose a higher-order accuracy model for the conservation law equation, the key step is to supply a series of partial differential equations in different time scales based on the Chapman–Enskog analysis [41–44]. The Chapman–Enskog analysis consists of the Chapman–Enskog expansion and multi-scale time and multi-scale space techniques. We can obtain the macroscopic equations by combining the equations in the different scales [41,42]. In Ref. [43], we introduce four partial differential equations to construct a LBM model. However, the four partial differential equations are not competent for constructing the LBM model with the higher-order dispersion term and dissipation term. A detailed analysis for the generalized lattice Boltzmann equation with multiple-relaxation-time collision operators is proposed by Ginzburg [44]. In Ref. [44], the model works with the multiple-relaxation-time collision operators (MRT), two-relaxationtime collision operators (TRT) and the Bhatnagar–Gross–Krook collision operator (BGK), some types of partial differential equations are discussed. In Ref. [44], the fourth and sixth order analysis of the diffusion term and third-order analysis of the convection term are given. In the frame of the TRT models in 2D and 3D, Ginzburg presents the third order analysis of the convection term and the fourth and sixth order analysis of the diffusion term. The approach is the only possible one, at least beyond the BGK model. The leading-order correction vanishes when,

   4 1 1 1 1 1 þ þ ¼ ; 3 2 ke 2 kD 9 where ke and kD are two eigenvalues of the TRT model. In the case of the BGK, k1e ¼ k1D ¼ s the solution gives s ¼ 12 þ 121=2 . This is similar to C3 = 0 in our model. In this paper, we furnished more higher-order multi-scale time expansion than the preceding work in Ref. [29]. It is difficult to obtain the lattice Boltzmann equation in any time scale. Therefore, we only provide a series of six partial differential equations in this paper. In our former works the equilibrium function for the Navier–Stokes equations up to fourth order moments in the frame of the two and three-dimensional problems [45,46], in the paper, we proposed the moments up to sixth order. Two open issues on the lattice Boltzmann method are considered in this paper. The first issue concerns the approach of how to obtain a higher-order truncation error for a lattice Boltzmann model. If the multi-scale technique and the Chapman– Enskog expansion are used, then how many series of partial differential equations are enough, and how many conservation laws in different time scales are enough? The second issue concerns numerical convergence, numerical stability of lattice Boltzmann scheme and its dependence on the Knudsen number e and relaxation time factor s and how to construct the third-order dispersion term and the fourth-order dissipation term. This paper is organized as follows: in the next section, a series of partial differential equations are given. In Section 3, we propose a lattice Boltzmann model for one-dimensional conservation law equation by using higher-order moment method. In Section 4, some numerical examples are given, and in Section 5 some conclusions are discussed. 2. A series of partial differential equations in different time scales in the lattice Boltzmann model 2.1. The lattice Boltzmann equation Consider a D-dimensional lattice meshes, the particles velocity can be discrete into b directions. In the lattice, each node has b nearest neighbors connected by b links. At each node, it allows rest particles to exist, thus, the particles velocity is discrete into b + 1 velocities actually. The distribution function fa (x, t) is defined as the one-particle distribution function with velocity ea at time t, position x. We assume fa(x, t) possesses the equilibrium distribution faeq ðx; tÞ, and it meets the condition

X

faeq ðx; tÞ ¼

a

X

fa ðx; tÞ:

ð2Þ

a

The distribution function fa(x, t) satisfies the lattice Boltzmann equation

fa ðx þ ea ; t þ 1Þ ¼ fa ðx; tÞ þ Xa ðx; tÞ;

ð3Þ

where s is the single-relaxation time factor, Xa(x,t) is the collision term. In Eq. (3), collision term Xa(x, t) can be expressed as multiple-relaxation-time collision operators, two-relaxation-time collision operators and the Bhatnagar–Gross–Krook collision operator [44]. In those models, the BGK collision operator is the simplest one [2]. The BGK collision operator is

Xa ðx; tÞ ¼ 

 1 fa ðx; tÞ  faeq ðx; tÞ ;

s

ð4Þ

Y. Dong et al. / Applied Mathematical Modelling 34 (2010) 481–494

483

where faeq ðx; tÞ is the local equilibrium distribution function at position x, and time t, with velocity ea, s is the single-relaxation-time. 2.2. A series of partial differential equations in different time scales The Knudsen number e is defined as e ¼ L‘, where ‘is the mean free path, and L is the characteristic length. We assume that the Knudsen number e is equal to the time step Dt[29], thus, the lattice Boltzmann equation (3) can be written as

fa ðx þ eea ; t þ eÞ ¼ fa ðx; tÞ 

1

s

 fa ðx; tÞ  faeq ðx; tÞ :

ð5Þ

Using the Taylor expansion to Eq. (5),

fa ðx þ eea ; t þ eÞ  fa ðx; tÞ ¼

1 X n¼1



en @ n! @t

þ ea

@ @x

n

fa ðx; tÞ:

ð6Þ

Retaining terms up to O(e7), we have

   2  3 @ @ e2 @ @ e3 @ @ fa ðx þ eea ; t þ eÞ  fa ðx; t Þ ¼ e fa ðx; tÞ þ fa ðx; tÞ þ ea fa ðx; tÞ þ þ ea þ ea @t @x @x @x 2 @t 6 @t      6 4 5 e4 @ @ e5 @ @ e6 @ @ þ fa ðx; tÞ þ fa ðx; tÞ þ fa ðx; tÞ þ Oðe7 Þ: þ ea þ ea þ ea @x @x @x 24 @t 120 @t 720 @t ð7Þ Next step is that Chapman–Enskog expansion [47] is applied to fa(x, t) under the assumption of small Knudsen number e, it is 6 X

fa ¼ faeq þ

ei fai þ Oðe7 Þ:

ð8Þ

i¼1

To discuss changes in different time scales, introduced as t0, . . . , t6,

ti ¼ ei t;

i ¼ 0; . . . ; 6

ð9Þ

and 6 X @ @ ei þ Oðe7 Þ: ¼ @t @t i i¼0

ð10Þ

Substituting Eqs. (8)–(10) into Eq. (7),the equation to order e is

1 Dfað0Þ ¼  fað1Þ ;

ð11Þ

s

ð0Þ

@ where fa  faeq , partial differential operator D  @t@0 þ ea @x . The equation to order e2 is

  @ ð0Þ 1 1 fa þ  s D2 fað0Þ ¼  fað2Þ : @t 1 2 s

ð12Þ

The equation to order e3 is



s2  s þ

   1 3 ð0Þ 1 @ ð0Þ @ ð0Þ 1 D fa þ 2  s D f þ f ¼  fað3Þ : 6 2 @t 1 a @t2 a s

ð13Þ

The equation to order e4 is

        2 3 7 1 1 @ ð0Þ 1 @ ð0Þ @ ð0Þ 1 @ ð0Þ 1 s3 þ s2  sþ D4 fað0Þ þ 3 s2  s þ D2 fa þ 2 fa þ fa þ f ¼  fað4Þ : s D s 2 12 24 6 @t1 2 @t2 @t 3 2 s @t21 a ð14Þ 5

The equation to order e is



     5 1 1 3 7 1 @ ð0Þ 1 2 @ ð0Þ D5 fað0Þ þ 4 s3 þ s2  s þ D3 f a þ 3 s2  s þ D f 4 4 120 2 12 24 @t1 6 @t 2 a       1 @ ð0Þ @ ð0Þ 1 @2 1 @2 1 þ2 f þ f þ 3 s2  s þ D 2 fað0Þ þ 2  s fað0Þ ¼  fað5Þ : s D 2 @t 3 a @t 4 a 6 2 @t s @t 1 1 @t 2

s4  2s3 þ s2  s þ

ð15Þ

484

Y. Dong et al. / Applied Mathematical Modelling 34 (2010) 481–494

The equation to order e6 is

    5 13 3 3 2 31 1 5 1 1 @ ð0Þ s5 þ s4  s þ s  sþ D6 fað0Þ þ 5 s4  2s3 þ s2  s þ D4 f 2 6 4 360 720 4 4 120 @t 1 a       3 7 1 @ ð0Þ 1 2 @ ð0Þ 1 @ ð0Þ @ ð0Þ þ 4 s3 þ s2  sþ D3 f þ 3 s2  s þ D f þ2 f þ f s D 2 12 24 @t 2 a 6 @t3 a 2 @t4 a @t 5 a       3 7 1 @2 1 @2 1 @2 þ 6 s3 þ s2  sþ D2 2 fað0Þ þ 6 s2  s þ D fað0Þ þ 2 f ð0Þ s 2 12 24 6 2 @t 1 @t 2 @t1 @t3 a @t 1     2 1 @ 3 ð0Þ 1 @ ð0Þ 1 þ s2  s þ fa þ f ¼  fað6Þ : s 3 6 @t1 2 s @t 22 a

ð16Þ

Eqs. (11)–(16) is so-called a series of partial differential equations in different time scales. It is suitable for one-dimensional, two-dimensional and three-dimensional cases. In Ref. [29], a series of four partial differential equations (11)–(14) are given. We find that Eqs. (11)–(14) are not enough to find higher-order moment. By increasing the number of the series of partial differential equations, we can obtain the model with higher-order truncation error, the third- order dispersion term and the fourth-order dissipation term. Even the fifth-order dispersion term and the sixth-order dissipation term can be obtained. We also find six polynomials of the relaxation time factor s in Eqs. (11)–(16), they are

C 1 ¼ 1; 1 C 2 ¼  s; 2

ð17Þ ð18Þ

1 1 ¼ C 22  ; 6 12 3 7 1 1 s þ ¼ C 32  C 2 ; C 4 ¼ s3 þ s2  2 12 24 6 5 2 1 1 1 1 4 3 ¼ C 42  C 22 þ ; C 5 ¼ s  2s þ s  s þ 4 4 120 4 120 5 13 3 3 2 31 1 1 17 C 6 ¼ s5 þ s4  s þ s  sþ ¼ C 52  C 32 þ C2: 2 6 4 360 720 3 720

C 3 ¼ s2  s þ

ð19Þ ð20Þ ð21Þ ð22Þ

Polynomials ((17)–(22)) are the first six Bernoulli polynomials, which are in full agreement with the results in the literature [48]. They can be used to indicate coefficients of the dispersion term and the dissipation term to the modified conservation law equation. These polynomials (Eqs. (17)–(22)) have the character: when s > 1.0,C2, C4 and C6 are negative numbers, but C3 and C5 are positive numbers. In Fig. 1, we point out the relations between Ci and the relaxation time factor s. 2.3. Definitions for moments of equilibrium distribution functions The macroscopic quantity u(x, t) is defined by

uðx; tÞ ¼

X

fa ðx; tÞ:

ð23Þ

a

Fig. 1. Relations between coefficients Ci ands.

Y. Dong et al. / Applied Mathematical Modelling 34 (2010) 481–494

485

According to Eq. (2), we have

X

faeq ðx; tÞ ¼ uðx; tÞ:

ð24Þ

a

Bases on Chapman–Enskog expansion, thus

X

faðnÞ ðx; tÞ ¼ 0;

n P 1:

ð25Þ

a

Combining these equations and higher-order moments, the conservation laws in different time scales are given in Appendix A. Some moments of the equilibrium distribution function are defined as following:

uðx; tÞ ¼

X a

ð0Þ

mj ðx; tÞ ¼

fað0Þ ðx; tÞ;

X

ð26Þ

fað0Þ ðx; tÞeaj ;

ð27Þ

fað0Þ ðx; tÞeai eaj ;

ð28Þ

fað0Þ ðx; tÞeai eaj eak ;

ð29Þ

a

X

pð0Þ ij ðx; tÞ ¼

a

X

ð0Þ

Pijk ðx; tÞ ¼

a ð0Þ

Q ijkm ðx; tÞ ¼

X

fað0Þ ðx; tÞeai eaj eak eam ;

ð30Þ

a

X

ð0Þ

Rijkmn ðx; tÞ ¼

fað0Þ ðx; tÞeai eaj eak eam ean ;

ð31Þ

fað0Þ ðx; tÞeai eaj eak eam ean eal :

ð32Þ

a

X

ð0Þ

Sijkmnl ðx; tÞ ¼

a

Especially, for one-dimensional model, these moments are defined as

uðx; tÞ ¼

X

fa ðx; tÞ;

ð33Þ

a

m0 ðx; tÞ ¼

X

fa ðx; tÞea ;

ð34Þ

fa ðx; tÞe2a ;

ð35Þ

fa ðx; tÞe3a ;

ð36Þ

a

X

p0 ðx; tÞ ¼

a

X

P0 ðx; tÞ ¼

a

Q 0 ðx; tÞ ¼

X

fa ðx; tÞe4a ;

ð37Þ

fa ðx; tÞe5a ;

ð38Þ

fa ðx; tÞe6a :

ð39Þ

a

R0 ðx; tÞ ¼

X a

S0 ðx; tÞ ¼

X a

For some velocity sets in common use, these equilibrium distribution functions can be solved. In the Appendix B, the equilibrium distribution functions for the 3-bit model, 5-bit model, and 7-bit model are given respectively. According to the definition of these moments, we obtain some conservation laws in different time scales, see Appendix C. 3. The conservation law equation with fourth-order truncation error 3.1. Conservation law equation Let us consider the one-dimensional conservation law equation (1), the flux F(u) is a nonlinear function of u. It means that mass diffusion vanishes [49–52], that is to say, the diffusion coefficient l equals to zero.

@u @FðuÞ @2u þ ¼l 2; @t @x @x

l ¼ 0:

ð40Þ

In order to recover the conservation law equation, the flux F(u) meets

FðuÞ ¼ m0 ¼

X

fað0Þ ea :

ð41Þ

a

Therefore, the conservation law in the time scale t0 is

@u @FðuÞ þ ¼ 0: @t 0 @x

ð42Þ

486

Y. Dong et al. / Applied Mathematical Modelling 34 (2010) 481–494

Eq. (42) can be written as

@u @FðuÞ @FðuÞ @u ¼ ¼ : @t0 @x @u @x

ð43Þ

3.2. Looking for p0 and P0 Summing the Eqs. (11)–(14) over a and making (11) + (12)*e + (13)*e2 + (14)*e3, we obtain

X X X 2 @u @FðuÞ @2u X @f D fað0Þ þ e2 C 3 D3 fað0Þ þ e3 C 2 2 þ 3e3 C 3 D2 a þ e3 C 4 D4 fað0Þ ¼ Oðe4 Þ: þ þ eC 2 @t @x @t @t1 1 a a a a ð0Þ

ð44Þ

From Eq. (43), we have

 2 @FðuÞ @u @FðuÞ @u þ ¼ 0: @u @t0 @u @x

ð45Þ

If p0 and P0 meet the following conditions:

 2 @ p0 @FðuÞ ¼ @u @u

ð46Þ

 3 @P0 @FðuÞ ; ¼ @u @u

ð47Þ

and

then

X

ð48Þ

D3 fað0Þ

ð49Þ

a

X

    @ @u @FðuÞ @ @FðuÞ @ p0 þ ¼ 0; þ þ @t 0 @t 0 @x @x @t 0 @x !     @ 2 @u @FðuÞ @2 @FðuÞ @ p0 @ 2 @ p0 @P0 ¼ 0: ¼ 2 þ þ þ þ2 þ 2 @x @x@t 0 @t 0 @x @x @t0 @x @t 0 @t 0

D2 fað0Þ ¼

a

The Eq. (44) becomes

X @u @FðuÞ @2u X 3 @fa 3e C 3 D2 þ e3 C 4 D4 fað0Þ ¼ Oðe4 Þ: þ þ e3 C 2 2 þ @t @x @t @t 1 1 a a ð0Þ

ð50Þ

3.3. Finding the e3 terms The e3 terms are

2

e3 C 2 @@t2u,

P

1

X @u ¼ C 2 D2 fað0Þ ¼ 0; @t1 a

3

a 3e

C 3 D2

ð0Þ

@fa @t 1

, and

P

ae

3

ð0Þ

C 4 D4 fa . According to Eqs. (12) and (48), we have

@2u ¼ 0: @t21

ð51Þ

Other term is

X

"

  @ 3 X ð0Þ @3 @F @ p0 @3 @ p0 @P 0 D f þ 3 þ þ þ3 a 3 2 2 @x @t 0 @x @t 0 @x @t0 a @t0 @x @t0 ! 0 0 3 @ @P @Q : ¼ e3 C 4 3 þ @x @t 0 @x

e3 C 4 D4 fað0Þ ¼ e3 C 4

a

! þ

@ 3 @P0 @Q 0 þ @x3 @t 0 @x

!#

ð52Þ

If we select that Q0 meets

 4 @Q 0 @FðuÞ ; ¼ @u @u

ð53Þ

then

X a

3

4 ð0Þ

3

e C 4 D fa ¼ C 4 e

@ 3 @P0 @Q 0 þ @x3 @t 0 @x

! ¼ 0:

ð54Þ

Y. Dong et al. / Applied Mathematical Modelling 34 (2010) 481–494

487

The Eq. (50) becomes

@u @FðuÞ þ ¼ Oðe4 Þ: @t @x

ð55Þ

Summing the Eqs. (11)–(16) over a and making (11) + (12)*e + (13)*e2 + (14)*e3 + (15)*e4 + (16)*e5, the modified conservation law equation is

@u @FðuÞ þ ¼ E4 þ E5 þ Oðe6 Þ; @t @x where E4, E5 are

E4 ¼ e4

X a

X

ð56Þ

" C 5 D5 fað0Þ þ 4C 4 D3 "

# @ ð0Þ @ ð0Þ @ ð0Þ @2 @2 fa þ 3C 3 D2 fa þ 2C 2 D fa þ 3C 3 D 2 fað0Þ þ 2C 2 fað0Þ ; @t1 @t 2 @t 3 @t 1 @t 2 @t 1

@ ð0Þ @ ð0Þ @ ð0Þ @ ð0Þ fa þ 4C 4 D3 fa þ 3C 3 D2 fa þ 2C 2 D fa @t @t @t @t 1 2 3 4 a # @2 @2 @2 @3 @2 þ 6C 4 D2 2 fað0Þ þ 6C 3 D fað0Þ þ 2C 2 fað0Þ þ C 3 3 fað0Þ þ C 2 2 fað0Þ : @t 1 @t 2 @t 1 @t3 @t 1 @t 1 @t 2

E5 ¼ e5

ð57Þ

C 6 D6 fað0Þ þ 5C 5 D4

ð58Þ

3.4. The dispersion term According to

E4 ¼ 

P

X

ð0Þ

a Dfa

¼ 0,

P

2 ð0Þ

a D fa

e4 C 5 D5 fað0Þ ¼ e4 C 5

a

¼ 0,

P

3 ð0Þ

a D fa

¼ 0, and

P

!

4 ð0Þ

a D fa

¼ 0, the dispersion term is

@ 4 @Q 0 @R0 : þ @x4 @t0 @x

ð59Þ

If assume that

 5 @R0 @FðuÞ þ k; ¼ @u @u

ð60Þ

then the fifth-order dispersion term is

E4 ¼ 

X

e4 C 5 D5 fað0Þ ¼ e4 C 5

a

@ 4 @Q 0 @R0 þ @x4 @t0 @x

! ¼ e4 kC 5

@5u ; @x5

ð61Þ

where k is a parameter. 3.5. The dissipation term P ð0Þ In the dissipation term, some terms are zero except for the term  a e5 C 6 D6 fa . It is



X

"

5

6 ð0Þ

!#   @5 @u @ 5 @R0 @S0 þ 5 : 5 k þ @x @t 0 @x4 @x @t0 @x

5

e C 6 D fa ¼ e C 6

a

ð62Þ

Therefore



X

e5 C 6 D6 fað0Þ ¼ e5 C 6

a

! @5 @F @R0 @S0 : 5k þ þ @x @t 0 @x5 @x

ð63Þ

If we select that

5k

@F @R0 @S0 @u þ þ ¼l ; @x @t 0 @x @x

ð64Þ

then S0 meets

 6 @S0 @F @F : þ ¼ l þ 6k @u @u @u

ð65Þ

Lastly, we obtain the conservation law equation with the fifth-order dispersion term and the sixth-order dissipation term

@u @FðuÞ @5u @6u þ ¼ e4 kC 5 5  e5 lC 6 6 þ Oðe6 Þ: @t @x @x @x

ð66Þ

488

Y. Dong et al. / Applied Mathematical Modelling 34 (2010) 481–494

In summary, higher-order moments p0, P0, Q0, R0 and S0 meet Eqs. (46), (47), (53), (60) and (65) respectively. When F(u) is a specific function, the higher-order moments above should be solved. 3.6. The Hirt’s heuristic stability conditions 5

6

In Eq. (66), the fifth dispersion is e4 kC 5 @@x5u, the sixth dissipation term is e5 lC 6 @@x6u. According to Hirt’s heuristic stability theory [53,54], coefficients e4kC5 > 0, e5lC6 < 0 are necessary stability conditions. Therefore, the stability of the lattice Boltzmann scheme Eq. (3) is controlled by the condition e4kC5 > 0 and e5lC6 < 0. Select k > 0, l > 0, and s > 1, then e4kC5 > 0 and e5lC6 < 0. 4. Numerical examples Since our purpose is to design the lattice Boltzmann model for higher- order accuracy of truncation error by using higherorder method, we will present the numerical examples in this section which are on linear equations and on nonlinear conservation law equations. All the examples are standard problems widely used in the literature to test various schemes for these sample equations. In these numerical examples, we use ‘‘macroscopic variable” boundary conditions. Three kinds of boundary conditions are: (1) Dirichlet boundary condition, at two ends, u(0, t) = uL, u(M, t) = uL, where M is the lattice size. Then,

fatþ1 ð0; tÞ ¼ faeq ð0; tÞ;

fatþ1 ðM; tÞ ¼ faeq ðM; tÞ;

where faeq ð0; tÞ and faeq ðM; tÞ are equilibrium distribution functions at time t on x = 0 and x = M, calculated by using the relations of equilibrium distribution function and macroscopic variable u. @ @ uð0; tÞ ¼ 0, @x uðM; tÞ ¼ 0. We use u(0, t) = u(Dx, t), u(M, t) = u(M  Dx, t), (2) Neumann boundary condition, at two ends, @x where Dx is the spatial step. Then,

fatþ1 ð0; tÞ ¼ faeq ðDx; tÞ;

fatþ1 ðM; tÞ ¼ faeq ðM  Dx; tÞ;

where faeq ðDx; tÞ and faeq ðM  Dx; tÞ are equilibrium distribution functions at time t on x = Dx and x = M  Dx, calculated by using the relations of equilibrium distribution function and macroscopic variable u. (3) Periodic boundary condition: In this paper, periodic boundary condition is

uð0; tÞ ¼ uðM  2; tÞ; uð1; tÞ ¼ uðM  1; tÞ;

uð2; tÞ ¼ uðM; tÞ;

where x = 0,1, . . . ,M, then,

fatþ1 ð0; tÞ ¼ faeq ðM  2; tÞ;

fatþ1 ðM; tÞ ¼ faeq ð2; tÞ;

where faeq ðM  2; tÞ and faeq ð2; tÞ are equilibrium distribution functions at time t on x = M  2 and x = 2, calculated by using the relations of equilibrium distribution function and macroscopic variable u. Example 1

8 @u @u þ ¼ 0; 06x61 t>0 > < @t @x uðx; 0Þ ¼ 1; 0:35 6 x 6 0:75 otherwise uðx; 0Þ ¼ 0; > : uð0; tÞ ¼ uð1; tÞ ¼ 0; t > 0;

ð67Þ

Example 2

8 @u @u þ @x ¼ 0; 06x61 t>0 > @t > > > ( > < 2 12 ð1  ½ði  50Þ=15 Þ 35 6 i 6 65; where x ¼ i=100 uðx; 0Þ ¼ > > 0 otherwise > > > : uð0; tÞ ¼ uðm; tÞ; t > 0:

ð68Þ

Examples 1 and 2 are one-dimensional linear conservation equations [55]. In Figs. 2a and 3a, we give the comparisons of exact solution and the LBM results. We also present the absolute error curves in Figs. 2b and 3b. These absolute errors are controlled in the scope of (5.0  105, 5.0  105). Especially, we find that this model has higher resolution in the discontinuous regions. Since this model possesses fourth-order accuracy of truncation error, the discontinuity should spread over one grid mesh only.

Y. Dong et al. / Applied Mathematical Modelling 34 (2010) 481–494

489

Fig. 2. Comparison of exact solution and the LBM result of the Example 1. (a) Comparison of exact solution and the LBM result. (b) Error bars of the lattice Boltzmann model. Parameters are: c = 500.0, lattice size M = 100, time t = 10,000Dt = 0.2, k = 0.1, l = 0.001.

Fig. 3. Comparison of exact solution and the LBM result of the Example 2 at time t = 0.4. (a) Comparison of exact solution and the LBM result. (b) Error of the lattice Boltzmann model. Parameters are: c = 500.0, lattice size M = 100, time t = 20,000, Dt = 0.4, k = 0.1, l = 0.001.

Example 3

8 @u þ u @u ¼ 0; 0:5 6 x 6 0:5 > @x  > < @t 1 x < 0:5 uðx; 0Þ ¼ > 0 x P 0:5 > : uð0; tÞ ¼ uð1; tÞ ¼ 0; t > 0:

t>0 ð69Þ

We also apply our model to the one-dimensional nonlinear conservation law equations [52]. In Example 3, we obtain a numerical result with higher resolution (see Fig. 4a). The absolute error curve is plotted in Fig. 4b. We find the error is less than 0.002. The last example is Eq. (70), the smooth 1-periodic initial data is

uðx; 0Þ ¼ SinðpxÞ:

ð70Þ

The periodic boundary condition u(0,t) = u(1,t) is used. In order to obtain the exact numerical result, we use the fourthorder Runge–Kutta method to solve the Eq. (70). It is well-known that solution of Eq. (70) develops a shock discontinuity at t  0.31. In Fig. 5 we give two numerical results in time t = 0.15 and t = 0.4. At time t = 0.15 the shock begins to shapes, and at

490

Y. Dong et al. / Applied Mathematical Modelling 34 (2010) 481–494

Fig. 4. Comparison of exact solution and the LBM result of the Example 1 at time t = 0.2. (a) Comparison of exact solution and the LBM result. (b) Absolute error of the lattice Boltzmann model. Parameters are: c = 500.0, lattice size M = 100, time t = 10,000, Dt = 0.2, k = 0.1, l = 0.001.

Fig. 5. Comparison between exact solution and the LBM result of the Example 4. Solid line is the exact solution, circles are this LBM result. (a) t = 75, Dt = 0.15, (b) t = 200, Dt = 0.4. Parameters are: k = 0.1, l = 0.001, c = 5.0, lattice size M = 100.

Table 1 L1-norm of the errors for Numerical Solutions of the Example 4. In LBM solutions, time steps N = 40, N = 80, N = 160, N = 320 are corresponding to M = 54, M = 108, M = 216 and M = 432 for t = 0.15. Time steps N = 40, N = 80, N = 160, N = 320 are corresponding to M = 20, M = 40, M = 80 and M = 160 for t = 0.4. t = 0.15

t = 0.4

N

LxF

ORD

STG

LBM

LxF

ORD

STG

LBM

The Eq. (44) becomes 40 80 160 320

0.023702 0.12249 0.06246 0.003158

0.02620 0.000667 0.000169 0.000043

0.000859 0.000232 0.000061 0.000016

0.000437 0.0000748 0.000034 0.000005

0.044449 0.023486 0.011393 0.005235

0.003612 0.001291 0.000498 0.000209

0.000849 0.000277 0.000098 0.000038

0.000422 0.0000851 0.000045 0.000003

time t = 0.4 the shock has formed. Table 1 shows the L1-norm of the errors at the pre-shock time t = 0.15 and post-shock time t = 0.4 for some classical numerical methods and this lattice Boltzmann model [52], in which the L1-norm is defined as

491

Y. Dong et al. / Applied Mathematical Modelling 34 (2010) 481–494

Table 2 Convergence rates for Numerical Solutions of the Example 4. In LBM solutions, time steps N1 = 40, N2 = 80, N3 = 160, N4 = 320 are corresponding to M1 = 54, M2 = 108, M3 = 216 and M4 = 432 for t = 0.15. Time steps N1 = 40, N2 = 80, N3 = 160, N4 = 320 are corresponding to M1 = 20, M2 = 40, M3 = 80 and M4 = 160 for t = 0.4. t = 0.15

R1 R2 R3

t = 0.4 LxF

ORD

STG

LBM

LxF

ORD

STG

LBM

1.7  105 1.1  105 1.8  106

1.2  105 6.6  108 4.3  109

3.5  107 2.1  108 1.6  109

1.8  107 6.2  109 1.2  109

4.1  105 5.6  106 6.8  107

3.7  106 3.2  107 3.1  108

8.7  107 7.1  108 6.1  109

4.8  107 1.9  108 3.4  109

P L1 ¼ m i¼1 jEr i j, M is the lattice size. In Table 1, we find this lattice Boltzmann model has higher accuracy and resolution than the first-order Lax–Friedrich scheme (LxF), the second-order non-oscillatory central differencing non-staggered scheme (ORD), and the second-order non-oscillatory central differencing with staggered scheme (STG). Our method has two advantages, which are: (1) the widths of the shock waves are only one grid that are less than the corresponding classical high-resolution model. (2) The absolute errors are in the scope of (104, 104). These are the features we are seeking for. In order to show the convergence effect, we define the convergence rate ðiÞ

Ri ¼

L1 Mi

ðiþ1Þ

L

 M1

iþ1

Niþ1  Ni

;

ð71Þ

ðiÞ

where L1 ,Mi, and Ni denote L1-norm, the lattice size, the time steps for case i, respectively. In the Table 2, we find the LBM scheme possesses higher convergence rate than the LxF, ORD and STG schemes. Example 4

8 @u þ u @u ¼ 0; 06x61 t>0 > @x > < @t uðx; 0Þ ¼ SinðpxÞ 0 6 x 6 1 > > : uð0; tÞ ¼ uð1; tÞ t > 0

ð72Þ

5. Conclusion In 1991, Frisch pointed out the possibility that increases the requirements of higher-order moment may be used to construct lattice gas model for Navier–Stokes equation [56]. According to this idea, we have completed the higher-order moment method for the conservation law equation. The numerical results agree well with the results of classical schemes. Some key points are as follows: (1) We obtain a new series of partial differential equations in the different time scales, this is the fundamental of higherorder moment method. (2) The first six Bernoulli polynomials Eqs. (17)–(22) have been found. These polynomials have the character: when s > 1.0, C2, C4 and C6 are negative numbers, but C3 andC5 are positive numbers. (3) We obtain lattice Boltzmann model for the conservation law equation with the fourth-order accuracy of truncation error. The dispersion term and dissipation term can be controlled to positive and negative. This ideal and method of this paper can be spread into one-dimensional Euler equations and Navier–Stokes equations. Perhaps we could use multi-layer lattice Boltzmann model for the Euler equations and Navier–Stokes equations to yield higher-order accuracy of truncation error. These are our future works. Acknowledgements This work is Project 20092005 supported by Gradate Innovation Fund of Jilin University, the National Nature Science Foundation of China, (Grant No. 10072023, Grant No. 90305013), and the ChuangXin Foundation of Jilin University (No. 2004CX041). We would like to thank Prof. Qian Yuehong, Prof. Wang Jianping, Liu Yanhong, Wang Huimin, Yan Bo, Shi Xiubo, Li Tingting, Chen Zhenlong, Chen Yuanjian, and Yin Xianli for their many helpful suggestions.

492

Y. Dong et al. / Applied Mathematical Modelling 34 (2010) 481–494

Appendix A. In order to obtain the conservation laws in different time scales, we need supply those partial differentials on time ti to the equilibrium distribution function. Eqs. (11)–(16) are changed into the following forms:

  @u X @ ð0Þ ¼ ea fa ; @t 0 @x a i @u X h ¼ C 2 D2 fað0Þ ; @t 1 a   @u X @ ð0Þ ¼ C 3 D3 fað0Þ  2C 2 D fa ; @t 2 @t 1 a " # @u X @ ð0Þ @ 2 ð0Þ 4 ð0Þ 2 @ ð0Þ ¼ C 4 D fa  3C 3 D f  2C 2 D f  C 2 2 fa ; @t 3 @t 1 a @t 2 a @t 1 a " # @u X @ ð0Þ @ 2 ð0Þ @2 5 ð0Þ 3 @ 2 @ ð0Þ ð0Þ ð0Þ ; ¼ C 5 D fa  4C 4 D f  3C 3 D f  2C 2 D f  3C 3 D 2 fa 2C 2 f @t 4 @t 1 a @t2 a @t3 a @t 1 @t 2 a @t 1 a " @u X @ ð0Þ @ ð0Þ @ ð0Þ @ ð0Þ @2 ¼ C 6 D6 fað0Þ  5C 5 D4 fa  4C 4 D3 fa  3C 3 D2 fa  2C 2 D fa  6C 4 D2 2 fað0Þ @t 5 @t @t @t @t @t 1 1 2 3 4 a # 2 2 3 2 @ @ @ @ 6C 3 D f ð0Þ  2C 2 f ð0Þ  C 3 3 fað0Þ  C 2 2 fað0Þ : @t 1 @t 2 a @t 1 @t3 a @t 1 @t 2

ðA:1Þ ðA:2Þ ðA:3Þ ðA:4Þ ðA:5Þ

ðA:6Þ

Appendix B. For the 3-bit model, those velocities are e0 = 0, e1 = c, e2 = c. These equilibrium distribution functions can be solved by combining Eqs. (33)–(35), they are

ð0Þ

1 ðm0 c þ p0 Þ; 2c2 1 ¼ 2 ðm0 c þ p0 Þ; 2c ð0Þ ð0Þ ¼ u  f1  f2 :

f1 ¼ ð0Þ

f2

ð0Þ f0

ðB:1Þ ðB:2Þ ðB:3Þ

For the 5-bit model, those velocities are ea = (0, c, c, 2c, 2c), where, a = 0,1, . . . , 4. These equilibrium distribution functions can be solved by combining Eqs. (33)–(37), they are ð0Þ

1 0 3 4m c þ 4p0 c2  Q 0  P0 c ; 4 6c

1 ¼ 4 4m0 c3 þ 4p0 c2  Q 0 þ P0 c ; 6c

1 0 ¼ 2P c  2m0 c3 þ Q 0  p0 c2 ; 4 24c

1 ¼ 2P 0 c þ 2m0 c3 þ Q 0  p0 c2 ; 24c4 ð0Þ ð0Þ ð0Þ ð0Þ ¼ u  f1  f2  f3  f4 :

f1 ¼ ð0Þ

f2

ð0Þ

f5

ð0Þ

f4

ð0Þ f0

ðB:4Þ ðB:5Þ ðB:6Þ ðB:7Þ ðB:8Þ

For the 7-bit model, those velocities are ea = (0, c, c, 2c, 2c, 3c, 3c), where, a = 0, 1, . . . , 6. These equilibrium distribution functions can be solved by combining Eqs. (33)–(39), they are ð0Þ

1 ð540p0  195Q 0 c2 þ 15S0 þ 540m0 c5  195P0 c3 þ 15R0 cÞ; 720c6 1 ¼ ð540p0  195Q 0 c2 þ 15S0  540m0 c5 þ 195P0 c3  15R0 cÞ; 720c6 1 ¼ ð6S0  12R0 c þ 60Q 0 c2 þ 120P0 c3  54p0 c4  108m0 c5 Þ; 720c6 1 ¼ ð6S0 þ 12R0 c þ 60Q 0 c2  120P0 c3  54p0 c4 þ 108m0 c5 Þ; 720c6 1 ¼ ðS0 þ 3R0 c  5Q 0 c2  15P0 c3 þ 4p0 c4 þ 12m0 c5 Þ; 720c6 1 ¼ ðS0  3R0 c  5Q 0 c2 þ 15P0 c3 þ 4p0 c4  12m0 c5 Þ; 720c6 ð0Þ ð0Þ ð0Þ ð0Þ ð0Þ ð0Þ ¼ u  f1  f2  f3  f4  f5  f6 :

f1 ¼ ð0Þ

f2

ð0Þ

f3

ð0Þ

f4

ð0Þ

f5

ð0Þ

f6

ð0Þ f0

ðB:9Þ ðB:10Þ ðB:11Þ ðB:12Þ ðB:13Þ ðB:14Þ ðB:15Þ

Y. Dong et al. / Applied Mathematical Modelling 34 (2010) 481–494

493

Appendix C. According to the definition of higher-order moments Eqs. (33)–(39), Eqs. (A.1)–(A.6) become the following forms:

@u @m0 þ ¼ 0; @t 0 @x   @u @ @m0 @ p0 þ C2 þ ¼ 0; @t 1 @x @t 0 @x @u @ @ 2 m0 @ 2 p0 @ 2 P 0 þ C3 þ 3 þ 2 @t 2 @x @x@t 0 @x2 @t 20

ðC:1Þ ðC:2Þ ! ¼ 0;

@u @ @ 3 m0 @ 3 p0 @ 3 P0 @3Q 0 þ C4 þ 6 þ 4 þ 3 @t 3 @x @x2 @t 0 @ 3 x3 @t 30 @x@t 20 @u @ þ C5 @t 4 @x þ 4C 4

@ @t1

þ 3C 3

@ @t2

@u @ þ C6 @t 5 @x þ 5C 5

@ @t1

þ 4C 4

@ @t2

þ 6C 4

@2 @t21

ðC:3Þ !

@ @ 2 m0 @ 2 p0 @2u þ Þ þ C 2 2 ¼ 0; ð @x @t 0 @t1 @x@t1 @t 1 ! @ 4 m0 @ 4 p0 @ 4 P0 @4Q 0 @ 4 R0 þ 10 þ 10 2 2 þ 5 3 þ 4 4 3 @x @t @x4 @t 0 @x@t 0 @x @t0 0 ! @ @ 2 m0 @ 2 p0 @ 2 P 0 þ3 þ 2 2 @x @x@t 0 @x2 @t 0  0  @ @m @ p0 @2u þ 2C 2 þ ¼ 0; @x @t 0 @t 1 @t 2 @x ! @ 5 m0 @ 5 p0 @ 5 P0 @5Q 0 @ 5 R0 @ 5 S0 þ 15 þ 20 þ 15 þ 6 þ 5 @x4 @t 0 @x5 @t 50 @x@t 40 @x2 @t 30 @x3 @t20 ! @ @ 3 m0 @ 3 p0 @ 3 P0 @3Q 0 þ6 þ4 2 þ 3 3 3 2 @x @x @t 0 @ x3 @t 0 @x@t 0 !   2 0 2 0 2 0 @ @ m @ p @ P @ @ @m0 @ p0 þ 3C þ 3 þ þ 2 3 @x @t3 @x @t 0 @x@t 0 @x2 @x @t 20  0  2 3 2 0 @ @m @p @ u @ u @ u þ 2C 2 þ þ C 3 3 þ C 2 2 ¼ 0: @x @t 0 @t 1 @t 3 @x @t 1 @t 2 þ 3C 3

ðC:4Þ

ðC:5Þ

ðC:6Þ

References [1] U. Frisch, B. Hasslacher, Y. Pomeau, Lattice gas automata for the Navier–Stokes equations, Phys. Rev. Lett. 56 (1986) 1505–1508. [2] Y.H. Qian, D. d’humieres, P. Lallemand, Lattice BGK model for Navier–Stokes equations, Europhys. Lett. 17 (6) (1992) 479–484. [3] H.D. Chen, S.Y. Chen, M.H. Matthaeus, Recovery of the Navier–Stokes equations using a lattice Boltzmann gas method, Phys. Rev. A 45 (1992) 5339– 5342. [4] R. Benzi, S. Succi, M. Vergassola, The lattice Boltzmann equations: theory and applications, Phys. Rep. 222 (1992) 147–197. [5] S.Y. Chen, G.D. Doolen, Lattice Boltzmann method for fluid flows, Annu. Fluid Mech. 3 (1998) 314–322. [6] L.-S. Luo, Theory of the lattice Boltzmann method: lattice Boltzmann method for nonideal gases, Phys. Rev. E 62 (2000) 4982. [7] X.W. Shan, H.D. Chen, Lattice Boltzmann model of simulating flows with multiple phases and components, Phys. Rev. E 47 (1993) 1815. [8] X.W. Shan, H.D. Chen, Simulation of non-ideal gases liquid–gas phase transitions by the lattice Boltzmann equation, Phys. Rev. E 49 (1994) 2941. [9] M. Swift, W. Osborn, J. Yeomans, Lattice Boltzmann simulation of nonideal fluids, Phys. Rev. Lett. 75 (1995) 830. [10] K. Gustensen, D.H. Rothman, S. Zaleski, et al, Lattice Boltzmann model of immiscible fluids, Phys. Rev. A 43 (1991) 4320–4327. [11] K.N. Premnath, J. Abraham, Three-dimensional multi-relaxation lattice Boltzmann models for multiphase flows, J. Comput. Phys. 224 (2007) 539–559. [12] S.P. Dawson, S.Y. Chen, G. Doolen, lattice Boltzmann computations for reaction–diffusion equations, J. Chem. Phys. 98 (1993) 1514–1523. [13] D.J. Holdych, R.O. Georgiadis Buckius, Magration of a van der Waals bubble: lattice Boltzmann formulation, Phys. Fluids 13 (2001) 817. [14] M.R. Swift, E. Orlandini, W.R. Osborn, et al, Lattice Boltzmann simulations of liquid–gas and binary systems, Phys. Rev. E 54 (1996) 5041. [15] R.S. Maier, R.S. Bernard, D.W. Grunau, Boundary conditions for the lattice Boltzmann method, Phys. Fluids 6 (1996) 1788. [16] S. Succi, E. Foti, F.J. Higuera, 3-Dimensional flows in complex geometries with the lattice Boltzmann method, Europhys. Lett. 10 (1989) 433. [17] A. Ladd, Numerical simulations of particle suspensions via a discretized Boltzmann equation. Part 2. Numerical results, J. Fluids Mech. 271 (1994) 311. [18] O. Filippova, D. Hanel, Lattice Boltzmann simulation of gas-particle flow in filtes, Comput. Fluids 26 (1997) 697–712. [19] G.W. Yan, Y.S. Chen, S.X. Hu, Simple lattice Boltzmann model for simulating flows with shock wave, Phys. Rev. E 59 (1999) 454–459. [20] C.H. Sun, Lattice-Boltzmann model for high speed flows, Phys. Rev. E 58 (1998) 7283–7287. [21] F.J. Alexander, H. Chen, S. Chen, et al, Lattice Boltzmann model for compressible fluids, Phys. Rev. A 46 (1992) 1967–1970. [22] M. De Cicco, S. Succi, G. Balla, Nonlinear stability of compressible thermal lattice BGK model, SIAM J. Sci. Comput. 21 (1999) 366–377. [23] R.J. Mason, A compressible lattice Boltzmann model, Bull. Am. Phys. Soc. 45 (2000) 168–170. [24] R.J. Mason, A multi-speed compressible lattice Boltzmann model, J. Stat. Phys. 107 (2002) 385–400. [25] G.W. Yan, Y.F. Dong, Y.H. Liu, An implicit Lagrangian lattice Boltzmann method for the compressible flows, Int. J. Numer. Methods Fluids 51 (2006) 1407–1418. [26] T. Kataoka, M. Tsutahara, Lattice Boltzmann method for the compressible Euler equations, Phys. Rev. E 69 (5) (2004). Art No. 056702. [27] T. Kataoka, M. Tsutahara, Lattice Boltzmann method for the compressible Navier–Stokes equations with flexible specific-heat ratio, Phys. Rev. E 69 (3) (2004). Art No. 035701. [28] G.W. Yan, J.Y. Zhang, Y.H. Liu, et al, A multi-energy-level lattice Boltzmann model for the compressible Navier–Stokes equations, Int. J. Numer. Methods Fluids 55 (2007) 41–56. [29] G.W. Yan, A lattice Boltzmann equation for waves, J. Comput. Phys. 161 (2000) 61–69.

494 [30] [31] [32] [33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48] [49] [50] [51] [52] [53] [54] [55] [56]

Y. Dong et al. / Applied Mathematical Modelling 34 (2010) 481–494 G.W. Yan, Studies of Burgers equation using a lattice Boltzmann method, Acta Mech. Sinica 31 (1999) 143–151 (in Chinese). G.W. Yan, M. Song, Recovery of the solitons using a lattice Boltzmann model, Chinese Phys. Lett. 16 (1999) 109–110. G.W. Yan, L. Yuan, Lattice Bhatnagar–Gross–Krook model for the Lorenz attractor, Physica D 154 (2001) 43–50. J.G. Zhou, Lattice Boltzmann Methods for Shallow Water Flows, Springer-Verlag, Berlin, Heidelberg, Germany, 2000. I. Ginzburg, Variably saturated flow described with the anisotropic lattice Boltzmann methods, J. Comput. Fluids 25 (2006) 831–848. Z.H. Chai, B.C. Shi, A novel lattice Boltzmann model for the Poisson equation, Appl. Math. Modell. 32 (2008) 2050–2058. S. Succi, Numerical solution of the Schrödinger equation using discrete kinetic theory, Phys. Rev. E 53 (1996) 1969–1975. S. Succi, Lattice quantum mechanics: an application to Bose–Einstein condensation, Int. J. Mod. Phys. C 9 (1998) 1577–1585. S. Palpacelli, S. Succi, Numerical validation of the quantum lattice Boltzmann scheme in two and three dimension, Phys. Rev. E 75 (2007) 066704. S. Palpacelli, S. Succi, R. Spiggler, Ground-state computation of Bose–Einstein condensates by an imaginary-time quantum lattice Boltzmann scheme, Phys. Rev. E 76 (2007) 036712. L.H. Zhong, S.D. Feng, P. Dong, S.T. Gao, Lattice Boltzmann schemes for the nonlinear Schrödinger equation, Phys. Rev. E 74 (2006) 036704. U. Frisch, D. d’Humières, B. Hasslacher, P. Lallemand, Y. Pomeau, J.P. Rivet, Lattice gas hydrodynamics in two and three dimensions, Complex Sys. 1 (1987) 649–707. D. d’Humières, Generalized lattice-Boltzmann equations, AIAA rarefied gas dynamics: theory and simulations, Prog. Astronaut. Aeronaut. 59 (1992) 450–548. G.W. Yan, Y.S. Chen, S.X. Hu, A lattice Boltzmann method for KdV equation, Acta Mech. Sinica 14 (1998) 18–26. I. Ginzburg, Equilibrium-type and link-type lattice Boltzmann model for generic advection and anisotropic-dispersion equation, Adv. Water Resour. 28 (2005) 1171–1195. P. Lallemand, L.S. Luo, Theory of the lattice Boltzmann method: dispersion, dissipation, isotropy, Galilean invariance and stability, Phys. Rev. E 61 (2000) 6546–6562. D. d’Humières, I. Ginzburg, M. Krafczyk, P. Lallemand, L.S. Luo, Multiple-relaxation time lattice Boltzmann models in three-dimension, Phil. Trans. Royal Soc. Lond. A 360 (2002) 437–451. S. Chapman, T.G. Cowling, The Mathematical Theory of Non-uniform Gas, Cambridge University Press, Cambridge, 1970. D.J. Holdych, D.R. Noble, J.G. Georgiadis, et al, Truncation error analysis of lattice Boltzmann methods, J. Comput. Phys. 193 (2004) 595–619. C.W. Shu, Stanley Osher, efficient implementation of essentially non-oscillatory shock capturing schemes, J. Comput. Phys. 77 (1988) 439–471. C.W. Shu, Stanley Osher, efficient implementation of essentially non-oscillatory shock capturing schemes. II., J. Comput. Phys. 83 (1989) 32–78. A. Harten, ENO schemes with subcell resolution, J. Comput. Phys. 83 (1989) 148–184. H. Nessyahu, E. Tadmor, Non-oscillatory central differencing for hyperbolic conservation laws, J. Comput. Phys. 87 (1990) 408–419. C.W. Hirt, Heuristic stability theory for finite-difference equations, J. Comput. Phys. 2 (1968) 339–355. J. Wang, R.X. Liu, A new approach to design high-order schemes, J. Comput. Appl. Math. 134 (2001) 59–67. H. Yang, An artificial compression method for ENO schemes: the slope modification method, J. Comput. Phys. 89 (1990) 125–160. U. Frisch, Relation between the lattice Boltzmann equation and the Navier–Stokes Equation, Physica D 47 (1991) 231–232.