International Conference on Fascinating Advancement in Mechanical Engineering (FAME2008),
11‐13, December 2008

Performance of 1-D and 2-D Lattice
Boltzmann (LB) in Solution of the Shock
Tube Problem
M. Komeili1, M. Mirzaei2, M. Shabouei3

Abstract- In this paper we presented a lattice Boltzmann
with square grid for compressible flow problems. Triple
level velocity is considered for each cell. Migration step
use discrete velocity but continuous parameters are
utilized to calculate density, velocity, and energy. So, we
called this semi-discrete method. To evaluate the
performance of the method the well-known shock tube
problem is solved, using 1-D and 2-D version of the
lattice Boltzmann method. The results of these versions
are compared with each other and with the results of the
analytical solution.
Keywords___ distribution function, Lattice Boltzmann,
phase velocity, shock tube, sound velocity.
I. INRTODUCTION
T is about 20 years since Frisch first time
(1987) succeeded to use the Lattice Boltzmann
Equations to calculate the viscosity of Lattice Gas
Cellular Automata (LGCA) [1]. Then, this method
has been considered as a new solver of the fluid
flow in the various regimes. Ranging from low
speed [2-4] to high speed flows [5-8]. In 1987
Alexander introduced a modified LB method to
model the compressible flows [9]. The Sung model
modifies the LB in such a manner that the fluid
velocity is added to particles’ velocity [10-17]. Its
method regains the results of the Euler and Navierstockes solutions with first-order and second-order
precision, respectively.
In the present work, we intend to develop the Sung
method to model 1-D and 2-D compressible flow
fields and to evaluate the performance of the model
for solution of 1-D problems.

I

1.
2.
3.

M. Komeili is the Master student of Aerospace
Engineering in K.N. Toosi University, Tehran, Iran
(e-mail: matin_444@yahoo.com).
M. Mirzaei is Associate professor in K.N.Toosi
University, Tehran, Iran (email: mirzaei@kntu.ac.ir).
M. Sabouei is the Master student of Aerospace
Engineering in K.N. Toosi University, Tehran, Iran
(email:m.shabouei@gmail.com).

II. LATTICE BOLTZMANN EQUATIONS
The first step in Lattice Boltzmann models is
determination of a distribution function. It defines
the probability of finding a particle on the specific
position and the specific time with a definite
velocity.
While, the distribution function is identified
throughout the domain. We will be able to calculate
the Microscopic quantities (i.e. density, velocity,
pressure, temperature, and energy).
The relations between the microscopic quantities
and the distribution function are

   fdv

(1)

V   vfdv

(2)

E   efdv

(3)

Where, is the distribution function, the lower
cases relate to particles, and upper cases relate to
flow.
So, if we can find the distribution function in
each time step our problem has been solved. The
distribution function in the next time step is
calculated from the following equation:
f  x  vt , v, t  t   f  x , v, t    x, v, t 

(4)

Department of Mechanical Engineering, Mepco Schlenk Engineering College, Sivakasi, India
Equation (3) is called BGK equation Where,
the collision operator that can be calculated by:
 x, v, t   

1

 

f x, v, t   f

eq 

x, v, t 



is

If we equate relaxation time (  )
collision operator converts to:

(5)

f is equilibrium distribution function. This
quantity is only a function of microscopic
quantities. In the other word, the distribution
function and consequently the microscopic values
are obtained in the subsequent time if
f is
identified.

III. MODIFIED LATTICE BOLTZMANN EQUATION
A. Distribution Function
Conventional Lattice Boltzmann models cannot
be used for compressible flows due to their
limitation in prediction of maximum particles’
velocities. In our model in compressible flows, we
added flow velocity to particles velocity on each
node. So, each particle is translated from a node to
any node into the field.
In this case, distribution function depends on the
following variables:

x ; The position of particle.

r ; The migrating velocity which translates
particles from the starting points to their
destinations.

 ; The particle velocity.
 ; The particle energy.
t ; The time .
Since is a discontinuous variable, but and
are the continuous quantities, we called this method
“semi-discrete”.
In fact each node carries mass (m), momentum

