Table of Contents Author Guidelines Submit a Manuscript
Journal of Applied Mathematics
Volume 2012 (2012), Article ID 135173, 12 pages
Research Article

Flow Simulations Using Two Dimensional Thermal Lattice Boltzmann Method

Department of Mechanical Engineering and Mechanics, Lehigh University, 27 Memorial Drive West, Bethlehem, PA 18015, USA

Received 13 February 2012; Accepted 9 April 2012

Academic Editor: Shuyu Sun

Copyright © 2012 Saeed J. Almalowi and Alparslan Oztekin. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.


Lattice Boltzmann method is implemented to study hydrodynamically and thermally developing steady laminar flows in a channel. Numerical simulation of two-dimensional convective heat transfer problem is conducted using two-dimensional, nine directional D2Q9 thermal lattice Boltzmann arrangements. The velocity and temperature profiles in the developing region predicted by Lattice Boltzmann method are compared against those obtained by ANSYS-FLUENT. Velocity and temperature profiles as well as the skin friction and the Nusselt numbers agree very well with those predicted by the self-similar solutions of fully developed flows. It is clearly shown here that thermal lattice Boltzmann method is an effective computational fluid dynamics (CFD) tool to study nonisothermal flow problems.

1. Introduction

Historically, the Lattice Boltzman method (LBM) evolves from Lattice Gas Cellular Automata. In 1988, LBM is proposed to be used to simulate flows for the first time. The LBM is a branch of statistics of mechanics which is an ideal approach to simulate flows in simple or complex geometries. Recently, LBM has been modified to solve nonlinear partial differential equations to model complex fluid flows. Different approaches of the LBM have been discussed by several investigators [1]. However a successfully LBM simulation rests on the correct implementation of the boundary conditions, where unknown distribution function originated from the operation. As it is stated in several literatures, the implementation of the boundary conditions in LBM is the key to successfully model flow problems. Each type of boundary condition requires different technique and has different degree of accuracy as well.

LBM has been used to model flows in various geometries by several research groups, but only few investigators compare LBM model against other numerical methods for computational fluid dynamics (CFD). Recently Chen et al. and Begum and Basit [2, 3] had proposed models to be applied to simulate complex flows. It has been shown that LBM can easily be implemented to study single- and multiphase flows. The equations governing conservation of mass and momentum are satisfied at each lattice nodes. The LBM approximation is a linear discretized equation which has two terms: streaming and collision terms. The collision term in the LBM has been approximated to be a linear term using models introduced by Bhatngar, Gross, and Krook (BGK) [4]. Recent study by Shi et al. [5] has shown that the LBM is a promising tool to study microscopic flows. Present work is to illustrate that LBM can be an effective CFD tool and in order to demonstrate that 2D developing nonisothermal flows in a channel is studied by implementing LBM method. The results predicted by LBM have been compared against those obtained using ANSYS-FLUENT and those obtained by self-similar solutions in the developed region for validation.

2. Lattice Boltzmann Governing Equation

Ludwig Eduard Boltzmann (1844–1906), the Austrian physicist, had the greatest achievement in the development of statistical mechanics. This approach has been used to predict macroscopic properties of matter such as the viscosity, thermal conductivity, and diffusion coefficient from the microscopic properties of atoms and molecules [68]. The probability of finding particles within certain range of velocities at a certain range of locations replaces tagging each particle as in molecular dynamics simulation. The lattice Boltzmann transportation can be governed by distribution function which represents particles at location at time , and the particle will be displaced by in time with the application force on the liquid molecules [9]. The equation governing the distribution function has two terms, the streaming step and the collision term. Here, and are spatial coordinates, is the time, and is the lattice discrete velocity.

