Introduction to Lattice Boltzmann Method

In this module, we are going to learn how to solve diffusion and advection-diffusion (Navier-Stokes problem) in Fluid mechanic by Lattice Boltzmann Method (LBM). As we know there are two main approaches in simulating the transport equations continuum and discrete. As a continuum approach we can refer to the Finite difference method (FDM) which we learned in the numerical class. As we experienced in this class there are some criteria which should be met in each types of problems because of (nonlinearity, complex boundary conditions, complex geometry, etc.). If you remember following procedures in continuum approach are expected:

1- First the governing equations are identified (mainly partial differential equations).

2- The next step is to discretize the domain into volume, girds, or elements depending on the method of solution (FVM, FDM or FEM).

3- The scale is macroscopic. The velocity, pressure, temperature of all those particles represented by a nodal value, or averaged over a finite volume, or simply assumed linearly or bi-linearly varied from one node to another. Each volume or node or element contains a collection of particles (huge number, order of 10${}^{16}$).

In discrete approach (microscopic approach) the medium can be considered made of small particles (atom, molecule) and these particles collide with each other. This scale is micro-scale. Hence, we need to identify the inter-particle (inter-molecular) forces and solve ordinary differential equation of Newton's second law (momentum conservation).

The procedures are as follow:

1- At each time step, we need to identify location and velocity of each particle, i.e., trajectory of the particles.

2- Temperature and pressure are related to the kinetic energy of the particles (mass and velocity) and frequency of particles bombardment on the boundaries. This method is called molecular dynamics (MD) simulations.

3- The simulation time step should be less than the particles collision time, which is in the order of seconds (10${}^{-12}$ s). Hence, it is impossible to solve large size problems (order of cm) by MD method.

4- Statistical mechanics need to be used as a translator between the molecular world and the macroscopic world.

5- It is important to specify the appropriate inter-particle force function.

6- The main drawback or obstacle of using MD in simulation of a relatively large system is the computer resource and data reduction process, which will not be available for us in the seen future.

What about a middle man, sitting at the middle of both mentioned techniques, the lattice Boltzmann method?

The main idea of Boltzmann is to bridge the gap between micro-scale and macro-scale by not considering each particle behavior alone but behavior of a collection of particles as a unit. The concepts of LBM are based on Ludwig Boltzmann's kinetic theory of gases. Therefore at first it is strongly recommended to study kinetic theory in detail. The fundamental idea is that gases or fluids can be imagined as consisting of a large number of small particles moving with random motions. Here we only mention to some basic parts of LBM (because of our limitation) .It is highly recommended to review the below mentioned references and also the history of LBM.

LBM comprises two basic steps called collision and propagation of fictitious particles. The exchange of momentum and energy is achieved through particle streaming and billiard-like particle collision. Considering $r$ as the position of a particle relative to the fixed reference node and lets us at first introduce particle distribution function $f$ defined as the number of molecules at time $t$ positioned between $r$and $r+dr$ which have velocities between $c$ and $c+dc$.

The LBM simplifies Boltzmann's original idea of gas dynamics by reducing the number of particles and confining them to the nodes of a lattice. For a two dimensional model, a particle is restricted to stream in a possible of nine directions, including the one staying at rest. These velocities are referred to as the microscopic velocities and denoted by $e_{i} $ which i=0... 8 are known D2Q9 it means two dimensional involved nine velocity directions. For one dimensional the particles are denoted as D1Q3 means one dimensional with two particle placed in the left and right hand side of a certain particle zero. There are also some other consideration however in this module we limit our codes based on these two types of arrangement for 2D and 1D problems . It is necessary to note that in LBM every lattice space and also time step is considered as one.

Main Equations and steps in LBM

As discussed earlier the exchange of momentum and energy is achieved through particle streaming and particle collision. This process can be modelled by the Boltzmann transport equation, which is:

$ \frac{\partial f}{\partial t} +\mathop{u}\limits^{\to } .\nabla f=\Omega $