(m  ), and energy (m  ) to its destination.
The relation between microscopic quantities and
the distribution function for this case is as follow:






  
Y    f x , r ,  ,  , t d d
D

r s

(6)


T
Y   , v , E 


  m, m , m

(7)



T

to 1 the

(9)
 x, v, t    f x, v, t   f



eq 

x, v, t 



So, the well-known BGK will be:

 (10)
   
eq  
f x  r t , r ,  ,  , t  t  f x , r ,  ,  , t









B. Equilibrium Distribution Function
First, we introduce some quantity in order to
define equilibrium distribution function


 
c jvk  vk  c jv

(11)


 
c jvk  v  c jv

(12)


 
v k  v  v k

(13)

1
2

 jv  c jv2   

(14)





   2
1 2
v  2c jv .v  c jv  
2





In the above equations c jv ; j  1,2,..., bv is a
set which its arrays are vector, bv defines number
of direction on each node and depends on the type
of the lattice, index v clarifies number of velocity
level on each direction.

 jv is

any kind of energy

except kinetic one and it define by :





  1    1 e
2


D

Where




r vectors.

D  D1  D2 which   D1 and   D2 .
‘s’ is a set which contains all the

D

is




lattice

energy  e  E 

dimensions;

(15)

e

is

interior

1 2
v .
2 

cjv has symmetric properties. This means that:
(8)

International Conference on Fascinating Advancement in Mechanical Engineering (FAME2008),
11‐13, December 2008

bv

 c

jv

j 1

0

 c c
j 1

jv



bv

 c c c
jv

j 1

(20)



eq  
f  x , t    f x , c jvk ,  ,  , t d d  dvk
eq
jvk

bv

jv

(16)

jv

jv

bv
cv I
D

(17)

0

(18)



D



We define  jv as:
(21)

Where, is unit tensor.

Consider x an arbitrary node; v is flow velocity


on this node; vk is the vector which join x to the

nodes that surround; vk are vectors which connect

the end point of v
to the surrounded nodes
(Fig.1).




 jv  m, mc jv , m



T

jv



1 

2
 m, mv  c jv , m  v  c jv     
2



eq
f jvk

Substitution of

into equation (5), the

microscopic quantities are obtained from:

Y   f
D

r S

eq






  
x , r ,  ,  , t d d 

d

(22)

 jv

vk

k ,v , j

Where d vk is determined as:

d vk   k d v

(23)

k


Fig.1. migrating velocity of particles in hexagonal lattice.

k

Equilibrium distribution function is defined:

We have chosen triple velocity level for each
direction (cell), which are d 0 , d1 , and d 2 . d 0 is an



(19)



 
 
dvk   cjv    jv  r  c jvk

eq  
f x, r , , , t  
 
r  c jvk
0





Where

  x  has following properties

  x   0 for x  0







  
f eq x , r ,  ,  , t



is identified on

 
r ,  ,   S  D1  D2 . However, this quantity is

 

nonzero if r ,  ,   c jvk  c jv   jv 



arbitrary portion of the total density, we choose a
value between 0.4 and 0.55 .
d1 and d2 are calculated from the following
relations:

c2 2   b0 d 0   D  1e
d1 
b1 c2 2  c1 2

(24)

D  1e  c1 2   b0 d 0 
b1 c 2 2  c1 2

(25)



d2 

 g x  x dx  g 0
Notice

is density fraction and it equals to



.
eq

Under these circumstances f jvk will be







c1 and c2 are particles module velocity which
should be defined such that d1 and d2 have nonnegative values.

Department of Mechanical Engineering, Mepco Schlenk Engineering College, Sivakasi, India

 D  1e
c1  int 
   b0 d 0
c2  c1  1












(26)

Pressure is obtained from the following equation

P    1e

(27)

IV. TYPES OF LATTICE
A. 1-D Lattice
Fig.2 shows a 1-D lattice. In this lattice, a
particle is able to move only from a node to its
neighbors located in its right side or left- side.

Fig.3.migrating velocity in square lattice

In this figure the migrating velocities cjv are:

  c21
  c31
  c41
 1
c11
  c22
  c32
  c42
 2
c12

Fig.2.migrating velocity in 1‐D lattice