The collision takes place between the molecules; there will be a net difference between the numbers of molecules in the interval . The rate of change of the distribution function is expressed as Here, denotes external forces applied, and is the source or the collision term. With the absence of the external forces, (2.1) becomes Equation (2.2) is known as the lattice Boltzmann governing equation. The right hand side of equation is called a source and is approximated by BGK as Here is the relaxation frequency and the is the relaxation time is the equilibrium value of distribution function and is written as where is the discrete velocities vector, is the bulk fluid velocity, is the lattice sound speed and is the weight factor, one has Equation (2.3) becomes

3. Lattice Boltzmann Arrangements (D2Q9)

Lattice Boltzmann is relatively recent technique that has been shown to be as accurate as traditional CFD methods having ability to integrate arbitrarily complex geometries. LBM can be used for different arrangements such as D1Q2, D2Q4, D2Q9, D3Q15, D3Q19, or D3Q27 [10]. However, in this paper, we only use D2Q9 which implies the two-dimensional and nine velocity components as shown in Figure 1.

Figure 1: D2Q9 LBM arrangement.

Each distribution function has position , velocity , and weight factor .

4. Momentum Lattice Boltzmann Model

The momentum LBM represents the particles velocity [4]. For instance, for D2Q9 lattice arrangements, the particle at the origin is at rest and the remaining particles move in different directions with different speed. Each velocity vector denotes a lattice per unit step. These velocities are exceptionally convenient in that all and components are either 0 or ±1. Mass of particle is taken as unity uniformly throughout the flow domain. The macroscopic fluid density, , is governed by conservation of mass The bulk fluid velocity is the average of microscopic lattice-directional velocity and the directional density and is governed by conservation of momentum Here are and components of the lattice directional velocity. Conservation of mass and momentum is also applied at each boundary, at the inlet and the outlet as shown in Figure 2 for the D2Q9 lattice arrangement for nodes placed on boundaries.

Figure 2: D2Q9 momentum lattice arrangements at boundaries and inside the flow domain.

The uniform inlet velocity is , the length of channel is , and the gap between plates is . is measured in units of and are measured in units of and and , respectively. Scaled inlet velocity becomes , and the flow domain becomes , and .

5. Thermal Lattice Boltzmann Model

There has been rapid progress in developing the construction of stable thermal lattice Boltzmann equation models to study heat transfer problems. McNamara and Zanetti successfully applied multispeed thermal fluid lattice Boltzmann method to solve heat transfer problems [6]. At the outlet, bounce back or extrapolation boundary conditions are considered as the thermal and flow boundary conditions. Bounce back type boundary conditions are proven to provide more accurate numerical approximations [11] and are used by the present work. The temperature at each wall is specified; however, the temperatures which are pointing to the flow domain are unknowns, and they can be evaluated from streaming and collision steps. The thermal lattice arrangement is illustrated in Figure 3.

Figure 3: D2Q9 thermal lattice arrangements at boundaries and inside the flow domain.

The rate of change of the thermal distribution function is written as With normalized temperature , the equilibrium value of the thermal distribution function is given by For simplicity the relaxation frequency of the thermal and momentum distribution function is selected as the same. Hence the kinematic viscosity, , and the thermal diffusivity, , are the same and are expressed by Here, is the wall temperature, and is the temperature distribution at the inlet. The Prandtl number . The temperature of the fluid is governed by the conservation of energy

6. Results and Discussion

The discretized equations (2.6) and (5.1) for momentum and thermal distribution and for nodes in and directions, respectively, are solved by employing Gauss-Siedel iterations. For boundary nodes the detailed discretized equations for and are described in Appendices A and B. The results are presented for steady incompressible two-dimensional laminar flows in an entrance region of a channel. Flow develops hydrodynamically and thermally in a 1?m long and 0.02?m height channel with the aspect ratio AR of 50. At the inlet the flow is uniform (, and the temperature of the fluid satisfies . Boundary conditions imposed on the velocity field at and 1 are no-slip and no-penetration, and the thermal boundary conditions applied on each surface are . The physical properties are considered to be constant and are determined for water at 300?K—( and ). For the example illustrated in this paper, the flow rate considered is 0.4?kg/s and the corresponding Reynolds number .

Spectral convergence is checked for LBM to ensure that the results predicted by the LBM are not dependent on the number of nodes selected for the numerical simulations. Nodes are placed uniformly in the direction of nodes) and nodes). The convergence test is displayed in Figure 4 as the velocity and temperature profile at plotted for various and . It is shown that the mesh provides satisfactory spectral convergence and numerical accuracy and is thereby chosen for the numerical simulation results predicted by LBM.