Microscopic lattice nodes are defined as follow:

D2Q9:

$$\mathop{e_{i} }\limits^{\to } \left\{\begin{array}{l} {(0,0)\, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, i=0\, \, \, \, \, \, \, \, \, \, \, \, \, \, } \\ {(1,0),(0,1),(-1,0),(0,-1)\, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, i=1,2,3,4} \\ {(1,1),(-1,1),(-1,-1),(1,-1)\, \, \, \, \, \, \, \, \, \, \, i=5,6,7,8\, \, } \end{array}\right. $$

D2Q5:

$$\mathop{e_{i} }\limits^{\to } \left\{\begin{array}{l} {(0,0)\, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, i=0\, \, \, \, \, \, \, \, \, \, \, \, \, \, } \\ {(1,0),(0,1),(-1,0),(0,-1)\, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, i=1,2,3,4} \end{array}\right. $$

D1Q3:

$$\mathop{e_{i} }\limits^{\to } \left\{\begin{array}{l} {(0,0)\, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, i=0\, \, \, \, \, \, \, \, \, \, \, \, \, \, } \\ {(1,0),(0,-1)\, \, \, \, \, \, \, \, \, \, \, \, \, \, \, \, i=1,2} \end{array}\right. $$

For each particle on the lattice, we associate a discrete probability distribution function $f$ which describes the probability of streaming in one particular direction. Therefore the macroscopic fluid density can be defined as a summation of microscopic particle distribution function:

$$\rho (\mathop{x}\limits^{\to } ,t)=\sum _{i=0}^{i=8}f_{i} (\mathop{x}\limits^{\to } ,t)$$

Accordingly, the macroscopic velocity $u(\mathop{x}\limits^{\to } ,t)$is an average of microscopic velocities weighted by the distribution functions $f$ :

$$u(\mathop{x}\limits^{\to } ,t)=\frac{1}{\rho } \sum _{i=0}^{i=8}cf_{i} (\mathop{x}\limits^{\to } ,t)\, \mathop{e_{i} }\limits^{\to } $$

The key steps in LBM are the streaming and collision processes which are given by:

$$\underbrace{\, \, \, f_{i} (\mathop{\mathop{x}\limits^{\to } +c\mathop{e_{i} }\limits^{\to } }\limits^{} \, \Delta t,t+\Delta t)-f_{i} (\mathop{\mathop{x}\limits^{\to } ,}\limits^{} \, t)\, \, \, \, }_{Strea\min g} =-\underbrace{\, \frac{[f_{i} (\mathop{\mathop{x}\limits^{\to } ,}\limits^{} \, t)-f_{i}^{eq} (\mathop{\mathop{x}\limits^{\to } ,}\limits^{} \, t)]}{\tau } \, \, }_{Collision} $$

In the actual implementation of the model, streaming and collision are computed separately, and special attention is given to these when dealing with boundary lattice nodes. In the collision term $f_{i}^{eq} (\mathop{x}\limits^{\to } ,t)$is the equilibrium distribution, and $\tau $is considered as the relaxation time towards local equilibrium. For simulating single phase flows, it suffices to use Bhatnagar-Gross-Krook (BGK) collision, whose equilibrium distribution $f_{i}^{eq} $ is defined by: $f_{i}^{eq} (\mathop{x}\limits^{\to } ,t)=\omega _{i} \, \rho +\rho \, s_{i} (\mathop{u}\limits^{\to } (\mathop{x}\limits^{\to } ,t))$ $$s_{i} (\mathop{\mathop{u}\limits^{\to } )=}\limits^{} \omega _{i} \, \left[3\frac{\mathop{e_{i} }\limits^{\to } .\mathop{u}\limits^{\to } }{c} +\frac{9}{2} \frac{(\mathop{e_{i} }\limits^{\to } .\mathop{u}\limits^{\to } )^{2} }{c^{2} } -\frac{3}{2} \frac{\mathop{u}\limits^{\to } .\mathop{u}\limits^{\to } }{c^{2} } \right]$$ $$\omega _{i} =\left\{\begin{array}{l} {4/9\, \, \, \, \, \, \, \, \, \, i=0} \\ {1/9\, \, \, \, \, \, \, \, \, \, \, i=1,2,3,4} \\ {1/36\, \, \, \, \, \, \, \, i=5,6,7,8} \end{array}\right. \, \, $$ $$c=\frac{\Delta x}{\Delta t} $$ Coefficient $c$ is the lattice speed. The fluid kinematic viscosity (or diffusion coefficient in the problem) in the D2Q9 model is related to the relaxation time by: $$\nu =\frac{2\tau -1}{6} \frac{\Delta x^{2} }{\Delta t} $$

The algorithm can be summarized as follows:

1- Initializing distribution function

2- Collision step: calculate the updated distribution function

3- Streaming step

4- Boundary condition

5- Obtaining velocity, density and other macroscopic properties based on distribution function

6- Repeat again until the satisfaction of criterion of steady state

we need to notice that relaxation time should be lower than 0.5 to satisfy low Mach number in LBM else instability problems arise during solution and also during the streaming and collision step, the boundary nodes require some special treatments on the distribution functions in order to satisfy the imposed macroscopic boundary conditions .

Boundary Condition

In these codes we used three types of boundaries (No-slip, Dirichlet and Neumann). For other boundaries it is recommended to read the book or related papers mentioned below.

No-slip boundary condition (Bounce-back BCs):

By the so-called bounce-back we mean that when a fluid particle (discrete distribution function) reaches a boundary node, the particle will scatter back to the fluid along with its incoming direction. IN our codes we use bounce back type which is on the grid and it is first order but we can also use mid-grid bounce back with second order accuracy. There are also different types of bounce back to guarantee no-slip condition. These types are important where we deal with immersed body such as a circle in our domain which needs special treatment to not allow fluid penetrates solid boundary.

First-order Zou-He BCs (Bounce back model)

We can also use this method for boundary with movement such as Lid-driven which will be discussed later. Please take attention how we derive these equations.

$$\begin{array}{l} {f_{1} +f_{5} +f_{8} =\rho -(f_{0} +f_{2} +f_{4} +f_{3} +f_{6} +f_{7} )} \\ {f_{1} +f_{5} +f_{8} =\rho \, u+(f_{3} +f_{6} +f_{7} )} \\ {f_{5} -f_{8} =\rho \, v-f_{2} +f_{4} -f_{6} +f_{7} } \\ {} \end{array}$$

Now we can derive density as: $$\rho =\frac{1}{1-u} \, \left[f_{0} +f_{2} +f_{4} +2(f_{3} +f_{6} +f_{7} )\right]$$

But we know that:

$$f_{1} -f_{1}^{eq} =f_{3} -f_{3}^{eq} $$

Therefore:

$$\begin{array}{l} {f_{1} =\frac{2}{3} \rho \, v+f_{3} } \\ {f_{5} =\frac{1}{6} \rho \, u+\frac{1}{2} \rho \, v-\frac{1}{2} (f_{2} -f_{4} )+f_{7} } \\ {f_{8} =\frac{1}{6} \rho \, u-\frac{1}{2} \rho \, v+\frac{1}{2} (f_{2} -f_{4} )+f_{6} } \\ {} \end{array}$$

Other Boundaries will be discussed while solving the corresponding problems.

References :

1- A.A.Mohamad , Lattice Boltzmann Method , Fundumental and Engineering Applications with Computer Codes ,Springer, DOI 10.1007/978-0-85729-455-5

2- Yuanxun Bill Bao & Justin Meskas ,Lattice Boltzmann Method for Fluid Simulations April 14,2011Lattice

3- Timm Kr¨uger,LBM Workshop (Edmonton, Canada), August 22–26, 2011