In this figure the migrating velocities are:

Notice that in this case we use dual velocity level
and assume d 0  0 . So there are two levels of

 0
c10
  c22
  2 , c11
  c21
 1
c12

velocity and four directions:

k  1 to k  4

Since we have three levels of velocity and two
directions, the subscripts are defined as:

k  1 to k  2

j  1 to

In addition b1  b2  4
The contribution of each node which surrounds
the end point of vector
(Fig. 4) to the total
is defined as:
density
 1   u 3 v 3

Moreover,

b0  1 , b1  b2  2
The contribution of each node to the total density
is calculated as:

 2   v1

j4

v  1 to v  2

j2

v  0 to v  2

1   v2

j  1 to

 2   u 4 v 4

(29)

 3   u 1v1
 4   u 2 v 2

(28)
3

2

4

1

Notice that in this case the velocity has just one
direction.
B. Square Lattice
If we employ the square lattice (Fig.3) the
particles may move in a 2-D space.

Fig.4.square lattice

International Conference on Fascinating Advancement in Mechanical Engineering (FAME2008),
11‐13, December 2008
V. NUMERICAL SIMULATIONS
A. Interior Nodes
Assuming, leads omission of collision step so:

 

f i out  x  c t , t  t   f i in  x , t 

(30)

D. Results

First, initial values are set for microscopic flow
parameters i.e. Velocity, density, and energy, then
calculate translated mass, momentum, and energy
from each node (based on velocity level and
direction) by the following equations:

m jvk  1

(31)



 jvk  c jv



1 2
2

(32)


2



 jvk  v  2cjv c 

shock tube problem. This lattice sizes guarantee
the grid independency of the results. The results of
1-D and 2-D solutions were compared with each
other and with those of the analytical solution [18]
in Fig.5.

(33)

By multiplying The above quantities by d jv and
summation over all j , v, and k and for each node,
the microscopic quantities in the next step are
obtained.
B. Boundary Condition
For 1-D lattice, there is no boundary condition.
We just add two additional points to both right side
and left side of the domain.
For 2-D lattice periodic boundary condition is
employed in lateral direction in addition for
longitudinal direction two additional columns are
use on the limits of the domain.

The results of the shock tube problem [19,20]
are shown in figures 5. Fig 5-a represents the
variation of the gas density along the tube. In this
figure the results of the 1-D and the 2-D solutions
are compared with that of the analytical solution. It
can be seen that the 1-D and 2-D solutions leads to
the same results and in good agreement with
analytical solution. The same trend can be seen for
variation if internal energy, velocity and pressure in
Figs 5-b to 5-d. notice that the CPU time for 2-D
solution is 15 times of that the 1-D solution.
I.

CONCLUSION

We modified the lattice Boltzmann method with
adding the flow velocity to the migrating velocities.
In this manner the particles can translate to any
node in the domain. Accordingly, high Mach
numbers flow can be easily simulated. The 1-D
lattice is introduced in this paper for the first time.
The shock tube problem and other 1-D flows were
previously solved using 2-D lattice methods ([14])
but for the first time we develop the 1-D lattice
method for such problems. The results of this
method are accurate and moreover, its main
advantage over 2-D lattice method is its
computational efficiency. For instance its CPU time

C. Shock tube
To evaluate the validity of our models we
selected Sod test [12]. The initial conditions are:

 L  1 .0 ;

 R  0.125

p L  1 .0 ;

p R  0 .1

u L  1 .0 ;

u R  0.125

Index L and R show the initial values on the left
side and the right side respectively.
A 400 node lattice is adopted for 1-D sock tube
problem and a 400×4 node lattice is use for 2-D

Fig.5. The shock‐tube problem: the profiles of density,
internal energy, velocity, and pressure at the 75th
iteration time. Analytical solution is based on [18]

Department of Mechanical Engineering, Mepco Schlenk Engineering College, Sivakasi, India

For shock tube problem was 15 times less than
that of the 2-D lattice. The present method may be
expanded for more complicated flows (e.g.
employing of different gases on the two sides of the
tube, consideration of chemical interaction, and
reflection of shock on walls).