Figure 4: (a) Velocity and (b) temperature profiles at plotted for various and .

The velocity and temperature profiles are displayed at various cross-sections in the developing region. The profiles that are obtained by LBM are compared against those obtained by ANSYS-FLUENT at the same conditions. The boundary conditions at the inlet and the outlet and on the surface of the plates are selected as the same for both methods. The velocity and temperature field are considered to be converged when the error tolerance is less than 10-3.

The velocity profiles predicted by LBM and FLUENT at various cross-sections are shown in Figure 5. The solid lines denote the prediction obtained by LBM while the symbols denote the predictions obtained by FLUENT. The velocity profiles at all cross-section predicted by LBM agree very well with those predicted by FLUENT, as shown in Figure 5(a). The development length for velocity field at in this flow is expected to be . The thermal field has the same development length as the hydrodynamic field since is selected to be unity. The nearly fully developed velocity profile obtained by both methods at also agrees very well with each other. They also agree well with the analytical solution, , obtained by the self-similar solution for fully developed laminar flow in a channel, as shown in Figure 5(b).

Figure 5: (a) Velocity profiles at various cross-sections predicted by LBM and FLUENT. (b) Velocity profile at predicted by LBM and FLUENT and the self-similar solution for the fully developed laminar channel flow.

The temperature profiles predicted by LBM and FLUENT at various cross-sections are shown in Figure 6. The solid lines denote the prediction obtained by LBM, and the symbols denote the predictions obtained by FLUENT. The temperature profiles predicted by LBM agree very well with those predicted by FLUENT.

Figure 6: Temperature profile at various cross-sections predicted by LBM and FLUENT.

Wall shear stress and heat transfer coefficient are predicted at various cross-sections in the developing region. The local value of skin friction, , and the Nusselt number, , are determined from the numerical solution as where is the bulk temperature of the fluid and is calculated at each cross-section as The skin friction and the Nusselt number predicted by LBM are plotted in Figure 7 as a function . tends to 24 as the fully developed region is approached. Similarly, Nusselt number tends to 7.54 as the thermally fully developed region is approached, as shown in Figure 7(b). These values are in perfect agreement with the fully developed values of and as well documented in the literature.

Figure 7: (a) Skin friction and (b) Nusselt number plotted as a function of .

7. Conclusion

Hydrodynamically and thermally developing laminar steady flow in a channel is considered as an example to illustrate that Lattice Boltzmann method is a promising computational fluid dynamics tool. D2Q9 lattice arrangement is used to predict both velocity and temperature field. Profiles obtained by LBM-D2Q9 and ANSYS-FLUENT agree very well. Away from the inlet as the fully developed region is approached, and the profiles tend to self-similar solutions of developed laminar channel flows. That is confirmed by prediction of the skin friction and Nusselts number. As full-developed region is approached away from the inlet both the skin friction coefficient and the Nusselt number tend to well-documented values in the literature. Extension of LBM method to three-dimensional unsteady complex multiphase flows is natural, and the implementation of LBM to tackle such problems is underway.


A. Discretized LB Equations for the Velocity Field at the Inlet, the Outlet, and at the Surface of Plates

With indices denoting 9 directions of LB D2Q9 arrangements, denoting the nodes placed in the -direction and denoting the nodes placed in the -directions, the distribution functions for velocity, , and temperature field, , are represented by three-dimensional arrays , and for to 9, and .

