Lattice Boltzmann Method-1

Kinetic Theory and micro-macro relations

Pressure and Temperature are related to kinetic energy of particle as1,

\begin {equation}
p = \frac{2}{3}KE\cdot n
\end {equation}
\begin {equation}
T = \frac{2}{3} KE\cdot \frac{1}{k_{B}}
\end {equation}
Phase space and probability distribution function

A coordinate system comprising position and velocity can be used to express the dynamics of a system. Such a system is called phase space. In general, the dynamics of a point in the phase space is defined by three positions (say x, y, z) and three velocity components (u, v, w) or three momenta.

The number of particles that can be found at a location between x and x+dx , y and y+dy, and z and z+dz that have velocities between V_{x} and V_{x}+dV_{x} , between V_{y} and V_{y}+dV_{y} , and between V_{z} and V_{z}+dV_{z} is d^{6} N.

d^{6} N (r, v) = f (r, v)d^{3} r d^{3} v = f (r, v){x} .dV_{y}. dV_{z}

The density of particles in a box of size dx dy dz can be found by integrating f over all possible velocity ranges, i.e.,

n(r) = \int \int \int f (r, v)d^{3}r

Let us define a normalized distribution function as:

\overline{f} (r, v) = \frac{f (r, v)}{n(r)} \text{ }\text{ }\text{ }\text{where} \int \int \int \overline{f} (r, v)d^{3}r = 1

Hence the average value of a quantity Q in an infinitesimal volume dx dy dz (usually referred to as d^{3}r ) is

\begin {equation}
\overline{Q} =  \int \int \int Q(v)\overline{f} (r, v)d^{3}v
\end {equation}

The above integration is over all possible velocities.

Boltzmann Transport Equation

f (r, c, t) gives the number of molecules at time t positioned between r and r + dr that have velocities between c and c + dc. An external force F that acts on a gas molecule of unit mass will change the velocity of the molecule from c to c + Fdt and its position from r to r + cdt.

if no collisions take place between the molecules, The number of molecules f (r, c, t) before the external force is applied is equal to the number of molecules after the disturbance, f (r + cdt, c + Fdt, t + dt)

f (r + cdt, c + Fdt, t + dt)dr dc − f (r, c, t)dr dc = 0

However, if collisions take place between the molecules, there will be a net difference between the numbers of molecules in the interval dr dc. The rate of change between the final status and the initial status of the distribution function is called the collision operator, . Hence the equation for the evolution of the number of molecules can be written as

f (r + cdt, c + Fdt, t + dt)dr dc − f (r, c, t)dr dc = \Omega ( f )dr dcdt 


\begin {equation}
\frac{\partial f}{\partial t}+c \cdot \frac{\partial f}{\partial r}+\frac{F}{m} \cdot\frac{\partial f}{\partial c} = \Omega
\end {equation}

in the absence of external force, we get force free Boltzmann equation

\begin {equation}
\frac{\partial f}{\partial t}+c \cdot \nabla f = \Omega
\end {equation}
BGK model

It is difficult to solve the Boltzmann equation, because the collision term is very complicated. The outcome of two body collisions is not likely to influence significantly the values of many measured quantities (Cercignani 1990). Hence it is possible to approximate the collision operator with a simple operator without introducing a significant error in the outcome of the solution. In 1954, Bhatnagar, Gross, and Krook (BGK) introduced a simplified model for the collision operator

\begin {equation}
\Omega = \frac{ f^{eq} − f }{\tau}
\end {equation}

τ is relaxation factor. The local equilibrium distribution function is denoted by f^{eq} , which is the Maxwell–Boltzmann distribution function:

\begin {equation}
f(c) = 4\pi c^{2}\Big( \frac{m}{2\pi kT}\Big)^{\frac{3}{2}} e^{-\frac{mc^{2}}{2kT}}
\end {equation}

this function increases parabolically from zero for low speeds (c), reaches a maximum value,as c increases up to a certain value, and then decreases exponentially. As the temperature increases, the position of the maximum of f shifts to the right. The total area under the curve is always one, by definition.

\begin {equation}
\frac{\partial f}{\partial t}+c \cdot \nabla f = \frac{ f^{eq} − f }{\tau}
\end {equation}

In the lattice Boltzmann method, the above equation is discretized and assumed to be valid along specific directions, called linkages. Hence the discrete Boltzmann equation can be written along a specified direction as

\begin {equation}
\frac{\partial f_{i}}{\partial t}+c_{i} \cdot \nabla f_{i} = \frac{ f^{eq}_{i} − f_{i} }{\tau}
\end {equation}

The equation is a linear partial differential equation with left-hand side representing advection (streaming) and right-hand-side representing the collision process, the source term.

Lattice Arrangements
1D lattice arrangements
D3Q19: 6 face centers, 12 edge centers and at center. The weighting factors are: for f_{0} it is 12/36, for f_{1} to f_{6} it is 2/36, and for f_{7} to f_{18} it is 1/36.

D3Q15: 6 face centers, 8 corners and one at center. D3Q27: 6 face centers, 8 corners, 12 edge centers and one at center

Weighting factors
  • The distribution function residing on the site has more weight than other distribution functions.
  • A distribution function streaming a short distance should have more weight than a function streaming to a longer distance.
  • The distribution functions streaming the same distance, regardless of their direction, should have the same weight.
  • Also, the sum of all weighting factors should be 1

noting that k = \frac{R}{N}, Maxwell distribution for a gas moving with bulk velocity u can be written as,

\begin {equation}
f^{eq} = \frac{\rho}{(2\pi RT)^{\frac{D}{2}}} e^{-\frac{(c-u)^{2}}{2RT}}
\end {equation}

expanding e^{x} = 1+x+\frac{x^{2}}{2!}+\frac{x^{3}}{3!}+... and ignoring terms of the order of u^{3} to get second order accurate estimate

\begin {equation}
f^{eq}_{i} = \rho W_{i}  \Big[1+\frac{2c \cdot u -u \cdot u  }{2c^{2}_{s}} + \frac{(c \cdot u)^{2}  }{2c^{4}_{s}}\Big]
\end {equation}

for stationery medium, u=0 and f^{eq}_{i} = \rho W_{i}

Modelling external force(Guo et al 2004)

Calculate the updated momentum as

\rho u = \Sigma f_{\alpha}c_{\alpha} + F \frac{\Delta t}{2}

Use the updated value of u for calculating f^{eq} eqn (11)

Calculate the post collision value of PDF f^{*}_{\alpha} as

f^{*}_{i} = f^{*}_{i} + (\Omega _{i}+S_{i})\Delta t


\Omega_{i} = - \frac{(f_{i} - f^{eq}_{i})}{\tau}
S_{i} =  (1- \frac{\Delta t} {2\tau})\omega_{i} \Big (     \frac{c  -u   }{c^{2}_{s}} + \frac{(c \cdot u)c  }{2c^{4}_{s}}   \Big ) \cdot F


  1. Mohamad, A. A. Lattice boltzmann method. Vol. 70. London: Springer, 2011.
  2. Krüger, Timm, et al. “The lattice Boltzmann method.” Springer International Publishing 10.978-3 (2017):4-15.

Leave a Comment

Your email address will not be published.