REFERENCES
[1]
[2]
[3]
[4]
[5]

[6]

[7]

[8]
[9]
[10]
[11]
[12]

[13]

[14]
[15]

[16]

Frisch, U., Hasslacher, B., and Pomeau, Y.,
“Lattice Gas Automata for the Navier-Stokes
Equation” Phys. Rev.Lett., 56:1505, 1986.
Buick, M., “Lattice Boltzmann Method in
Interfacial Wave Modeling”, PHD thesis at
University of Edinburg, 1990.
Higuera, F., Jemenez, J., And Succi., S., “Lattice Gas
Dynamics with Enhanced Collision”, Europhys. Lett,
9:345, 1989.
Alexander, W., “Theory and Applications of the
Lattice Boltzmann Method”, PHD thesis at university
of Oxford, 1997.
Palmer, B.J., Rector, D.R., “Lattice Boltzmann
Algorithm for Simulating Thermal Flow in
Compressible Fluids”, Journal of Computational
Physics, 161, 11 1-20, 2000.
Vahala, G., Pavlo, P., Vahala, L., Martys, N.S.,
“Thermal Lattice Boltzmann Model(TLBM) for
compressible flow”, International Journal of Modern
Physics C, Vol. 9, No. 8, pp 1247-1261, 1998.
Hinton, F.L., Rosenbluth, M.N., Wong, S.K., LinLio, Y.R., Miller, R.L., “Modified Lattice Boltzmann
method for compressible fluid simulations”, Physical
Review E, Vol. 63, 061212, 2001.
Yang, T., “The Development of a Compressible
Lattice Boltzmann Method for Turbo Machine
Application”, PHD theses at Purdue university.
Alexander, F.J., Chen, H., Doolen, G.D., “Lattice
Boltzmann Model for Compressible Fluids”, Physical
Review A, Vol. 46, No. 4, pp 1967-1970, 1991.
Sun, Ch., “Adaptive Lattice Boltzmann Model for
Compressible Flows”, Tsinghue science and
technology, Vol.5, No.1, 2000.
Sun, Ch., Hsu, A., “Three-Dimensional Lattice
Boltzmann Model for Compressible flows”, Phys.
Rev. E, Vol. 68, 2003.
Sun, Ch., Munn, L.L.,” Influence of erythrocyte
aggregation on leukocyte margination in post
capillary expansions: A lattice Boltzmann analysis”,
Physics A, 362, pp191–196, 2006.
Sun, Ch., Munn, L.L., “Particulate nature of blood
determines macroscopic theology: A 2-D lattice
Boltzmann analysis”, Biophysical Journal, Vol.
88,pp. 1635–1645,2005.
Sun, Ch., “Lattice-Boltzmann models for high speed
flows”, Phys. Rev. E, Vol. 58,No.6, pp. 7283-7287,
1998.
Sun, Ch., “Adaptive lattice Boltzmann model for
compressible flows: Viscous and conductive
properties”, Phys. Rev. E, Vol. 61, No.3, pp. 26452653, 1999.
Moradi, B., “Simulation of Compressible Fluid
Flows with Lattice Boltzmann Method”. M.sc thesis
at Isfahan University of technology, Iran, 2004.

[17] Shiravi.Kh., V.R., “ Parallel Simulation of
Compressible Fluid with Lattice Boltzmann
Method”. M.sc thesis at Isfahan University of
Technology, Iran, 2006.
[18] G. A. Sod, J. Comput. Phys. 27, 1 (1978).
[19] Shabouei M, Ebrahimi R, Mazaheri Body K.,
“Numerical Simulation of Converging Spherical
Shock Problem”, In Ninth International Congress of
Fluid Dynamics & Propulsion (ICFDP 9),
Alexandria, Egypt 2008.
[20] Shabouei M, Ebrahimi R, Mazaheri Body K.,
“Numerical Solution of Cylindrically Converging
Shock Waves”, In International Conference on
Fascinating
Advancement
in
Mechanical
Engineering (FAME2008), India, 2008, arXiv
preprint arXiv:1602.02680.

International Conference on Fascinating Advancement in Mechanical Engineering (FAME2008),
11‐13, December 2008