At the inlet to , the conservation of mass and momentum yield The boundary conditions at the surface of the lower plate imposed on the velocity field give At the top surface , the velocity boundary conditions yield The outlet conservation of mass and momentum is approximated by using second-order extrapolation which yields

B. Discretized LB Equations for the Temperature Field at the Inlet, the Outlet, and at the Surface of Plates

The inlet to conservation of energy gives With the thermal boundary condition at the surface of the lower plate gives With the thermal boundary condition at the surface of the upper plate gives The outlet conservation of energy is approximated by the second-order extrapolation which yields


:Density distribution function
:Local equilibrium density distribution function
:Temperature distribution function
:Local equilibrium temperature distribution function
:Lattice discrete velocity, and are and components
:Bulk velocity of the fluid, and are and components
:Normalized component of the fluid velocity
:Normalized temperature
:Dimensionless relaxation frequency
:Weight factor
:Wall temperature in (°C)
:Lattice sound speed
:Dimensionless collision relaxation time
:Position vector
:Dimensionless and coordinate
:Fluid speed at the inlet
:Temperature profile at the inlet in (°C),
:Density of fluid, kg/m3
:Kinematic viscosity of fluid, m2/sec
:Thermal diffusivity of fluid, m2/sec
:Prandtl number
:Reynolds number
:Bulk temperature in (°C).


One of the authors of this paper is grateful to Saudi Arabia government for its support toward his graduate education.


  1. X. D. Niu, C. Shu, and Y. T. Chew, “A thermal lattice Boltzmann model with diffuse scattering boundary condition for micro thermal flows,” Computers and Fluids, vol. 36, no. 2, pp. 273–281, 2007. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at Scopus
  2. S. Chen, D. Martínez, and R. Mei, “On boundary conditions in lattice Boltzmann methods,” Physics of Fluids, vol. 8, no. 9, pp. 2527–2536, 1996. View at Google Scholar
  3. R. Begum and M. A. Basit, “Lattice Boltzmann method and its applications to fluid flow problems,” European Journal of Scientific Research, vol. 22, no. 2, pp. 216–231, 2008. View at Google Scholar · View at Scopus
  4. Q. Zou and X. He, “On pressure and velocity boundary conditions for the lattice Boltzmann BGK model,” Physics of Fluids, vol. 9, no. 6, pp. 1591–1596, 1997. View at Google Scholar
  5. Y. Shi, T. S. Zhao, and Z. L. Guo, “Finite difference-based lattice Boltzmann simulation of natural convection heat transfer in a horizontal concentric annulus,” Computers and Fluids, vol. 35, no. 1, pp. 1–15, 2006. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at Scopus
  6. G. R. McNamara and G. Zanetti, “Use of the boltzmann equation to simulate lattice-gas automata,” Physical Review Letters, vol. 61, no. 20, pp. 2332–2335, 1988. View at Publisher · View at Google Scholar
  7. M. C. Sukop and D. T. Thorne Jr., Lattice Boltzmann Modeling an Introduction for Geosceintific and Engineers, Springer, Berlin, Germany, 2010.
  8. A. J. C. Ladd and R. Verberg, “Lattice-Boltzmann simulations of particle-fluid suspensions,” Journal of Statistical Physics, vol. 104, no. 5-6, pp. 1191–1251, 2001. View at Publisher · View at Google Scholar · View at Zentralblatt MATH
  9. M. Sbragaglia and S. Succi, “Analytical calculation of slip flow in lattice boltzmann models with kinetic boundary conditions,” Physics of Fluids, vol. 17, no. 9, Article ID 093602, 8 pages, 2005. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at Scopus
  10. S. Succi, The Lattice Boltzmann Equation for Fluid Dynamics and Beyond, Numerical Mathematics and Scientific Computation, The Clarendon Press/Oxford University Press, New York, NY, USA, 2001.
  11. Q. Liao and T.-C. Jen, “Application of Lattice Boltzmann method in fluid flow and heat transfer,” in Computational Fluid Dynamics Technologies and Applications, InTech.