**1 **

國立臺灣大學工學院應用力學研究所 博士論文

### Department of Applied Mechanics College of Engineering National Taiwan University

### Doctoral Dissertation

半古典晶格波滋曼方法

Semiclassical Lattice Boltzmann Method

洪立昕 Li-Hsin Hung

指導教授：楊照彥 博士

Advisor: Jaw-Yen Yang, Ph.D.

中華民國 99 年 10 月 October, 2010

### 終於走到這一步了，四年多前的熱情慢慢累積成一行行的程式、

### 一篇篇的文章，也差不多消耗殆盡了。很高興也很欣慰能在不算太久 的時間內走出台大大門。最感謝的是父母的支持，無論是精神上還是 經濟上，都是給我最多幫助的。妹妹來自美國的鼓勵也幫助我不少，

### 尤其妳比我更累的研究生活，更是鞭策我繼續向上的動力。謝謝珮君 和冠孜，妳們的陪伴帶給我豐富而精采的精神生活，方能戰勝這漫長 而艱辛的一役。

### 感謝指導教授 楊照彥博士四年多來的提攜，教授不論是嚴謹的 研究態度還是充沛的精力跑遍各國都另人景仰和欽羨。

### 論文定稿之際，承蒙本校 許文翰教授、黃美嬌教授，清華大學 林 昭安教授，中正大學 何正榮教授，成功大學 陳朝光教授、黃吉川教 授以及楊玉姿教授的指導，讓我看到許多不同思考層面，更進一步精 進論文內容，令我受益匪淺。

### 實驗室待了許久，相處過的學長學弟不可勝數，每個人都是我永

### 遠的好友。特別感謝育炘帶我進入這全新的研究領域、澤揚在博後階

### 段和我的討論，還有聖鑫及耀天在部分研究內容上的幫忙，讓我能夠

### 發更多文章，順利畢業。最後，感謝許多我沒提及但仍愛我的人，少

### 了你(妳)們，這一切也難以成真，謝謝。

## 摘要

### 本論文提出了半古典晶格波滋曼法，此方法是在波滋曼法的 基礎下，使用量子統計(Bose-Einstein Statistics 和 Fermi-Dirac Statistics)取代古典統計

### (Maxwell-Boltzmann Statistics)後展開所得的方法。此方 法可驗證在古典極限(Classical Limit)下可回復原本古典 晶格波滋曼法，亦可由數值驗證得知古典極限下可得傳統晶 格波滋曼法之結果。本文以此一新提出之模擬方法研究量子 氣體運動問題，模擬之問題忽略粒子間交互作用，但因為使 用量子統計之故，粒子在量子統計下的特性如庖立不相容原 理(Pauli Exclusion Principle)、巨觀上傳輸係數之修正 等等均有考慮在內。本文探討一維量子氣體之震波管問題、

### 二維圓柱流問題、二維微管問題以及三維頂蓋流問題作為驗 證半古典晶格波滋曼法之方式，模擬結果指出了量子統計和 古典統計結果的主要差異。此外，本文亦提出了另一種雙分 布函數熱晶格波滋曼法，可得出一般雙分布函數熱晶格波滋 曼法所得結果，也可作為未來更進一步推廣半古典晶格波滋 曼法之方向。

### 關鍵詞:晶格波滋曼法、量子統計、波滋曼方程式

### Unlike describing the physical phenomenon in coordinate or momentum spaces in quantum mechanics, semiclassical Boltzmann equation treats the system in phase space, and it is much easier to describe the dynamics of quantum gases. In this thesis, a class of semiclassical lattice Boltzmann methods is developed for solving quantum hydrodynamics and beyond.

### The present method is directly derived by projecting the Uehling-Uhlenbeck Boltzmann-BGK equations onto the tensor Hermite polynomials following Grad's moment expansion method. The intrinsic discrete nodes of the Gauss-Hermite quadrature provide the natural lattice velocities for the semiclassical lattice Boltzmann method. Formulations for the second-order and third order expansion of the semiclassical equilibrium distribution functions are derived and their corresponding hydrodynamics are studied. Gases of particles of arbitrary statistics can be considered. Simulations of one-dimensional compressible gas flow by using D1Q5 lattices, two dimensional microchannel flow, two dimensional flow over cylinder by using D2Q9 lattices and three dimensional lid driven cavity flow by using D3Q19 lattices are provided for validating this method. It is shown that the classical flow patterns such like vortex and vortices shedding in flow over cylinder simulations, temperature and pressure contours together with streamline patterns could be produced from the present method in classical limit. The results also indicate the distinct characteristics of the effects of quantum statistics when they are compared with fluid phenomena in classical statistics.

### Keywords: Lattice Boltzmann Method, Semiclassical, Quantum.

## Contents

LIST OF TABLES . . . v

LIST OF FIGURES . . . 1

1 Introduction 1 1.1 Overview. . . 1

1.2 Conventional LB method . . . 3

1.3 Basics of Semiclassical LB method. . . 5

1.4 Contents of the Dissertation . . . 8

2 Semiclassical Kinetic Theory 9 2.1 Overview. . . 9

2.2 Introduction to Quantum Gases . . . 10

2.3 Semiclassical Boltzmann Equation(UUB Equation) . . . 12

2.4 Semiclassical Hydrodynamic Equations . . . 14

3 Semiclassical Lattice Boltzmann Method 21 3.1 Overview. . . 21

3.2 The Derivations of Conventional LB method . . . 22

3.2.1 Time Discretization . . . 22

3.2.2 Space and Velocity Discretization . . . 23

3.3 The Derivations of Single Relaxation Time Semiclassical LB method . . . 28

3.3.1 Summary of Single Relaxation Time Semiclassical LB Method . . . 32

3.4 Generalized Semiclassical LB method . . . 35

4 Initial Conditions and Boundary Conditions for Semiclassical

LB method 43

4.1 Overview. . . 43

4.2 Initial Conditions . . . 45

4.3 Boundary Conditions . . . 45

4.3.1 Periodic Boundary Condition . . . 45

4.3.2 Bounce Back Boundary Condition . . . 46

4.3.3 Non-Equilibrium Extrapolation Boundary Condition . 49 4.3.4 Immersed Boundary Condition . . . 52

4.3.5 Issues on Microchannel Boundaries . . . 53

5 Numerical Results 63 5.1 Overview. . . 63

5.2 One Dimensional Shock Tube . . . 64

5.3 Two Dimensional Microchannel Flow . . . 65

5.4 Two Dimensional Flow over Cylinder . . . 69

5.5 Two Dimensional Natural Convection and Rayleigh-Benard Con- vection Flow . . . 71

5.6 Three Dimensional Lid Driven Cavity Flow. . . 73

6 Conclusions and Future Work 101 6.1 Accomplishments . . . 101

6.2 Future Work. . . 103

A Nomenclature 107 B Chapman-Enskog Analysis of semiclassical LB method 111 B.1 Chapman-Enskog Analysis of Single Relaxation Time semiclas- sical LB method . . . 111

Contents iii

C Derivations of Double Distribution Function LB method 117

C.1 Two Relaxation Times Kinetic Model . . . 118

C.2 Expansion of the Equilibrium Distribution Functions . . . 121

C.3 Discretization of Velocity Space . . . 122

C.4 Time Discretization . . . 125

C.5 Chapman-Enskog Analysis of Double Distribution Functions LB method . . . 127

Bibliography 133

## List of Tables

3.1 D1Q5 quadratures and weights . . . 41 3.2 D2Q9 quadratures and weights . . . 41 3.3 D3Q19 quadratures and weights . . . 41 5.1 Comparisons of drag coecient of dierent boundary condi-

tions in LB method . . . 76 5.2 Comparisons of the average Nusselt number and the maximum

velocity components across the cavity center. The data in parentheses are the locations of the maxima. . . 76

## List of Figures

*2.1 Bose and Fermi functions in dierent fugacity z. (a) Bose*
*Function g**1.5*^{Bose}*(z), g**2.5*^{Bose}*(z)*, and ^{g}_{g}^{Bose}^{2.5}^{Bose}^{(z)}

*1.5* *(z)*, (b) Fermi Function
*g*_{1.5}^{F ermi}*(z),g*^{F ermi}*2.5* *(z)*, and ^{g}_{g}^{F ermi}^{2.5}^{F ermi}^{(z)}

*1.5* *(z)*, (c) Comparing Bose ^{g}_{g}^{Bose}^{2.5}^{Bose}^{(z)}

*1.5* *(z)*,
and Fermi Function ^{g}_{g}^{F ermi}^{2.5}*F ermi*^{(z)}

*1.5* *(z)* . . . 19

3.1 Comparisons of streamlines of uniform ow over a circular cylin-
der in a BE quantum gas in dierent expansion order with
*z = 0.2* *and Re**∞* = 40*. (a) N = 3, (b) N = 2.* . . . 40

3.2 D1Q5 lattice structure . . . 41

3.3 D2Q9 lattice structure . . . 41

3.4 D3Q19 lattice structure . . . 41

4.1 Layout of the curved wall boundary in the regularly spaced lattices. . . 57

4.2 Periodic boundary conditions. . . 57

4.3 Layout of the curved wall boundary in the regularly spaced lattices. . . 58

*4.4 Bounce back boundary conditions (a) q < 1/2, (b) q ≥ 1/2.* . . 59

4.5 The inuences of dierent accommodation coecients on ve- locity prole based on two dierent boundary conditions. (a) Bounce-back specular reection boundary condition, (b) Dis- crete maxwellian boundary condition. . . 60

4.6 The inuences of dierent boundary conditions on velocity pro-
*le based on the same accommodation r = 0.8. (a) Kn=0.5,*
(b) Kn=2. . . 61

4.7 The velocity prole of dierent Knudsen number based on the
*same accommodation coecient r = 0.8 for two dierent bound-*
ary conditions. (a) Bounce-back specular reection boundary
condition, (b) Discrete maxwellian boundary condition. . . 62

5.1 Convergence of solution with rened grids, the cells number is 100, 200, and 400. (a) Velocity, (b) Density, (a) Pressure, (d) Temperature. . . 78 5.2 Shock tube ow in a quantum gas. The eect due to dier-

ent particle statistics; Bose-Einstein: BE, Fermi-Dirac: FD, Maxwell-Boltzmann: MB. (a) Density, (b) Velocity, (c) Pres- sure, (d) Temperature. . . 80 5.3 Shock tube ow in a quantum gas. The eect in high tem-

perature limit; Bose-Einstein: BE, Fermi-Dirac: FD, Maxwell- Boltzmann: MB. (a) Density, (b) Velocity, (c) Pressure, (d) Temperature. . . 82 5.4 Velocity proles in a channel ow of gases of arbitrary statistics

*(gas with z = 0.2). (a) Kn=0.05, (b) Kn=0.5, (c) Kn=1. (d)*
*Normalized mass ux in a channel ow as a function of Kn*
number. . . 83

*5.5 Drag and lift coecients in classical limit. (a) Re = 100, (b)*
*Re = 500*. . . 84
5.6 Drag and lift coecients in dierent boundary conditions (im-

mersed boundary condition and bounce back boundary condi-
*tion). (a) Re = 100, (b) Re = 200.* . . . 85

List of Figures ix

5.7 Streamlines of uniform ow over a circular cylinder in a quan-
*tum gas with z = 0.2 and Re** _{∞}* = 20. (a) BE gas, (b) MB gas,
(c) FD gas. . . 86
5.8 Streamlines of uniform ow over a circular cylinder in a quan-

*tum gas with z = 0.2 and Re** _{∞}* = 40. (a) BE gas, (b) MB gas,
(c) FD gas. . . 87
5.9 Streamlines of uniform ow over a circular cylinder in a quan-

*tum gas with z = 0.2 and Re**∞*= 200. (a) BE gas, (b) MB gas,
(c) FD gas. . . 88
5.10 Dierent contours and streamlines of uniform ow over a circu-

*lar cylinder in a quantum gas in a classical limit with Re**∞* =
200*, steps=94200. (a) X-direction velocity(u**x*), (b) y-direction
*velocity(u**y**), (c) fugacity(z), (d) density(ρ), (e) temperature(θ),*
(f) streamlines. . . 89
5.11 Streamlines for dierent Rayleigh numbers of natural convec-

*tion. (a) Ra = 10*^{3}*, (b) Ra = 10*^{4}*, (c) Ra = 10*^{5}. . . 90
5.12 Isotherm lines for dierent Rayleigh numbers of natural con-

*vection. (a) Ra = 10*^{3}*, (b) Ra = 10*^{4}*, (c) Ra = 10*^{5}. . . 90
5.13 Streamlines of Rayleigh-Benard convection for dierent Rayleigh

*numbers. (a) Ra = 4000, (b) Ra = 10000, (c) Ra = 50000.* . . 91
5.14 Isotherms of Rayleigh-Benard convection for dierent Rayleigh

*numbers. (a) Ra = 4000, (b) Ra = 10000, (c) Ra = 50000.* . . 92
*5.15 Grid renement of three dimensional lid driven simulation, Re =*

400 for two dierent centerlines. (a) (x,0,0), (b) (0,0,z). . . 93
*5.16 Grid renement of three dimensional lid driven simulation, Re =*

1000 for two dierent centerlines. (a) (x,0,0), (b) (0,0,z). . . . 93

*5.17 Streamlines of MB statistics in three mid-planes, Re = 100. (a)*
xy planes, (b) zy planes, (c) xz planes. . . 94
*5.18 Streamlines of MB statistics in three mid-planes, Re = 400. (a)*

xy planes, (b) yz planes, (c) xz planes. . . 94
*5.19 Streamlines of MB statistics in three mid-planes, Re = 1000.*

(a) xy planes, (b) yz planes, (c) xz planes. . . 94 5.20 Centerlines velocity of classical and semiclassical model. (a)

Classical, (b) semiclassical . . . 95
*5.21 Velocity vectors of BE statistics in three mid-planes, Re = 100,*

*z = 0.2*. (a) xy planes, (b) yz planes, (c) xz planes. . . 96
*5.22 Pressure contours of BE statistics in three mid-planes, Re =*

100*, z = 0.2. (a) xy planes, (b) yz planes, (c) xz planes.* . . . . 96
*5.23 Streamlines patters of BE statistics in three mid-planes, Re =*

100*, z = 0.2. (a) xy planes, (b) yz planes, (c) xz planes.* . . . . 96
*5.24 Velocity vectors of FD statistics in three mid-planes, Re = 100,*

*z = 0.2*. (a) xy planes, (b) yz planes, (c) xz planes. . . 97
*5.25 Pressure contours of FD statistics in three mid-planes, Re =*

100*, z = 0.2. (a) xy planes, (b) yz planes, (c) xz planes.* . . . . 97
*5.26 Streamlines patters of FD statistics in three mid-planes, Re =*

100*, z = 0.2. (a) xy planes, (b) yz planes, (c) xz planes.* . . . . 97
*5.27 Streamlines patters of FD statistics in three mid-planes, Re =*

400*, z = 0.2. (a) xy planes, (b) yz planes, (c) xz planes.* . . . . 98
*5.28 Streamlines patters of BE statistics in three mid-planes, Re =*

400*, z = 0.2. (a) xy planes, (b) yz planes, (c) xz planes.* . . . . 98
*5.29 Streamlines patters of FD statistics in three mid-planes, Re =*

1000*, z = 0.2. (a) xy planes, (b) yz planes, (c) xz planes.* . . . 98

List of Figures 1

*5.30 Streamlines patters of BE statistics in three mid-planes, Re =*
1000*, z = 0.2. (a) xy planes, (b) yz planes, (c) xz planes.* . . . 99
*5.31 Centerlines velocity of dierent quantum statistics, z = 0.1 and*

*z = 0.2* . . . 100

### Chapter 1

## Introduction

### Contents

1.1 Overview . . . 1

1.2 Conventional LB method . . . 3

1.3 Basics of Semiclassical LB method . . . 5

1.4 Contents of the Dissertation . . . 8

### 1.1 Overview

After more than twenty years developments on lattice Boltzmann (LB) method since its introduction [1][2], the LB method is not just the extension of its ancestor lattice-gas automaton [3], and even not just a hydrodynamical equa- tions solver. It has been another approach analyzing Boltzmann equation in compared with Chapman-Enskog expansion [4] and Grad moment method [5].

LB method is based on the kinetic equations for simulating uid ow, see [6][7][8]. Over the past two decades, signicant advances in the development of the LB method [2][9][10][11] based on classical Boltzmann equations with the relaxation time approximation of Bhatnagar, Gross and Krook (BGK) [12] have been achieved. The LB method has demonstrated its ability to simulate hydrodynamic systems, magnetohydrodynamic systems, multi-phase

and multi-component uids, multi-component ow through porous media, and complex uid systems, see [13]. The LB equations can also be directly derived in a priori manner from the continuous Boltzmann equations [14][15][16]. Most of the classical LB methods are accurate up to the second order, i.e., Navier- Stokes hydrodynamics and have not been extended beyond the level of the Navier-Stokes hydrodynamics. A systematical method [17][16] was proposed for kinetic theory representation of hydrodynamics beyond the Navier-Stokes equations using Grad's moment expansion method [5][18]. The use of Grad's moment expansion method in other kinetic equations such as quantum kinetic equations and Enskog equations can be found in [19][20].

Despite their great success, however, most of the existing LB methods are limited to hydrodynamics of classical particles. Modern development in nanoscale transport requires carriers of particles of arbitrary statistics, e.g., phonon Boltzmann transport in nanocomposite and carrier transport in semi- conductors. The extension and generalization of the successful classical LB method to quantum LB method for quantum particles is desirable. Analogous to the classical Boltzmann equations, a semiclassical Boltzmann equations for transport phenomenon in quantum gases has been developed by Uehling and Uhlenbeck [21][22]. Following the work of Uehling and Uhlenbeck based on the Chapman-Enskog procedure [4], the hydrodynamic equations of a trapped dilute Bose gas with damping have been derived [23]. In [19], the quantum Grad expansion using tensor Hermite polynomials has been applied to obtain the non-equilibrium density matrix which reduces to the classical Grad mo- ment expansion if the gas obeys the Boltzmann statistics. The full Boltzmann equations is mathematically dicult to handle due to the collision integral in dierent types of collisions. To avoid the complexity of the the collision term, the relaxation time model originally proposed by BGK model for the clas-

1.2. Conventional LB method 3

sical non-relativistic neutral and charged gases has been widely used. Also, BGK-type relaxation time models to capture the essential properties of carrier scattering mechanisms can be similarly devised for the UUB equations for var- ious carriers and have been widely used in carrier transports [24]. Recently, kinetic numerical methods for ideal quantum gas dynamics based on Bose- Einstein(BE) and Fermi-Dirac(FD) statistics have been presented [25][26]. A gas-kinetic method for the semiclassical Boltzmann-BGK equations for non- equilibrium transport has been devised [27].

In this thesis, a new semiclassical LB method for the UUB-BGK equations based on Grad's moment expansion method by projecting the UUB-BGK equations onto Hermite polynomial basis has been derived. The relations between the relaxation time, viscosity and thermal conductivity are obtained by applying the Chapman-Enskog method [4] to the UUB-BGK equations for providing the basis for determining relaxation time used in the present semiclassical LB method. Hydrodynamics based on moments up to second and third order expansions are presented. Computational examples to illustrate the methods are given and the eects due to quantum statistics are delineated.

In the rest of this dissertation, the conventional LB method and quantum LB method will be introduced, then comes the basics of semiclassical LB method which was rst proposed in [28].

### 1.2 Conventional LB method

This section will describe the basic algorithm of the LB method. First, an overview of the historical development of LB method including recently devel- opments and reviews will be given. Next, the method itself will be described.

Finally, a derivation of the necessary equations and the lattice structures will

be discussed. The conventional LB method originated from its predecessor, the lattice gas cellular automata models [3] has become a competitive nu- merical tool for simulating uid ows over a wide range of complex physical problems [2] [29] [30] [31] [32] [33]. The theoretical background of LB method is the kinetic theory and Boltzmann equation, which are connected with the macroscopic Navier-Stokes equation by the Chapman-Enskog expansion. The LB method can be regarded as a simplied kinetic scheme by using a nite set of discrete velocities and a simplied collision integral. Both the algorithm and the boundary conditions are easy to implement in LB method. As the LB method computes macroscopic behavior, such as the motion of a uid, with equations describing microscopic scales, it operates on a mesoscopic level in between those two extremes. When compared to the traditional computa- tional uid dynamics techniques, the advantages of LB method are mostly on its clear physics and simple algorithm. While conventional solvers directly dis- cretize the Navier-Stokes equations, the LB method is essentially a rst order explicit discretization of the Boltzmann equation in a discrete phase-space. It can also be shown, that the LB method approximates the Navier-Stokes equa- tions with good accuracy. A good review of LB method can be found, e.g., in [6][7][8]. In LB method, the simulation region is divided into a cartesian grid of cells, each of which only interacts with cells in its direct neighbor- hood. The LB method consists of two steps, the stream step, and the collide step. These are usually combined with no-slip boundary conditions for the domain boundaries or obstacles. The simplicity of the algorithm is especially evident when implementing it. Using a LB method, the particle movement is restricted to a limited number of directions. As shown in Fig. 3.3 and Fig. 3.4, a two-dimensional model with 9 velocities (commonly denoted as D2Q9), and three-dimensional model with 19 velocities (commonly denoted

1.3. Basics of Semiclassical LB method 5

as D3Q19) are provided.

### 1.3 Basics of Semiclassical LB method

Although Boltzmann Equation has been successfully applied to many eld like dilute gas dynamics, multi-scale simulations, however, in recent years, micro and nano technology has been emerged quickly, and the transport phenomena in semiconductors at low temperatures is very important. There have been a successful theory in statistical mechanics which can predict the transport coecients like shear viscosity and thermal conductivity of ideal quantum

uid such like electrons in the metal. The questions arise of whether quantum systems like that can be described similar to the one developed for the clas- sical counterpart. When solving these kinds of problems, classical Boltzmann Equation is not enough and it require quantum mechanical treatments. In quantum mechanics identical particles are absolutely indistinguishable from one another and N-particle system can be described by a wave function with permutation symmetry. In nature, it is found that particles with antisymmet- ric wave functions are called fermions which obey FD statistics and particles with symmetric wave functions are called bosons which obey BE statistics.

The statistical properties of fermion and boson systems are profoundly dier- ent at low temperature. However, in the classical limit, both quantum dis- tributions reduce to the Maxwell-Boltzmann(MB) distribution. Boltzmann equation describes the dynamic behavior of ideal gas by a single-particle dis- tribution function.

In general, there are three strategies to take for statistically treating a quantum system [34]. One is to use a kinetic equation governing the density matrix, another one is to use a kinetic equation with the Wigner distribution

function, the last one, which is also the method used in this thesis, is to assume a semiclassical kinetic equation such as the UUB [21][22] kinetic equation as a generalization of the classical one. In UUB equation, the collision term of the Boltzmann equation is rewritten with the quantum particle scattering form.

The equilibrium distribution to semiclassical Boltzmann equation is the BE or FD distributions. The particles obey BE distribution are called bosons and FD distribution for fermions, no third category has yet been found.

Considering of the recently successful developing on LB method and well deriving semiclassical Boltzmann equation, it is natural to extend the conven- tional LB method to the semiclassical LB method for dealing with dierent quantum statistics. Although UUB equation introduced above has considered the quantum statistics and would reveal quantum eects in specic problems.

However, dealing with UUB equation is still a challenge. In this dissertation, we proposed a new semiclassical LB method which extends well known LB method to solve semiclassical Boltzmann equation with BGK approximation.

In [28], a new semiclassical LB-BGK method had been developed, and this
method would describe quantum systems in dierent approaches. The follow-
ing works about rareed channel ow [35] and axisymmetric ow of quantum
statistics [36] present dierent applications of this new method. The idea
of extending the conventional LB method to semiclassical LB method is to
adopt the Grad's moment method to nd solutions to semiclassical Boltz-
**mann equation by expanding f(x, ζ, t) in terms of Hermite polynomials. And***in simulations, the N-th nite order truncated distribution function f** ^{N}* was
considered. the details will be shown in chapter 3.

It is worth mentioning here the dierences between the present method and quantum LB method which is proposed by Succi and Benzi [37][38]. In the quantum LB method, Dirac equation which is the most general equation de-

1.3. Basics of Semiclassical LB method 7

scribing single particle motion in compliance with quantum theory and special
relativity is solved by LB method. In [37], the procedure builds on a formal
analogy between the Dirac equation and a special discrete kinetic equation
*known as LB method, it was then shown that the non-relativistic Schr¨odinger*
equation ensues from the Dirac equation under an adiabatic assumption that
is formally similar to the one which takes the Boltzmann equation to the
Navier-Stokes equations in kinetic theory. In [38], it was further shown that
by a proper resort to operator splitting methods, the Dirac equation can be in-
tegrated as a sequence of three one dimensional LB equation evolving complex
valued distribution function. In these works, LB method with complex distri-
bution function is treated as a numerical tool for solving a complex equation.

By using multiscale technique and the Chapman-Enskog expansion on com- plex variables, the complex partial dierential equations could be recovered.

Recent years, this procedure is applied for solving many dierent equations, for examples, the one-dimensional nonlinear Dirac equation[39], the nonlin- ear convection-diusion equations in [40] and the complex Ginzburg-Landau equation in [41]. It should be also noticed that several quantum lattice gas cellular automata methods [37][42][43][44] have been recently presented which applying and extending the concept of classical lattice gas cellular automata models to treat the time evolution of wave functions for spinning particles and the Schrödiger equation or the Dirac equation directly. For a more detailed review, see [45]. However, present semiclassical LB method associated with the works presented in [25][26][27] are based on the semiclassical kinetic de- scription. i.e., the particle motion (velocity or momentum) and position are treated in classical mechanics manner while the particles can be of quantum statistics. The procedure and physical meanings are much dierent.

### 1.4 Contents of the Dissertation

This Dissertation is organized as following.

Chapter 1 gives an overview of this dissertation, briey describes the basic ideas of semiclassical LB method which combines the elements of semiclassi- cal kinetic theory and LB method. History and recent developments of LB method are also given.

Chapter 2 describes the semiclassical kinetic theory. First, some important concepts of quantum kinetic theory are introduced, then comes the semiclassi- cal Boltzmann equation which has been developed for several years since UUB equation [21][22]. Following the previous works, the corresponding semiclas- sical hydrodynamic equations are developed and compared with the classical hydrodynamic equations.

In chapter 3, the semiclassical LB method is derived based on expanding the BE and FD equilibrium distribution function onto the Hermite polyno- mials. The analysis of semiclassical LB method is also given. In chapter 4, the initial condition and the boundary conditions of the semiclassical LB method based on conventional LB method are introduced. Since LB method is a mature uid simulating method, extending the boundary conditions from conventional LB method to semiclassical LB method is straightforward.

The numerical validation and examples are given in chapter 5, including one dimensional shock tube simulation, two dimensional ow over cylinder, natural convection ow, and three dimensional lid-driven cavity simulation.

In chapter 6 the conclusions and future works are given. In Appendix 1, the Chapman-Enskog analysis of semiclassical LB method is given. Finally, a new thermal LB model is proposed for the future development of semiclassical LB method.

### Chapter 2

## Semiclassical Kinetic Theory

### Contents

2.1 Overview . . . 9 2.2 Introduction to Quantum Gases . . . 10 2.3 Semiclassical Boltzmann Equation(UUB Equation) . 12 2.4 Semiclassical Hydrodynamic Equations . . . 14

### 2.1 Overview

In quantum mechanics, the term semiclassical has dierent meanings and all refer to some approximations or situations that combine quantum and classical properties, for example, Wentzel-Kramers-Brillouin(WKB) approxi- mation or UUB equation. In this chapter, from the basic concepts to quantum gases, the collision term of semiclassical Boltzmann equation is derived. Then, the approximation method, BGK method is introduced, and the equilibrium distribution function is proposed based on Boltzmann H theory. Finally, semi- classical hydrodynamic equations are derived from moment integration of the semiclassical Boltzmann equation with BGK approximation and the results are compared with their classical counterparts.

### 2.2 Introduction to Quantum Gases

Comparing with classical mechanics statistics and quantum mechanics statis- tics, distinguishability is a very important concept. In classical mechanics statistics, particles in a system are distinguishable even if they are identical. It is possible to label particles and track their phase space trajectories with cer- tainty. Conversely, in quantum mechanics statistics the precise descriptions of particles motions are limited by Heisenberg uncertainty relation where the si- multaneous measurements on momentum and position are not allowed. Under these fuzzy descriptions, the particles become untrackable and the property of indistinguishability is inherent. The distinguishable property makes the conguration dierent and the corresponding microstate distinct as particle in- terchanges. The number of microstates are usually very large and microstates will continually change as time goes on. The time-average behavior of the system will be identical to the average behavior of all sort microstates cor- responding to the macrostates at time t. The average behaviors on all the collection of microstates in any systems are called ensembles of the systems.

In statistical theory, there are three dierent statistics which are microcanoni-
cal ensemble, canonical ensemble and grand canonical ensemble. In canonical
*ensemble the macrostates are characterized by (N, V, T ), and the partition*
*function is written as Q**N**(V, T ) =*∑

*E**e*^{−βE}*, wherein β = 1/k**b**T*. In the ideal
*quantum gas system, the energy eigenvalues E are composed of single particle*
*states in occupation set {n**ε**}, ε is the single-particle energy, and E can be*
*expressed in terms of ε as E =* ∑

*ε**n*_{ε}*ε. The values of the numbers n**ε* must
satisfy the condition

∑

*ε*

*n*_{ε}*= N* (2.1)

According to the discussions of distinguishable property given above, we

2.2. Introduction to Quantum Gases 11

can write down the statistical weight factor for FD, BE and MB statistics as:

*W**F D**{n**ε**} = {* 1 *if all ε = 0 or 1*

0 otherwise *},* (2.2)

*W*_{BE}*{n**ε**} = 1,* (2.3)

*W** _{M B}* =∏

*ε*

1

*n** _{ε}*!

*.*(2.4)

We can further write the partition function of a N particle system
*Q*_{N}*(V, T ) =*

∑*′*
*{n**ε**}*

*W{n**ε**}e*^{−β}^{∑}^{ε}^{n}^{ε}* ^{ε}* (2.5)

∑_{′}

means the summation over all distribution sets that satises (2.1). To
work out the partition function within ∑_{′}

is relatively complicated since the summation is under some constraints. Now substituting (2.3) and (2.2) into (2.5), one has

*Q*_{N}*(V, T ) =*

∑*′*
*{n**ε**}*

*e*^{−β}^{∑}^{ε}^{n}^{ε}^{ε}*.*

The dierences between BE and FD statistics arises from the values that
*the number n**ε* can take. However, the restriction (2.1) makes the canonical
ensemble cumbersome here. The grand partition function seems to provide
an ecient way

*ϑ(z, V, T )≡*

∑*∞*
*N =0*

*z*^{N}*Q*_{N}*(V, T ) =*

∑*∞*
*N =0*

[

∑*′*
*{n*^{ε}*}*

∏

*ε*

*(ze** ^{−βε}*)

^{n}*]*

^{ε}= ∑

*{n*0*,n*1*,··· }*

*(ze*^{−βε}^{0})^{n}^{0}*(ze*^{−βε}^{1})^{n}^{1}(*· · · )*

=

∏

*ε*

(1*− ze** ^{−βε}*)

^{−1}*BE case, ze*

^{−βε}*< 1*

∏

*ε*

*(1 + ze** ^{−βε}*) FD case

*The average occupation number in state ε*

*< n**ε**> =* 1
*ϑ*[*−*1

*β*(*∂ϑ*

*∂ε*)*z,T ,all* *other* *ε*]

= 1

*z*^{−1}*e*^{βε}*+ η*

*where η = −1 for bosons, η = −1 for fermions, η = 0 for MB gas. It*
is convenient to dene the thermodynamic potential q to nd the average
number and energy

*q(z, V, T )≡ lnϑ =* *P V*

*kT* *= η*∑

*ε*

*ln(1 + ηze** ^{βε}*) (2.6)

*N = z(*¯ *∂q*

*∂z*)* _{V,T}* =∑

*ε*

1

*z*^{−1}*e*^{βε}*+ η* (2.7)

*E =*¯ *−(∂q*

*∂β*)* _{z,V}* =∑

*ε*

*ε*

*z*^{−1}*e*^{βε}*+ η* (2.8)

### 2.3 Semiclassical Boltzmann Equation(UUB Equa- tion)

Boltzmann equation is successful in describing dilute classical monatomic gases. People showed the derivation of the Boltzmann equation from the Newtonian of motion for the many particle system under suitable assump- tions. Hereinafter, the semiclassical Botlzmann equation (also called UUB equation) is used for describing the ow of fermions and bosons. The gener- alized Boltzmann equation is

*∂*_{t}*f +* **p**

*m* *· ∇***x****f + F***· ∇***ξ*** f = R(f ),* (2.9)

*to obey the streaming mythology like the left hand side of the Boltzmann*

**where F = ma is the external force. The quantum particles are assumed**2.3. Semiclassical Boltzmann Equation(UUB Equation) 13

equation. It should be emphasized here that the major and only dierence
between semiclassical Boltzmann equation and classical Boltzmann equation
* is the collision term R(f). In the classical one, collision term is derived based*
on classical mechanics [46] while it is derived by quantum scattering theory
in the semiclassical one [23]. The collision term for semiclassical Boltzmann
equation is described below [21].

**R(f ) =**

∫
**dp**_{1}

∫

**dΩK(p, p**_{1}*, Ω) {[1 + ηf(p, t)][1 + ηf(p*

**1**

**)]f (p**

^{∗}

**, t)f (p**

^{∗}

_{1}*, t)*(2.10)

**−[1 + ηf(p**^{∗}**, t)][1 + ηf (p**^{∗}_{1}**)]f (p, t)f (p****1***, t)},*

*where the function K is the collision kernel, Ω is the solid angle. η = +1 de-*
*notes the case of BE statistics and η = −1 denotes the case of FD statistics.*

We observe that the Pauli exclusion principle is included in the case of fermion
and a population factor for boson is included. It is noted that the Boltzmann
equation of the classical statistics is included in (2.10) as a special case when
*η = 0*. The population in nal state is neglected. This shows the classical col-
lision term can be recovered when considering dilute, non-degenerate classical
gases.

Although semiclassical Boltzmann equation has been derived from consid-
ering quantum particles collisions, the complexity of the collision term is still
a obstacle for computation. One straightforward way to solve this problem
**is to approximate the complex collision term R(f) with BGK approximation**

*f**−f*^{eq,Q}

*τ* which was rst proposed in [12]. BGK approximation makes the di-
rect solution of the BGK-Boltzmann equation tractable. After introducing
the BGK approximation, (2.9) could be rewritten as

*∂*

*∂t f + ξ· ∇*

**x**

**f + a**· ∇

**ξ***f =−f*

*− f*

^{eq}*τ* *.* (2.11)

*In this BGK-Boltzmann equation, f** ^{eq}* could be quantum equilibrium distri-

*bution function f*

^{eq,Q}*or classical equilibrium distribution function f*

^{eq,C}*. f*

*means the probability of particles (electrons, phonons, photons, etc...) occu- pying a specic quantum state in equilibrium, such that, the physical meaning*

^{eq}*of f*

^{eq}*is the same as < n > which has been derived in last section. The gen-*eralized equilibrium distribution function is

*f*^{eq,Q}*(ε) =< n*_{ε}*>=* 1

*z*^{−1}*e*^{βε}*+ η* (2.12)

*where z ≡ e*^{−βµ}**is the fugacity and ε = m(ξ − u)**^{2}*/2*is the energy of particle.

*The parameter η determines which statistics f** ^{eq,Q}* belong to.

### 2.4 Semiclassical Hydrodynamic Equations

The semiclassical hydrodynamic equations are obtained by taking moments
**ψ = [1, mξ, mξ**^{2}*/2]*on the semiclassical Boltzmann equation of (2.9) with the
collision term (2.10), then integrating the resulting equations over all ξ.

*∂*_{t}

∫

**ψf dΞ +**∇**x**

∫

**ψξf dΞ****− a**

∫

*ψ∇***ξ****f dΞ =**

∫

* R(f )ψdΞ* (2.13)
The integrals of the collision terms in all three cases should preserve the
conservation property. That means the conservation of mass, momentum
and energy need to be satised all the time, which is called the compatibility
condition,

∫

**R(f )ψdΞ = 0,***ψ =*

1
**mξ**

1

2*m(ξ*_{x}^{2}*+ ξ*_{y}^{2}*+ ξ*_{z}^{2})

*,* (2.14)

**where dΞ = m**^{3}**dξ/h**^{3} is the innitesimal volume in momentum space. We
note that the collision in BGK model should also preserve the compatibility

2.4. Semiclassical Hydrodynamic Equations 15

condition,

1
*τ*

∫

*(f* *− f*^{eq}**)ψdΞdt = 0,***ψ =*

1
*ξ**i*
1

2*m(ξ*_{i}^{2}*+ ξ*_{j}^{2}*+ ξ*_{k}^{2})

*.* (2.15)

The denitions of the number density, number density ux, and energy density are given, respectively, by

**n(x, t) =**

∫
**f dΞ,****j(x, t) = m**

∫

**ξf dΞ,****E(x, t) = m**

∫ *ξ*^{2}
2**f dΞ.**

Other denitions of higher order moments such as internal energy density
**ε(x, t)**, stress tensor P*ij***(x, t)**, heat ux vector q*i** (x, t)* are also given,

**ε(x, t) =***m*
2

∫

**|ξ − u|**^{2}**f (ξ, x, t)dΞ**

*P*_{ij}**(x, t) = m**

∫

*(ξ*_{i}*− u**i**)(ξ*_{j}*− u**j***)f (ξ, x, t)dΞ***q*_{i}**(x, t) =***m*

2

∫

**|ξ − u|**^{2}*(ξ*_{i}*− u**i***)f (ξ, x, t)dΞ.**

We can have the familiar form of hydrodynamic equations,

*∂n*

*∂t* + *∂nu*_{i}

*∂x** _{i}* = 0 (2.16)

*ρdu*_{i}

*dt* + *∂P*_{ij}

*∂x*_{j}*= ρa**i* (2.17)

*ρdε*

*dt* *+ P*_{ij}*∂u*_{i}

*∂x** _{j}* =

*∂q*

_{i}*∂x*_{j}*,* (2.18)

where, _{dt}* ^{d}* =

_{∂t}

^{∂}*+ u*

_{i}

_{∂x}

^{∂}*i*. We have briey described the equilibrium distribu-
tion including MB, FD, and BE distributions. The next step is to derive the
corresponding hydrodynamic equations under these statistics. Although the

hydrodynamic equations are the same in both classical and semiclassical ver-
*sion, the transport coecients are dierent due to the dierences between f*^{eq,Q}*and f** ^{eq,C}*. After some straightforward algebraic manipulations, the equation

*of state of pressure p, the number density n and total energies E in classical*and quantum statistics as:

*p*_{c}*= nk*_{B}*T*
*n** _{c}* =

*z*

Λ^{3}
*E**c* = 3

2*nk**B**T +* *D*
2*mnu*^{2}*,*
and

*p*_{q}*= nk*_{B}*Tg*_{5/2}*(z)*
*g*_{3/2}*(z)*
*n**q* = *g*_{3/2}*(z)*

Λ^{3}
*E** _{q}* = 3

2*nk*_{B}*Tg*_{5/2}*(z)*
*g*_{3/2}*(z)*+ *D*

2*mnu*^{2}*.*

*g** _{ν}* is the generalized FD/BE function which is dened as

*g*_{ν}*(z)≡* 1
*Γ(ν)*

∫ _{∞}

0

*x*^{ν}^{−1}

*z*^{−1}*e*^{x}*+ ηdx,* (2.19)
*wherein Γ(ν) is the Gamma function. Notice that in classical limit, g**ν**(z)*
*will approach z no matter what the number ν is. We can check that the*
*pressure p**q**, number density n**q**and total energies E**q*will approach the classical
*counterpart in classical limit. We can also check the specic heat γ = C**p**/C** _{v}*,

*in ideal classical gas,γ*

*c*

*= 5/3*

*for monatomic gas, but γ*

*q*=

^{5}

_{3}

^{g}

^{5/2}

_{g}*2*

^{(z)g}

^{1/2}

^{(z)}*3/2**(z)*

*in ideal quantum gas, which depends on the fugacity z and not constant.*

*Obviously, γ**q**approaches γ**c*in classical limit. Moreover, more hydrodynamical

2.4. Semiclassical Hydrodynamic Equations 17

coecients like viscosity and thermal conductivity are given by

*µ**q* *= τ nk**B**Tg*_{5/2}*(z)*
*g*_{3/2}*(z)* *= τ P*
*κ*_{q}*= τ5k*_{B}

*2m*[*7g*_{7/2}*(z)*

*2g*_{3/2}*(z)−* *5g*_{5/2}*(z)*

*2g*_{3/2}^{2} *(z)]nk*_{B}*T,* (2.20)
and the classical coecients are also listed for comparing:

*µ*_{c}*= τ nk*_{B}*T = τ P*
*κ*_{c}*= τ5k*_{B}

*2mnk*_{B}*T = τ5k*_{B}

*2mP.* (2.21)

All the semiclassical coecients will approach their classical counterparts and
Prandtl number equals 1 in classical coecients. However, in semiclassical
*model, the Prandtl number depends on the fugacity z and is no more a con-*
stant number. Similar results could be found in linearized semiclassical Boltz-
mann equations [23]. Until now, the semiclassical Boltzmann Equation and
the semiclassical hydrodynamical equations have been introduced and com-
pared with their classical counterparts. In those descriptions there appears an
*important parameter z which is not shown in classical Boltzmann equation*
*or MB statistics. The physical meaning of the fugacity z is described below:*

*Recall number density in semiclassical representation n**q* = ^{g}^{3/2}_{Λ}3* ^{(z)}*, if we con-

*sider 0 < z < 1, then g*

*3/2*

*(z)*

*≃ z, we have z = n*

*q*

*(h*

^{2}

*/2πmk*

_{B}*T )*

^{3/2}*= n*

*Λ*

_{q}^{3}.

*From this equation, z can be interpreted as the ratio of Λ*

^{3}to average volume occupied by particles. In other words, it is the ratio of occupied length to

*particles thermal wavelength. When z is small, it means the order of spa-*tial dimension is larger than the thermal wavelength and we can neglect the

*degeneracy eect of particles (highly non-degenerate gas). For large z, the*degeneracy eect is important (highly degenerate gas) because the order of thermal wavelength is comparable to the spatial dimension. Moreover, the

*thermal wavelength is proportional to T*

*, that means when the number*

^{−1/2}density is in normal condition this eect will be obvious especially in the low
*temperature. In summary, when z → 0 the quantum distribution will coin-*
cide with the classical one, and the physical explanation is that the length
dimension of particle is larger than the particle de Broglie wavelength. The
*wave property will not be important. When z is considerably large, two length*
scales become comparable and one cannot omit the quantum eect anymore.

*So, we can think the fugacity z is the index of the degree of degeneracy. The*
*fugacity z has some restrictions in two dierent quantum distributions. In the*
*case of Boson, z should not exceed 1 because of the non-negative density, and*
*in the Fermion case there are no such restrictions on z.*

Finally, the actual correction values of the generalized Fermi function is shown
in Fig. 2.1. One nds the BE and FD curves overlap MB curve in z → 0
limit, that means the classical statistics could only work in the classical limit
*(z → 0) of quantum statistics. Moreover, z is restricted in 0 ≤ z ≤ 1 for Bose*
*gas, the function g**3/2**(z)* *increases monotonically with z and is bounded, its*
*largest value being g**3/2**(1) = 2.612*.

2.4. Semiclassical Hydrodynamic Equations 19

**Fugacity**

0.2 0.4 0.6 0.8

0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2

**Bose(2.5)**
**Bose(1.5)**
**Bose(2.5)/Bose(1.5)**
**(a)**

**Fugacity**

0.2 0.4 0.6 0.8

0.2 0.4 0.6 0.8 1

**Fermi(2.5)**
**Fermi(1.5)**
**Fermi(2.5)/Fermi(1.5)**
**(b)**

**Fugacity**

0.2 0.4 0.6 0.8

0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2

**Bose(2.5)/Bose(1.5)**
**Fermi(2.5)/Fermi(1.5)**
**(c)**

*Figure 2.1: Bose and Fermi functions in dierent fugacity z. (a) Bose Func-*
*tion g**1.5*^{Bose}*(z), g**2.5*^{Bose}*(z)*, and ^{g}_{g}^{2.5}^{Bose}^{Bose}^{(z)}

*1.5* *(z)**, (b) Fermi Function g*^{F ermi}*1.5* *(z),g**2.5*^{F ermi}*(z)*,
and ^{g}_{g}^{F ermi}^{2.5}^{F ermi}^{(z)}

*1.5* *(z)*, (c) Comparing Bose ^{g}_{g}^{Bose}^{2.5}^{Bose}^{(z)}

*1.5* *(z)*, and Fermi Function ^{g}_{g}^{F ermi}^{2.5}^{F ermi}^{(z)}

*1.5* *(z)*

### Chapter 3

## Semiclassical Lattice Boltzmann Method

### Contents

3.1 Overview . . . 21 3.2 The Derivations of Conventional LB method . . . 22 3.2.1 Time Discretization . . . 22 3.2.2 Space and Velocity Discretization . . . 23

3.3 The Derivations of Single Relaxation Time Semiclas- sical LB method . . . 28

3.3.1 Summary of Single Relaxation Time Semiclassical LB Method . . . 32 3.4 Generalized Semiclassical LB method . . . 35

### 3.1 Overview

In this chapter, the semiclassical LB method will be derived from semiclas- sical Boltzmann equation with BGK approximation (2.11). In general, there are two steps for deriving LB method from the continuous Boltzmann equa- tion. First, the time discretization is achieved by integrating the Boltzmann

equation along characteristic line. Second, the discretization of space is done by low Mach expansion or Hermite polynomials expansion. The derivations of semiclassical LB method is following the Hermite expansion procedures and the results are extended to multiple relaxation time version. The derived equations could be validated by reducing the SLB equation to classical one in classical limit.

### 3.2 The Derivations of Conventional LB method

Just like other numerical methods, discretization of time and space before simulating continuous models or equations on a discrete digital computer is necessary, the detailed procedures are listed below.

### 3.2.1 Time Discretization

Consider the semiclassical Boltzmann equation with BGK approximation (2.11), and neglect the forcing term here:

*∂*_{t}**f + ξ****· ∇f = R***BGK**(f ) =−*1

*τ[f− f** ^{eq}*] (3.1)
Rewrite the Boltzmann BGK Equation (3.1) in the form:

*D*_{t}*f +* 1
*τf =* 1

*τf** ^{eq}* (3.2)

*where D**t* *≡ ∂**t* **+ ξ***· ∇. Then, integrating (3.2) over a time step δ**t* along
characteristics:

**f (x + ξδ**_{t}**, ξ, t + δ**_{t}*) = e*^{−δ}^{t}^{/τ}* f (x, ξ, t) +* 1

*τe*

^{−δ}

^{t}

^{/τ}∫ _{δ}_{t}

0

*e*^{t}^{′}^{/τ}*f*^{eq}**(x + ξδ**_{t}^{′}**, ξ, t + δ**_{t}^{′}*)dt** ^{′}*
(3.3)

3.2. The Derivations of Conventional LB method 23

*By using Taylor expansion, and neglecting the high order terms O(δ**t*2), also
*with τ*^{∗}*= τ /δ** _{t}*, the standard LB equation is derived as:

**f (x + ξδ**_{t}**, ξ, t + δ*** _{t}*)

**− f(x, ξ, t) = R***BGK*

*(f ) =−*1

*τ*^{∗}**[f (x, ξ, t)**− f^{eq}* (x, ξ, t)],*
(3.4)
This equation indicates the basic procedures of LB method: streaming and
collision. Since we have got the time discretization LB equation, the follow-
ing processes will focus on discretizing the equilibrium distribution function

*f*

^{eq}*on velocity space .*

**(x, ξ, t)**### 3.2.2 Space and Velocity Discretization

In LB method, computational domain is discretized on a regular lattice with

xed grid spacing which makes this method simple. Although original LB method is empirically derived from lattice gas automata, however, modern LB method is derived from continuous Boltzmann equation and has a concrete physical interpretation compared with those come from empirical derivations.

The two dierent approaches to discretizing the computational domain are low mach number expansion and hermite polynomials expansion. In the following derivations, particles are in high temperature, low density and under classical statistics are assumed. That means, the distribution follows the classical MB statistics. And the MB equilibrium distribution function is:

*f** ^{eq,C}* =

*n*

*(2πθ)*^{D/2}*exp[−c*^{2}

*2θ*] (3.5)

*wherein n is the number density, D is the space dimension, θ = RT the*
**temperature and the peculiar velocity c = ξ − u. ξ is the molecular velocity*** and u the mean velocity. C means classical equilibrium distribution function*
here.

3.2.2.1 Low Mach Number Expansion

Expanding (3.5) by Taylor expansion, we can get the discrete equilibrium distribution function as:

*f*^{eq,C}*≈* *n*

*(2πθ)*^{D/2}*exp(−ξ*^{2}

*2θ*)[1 +**ξ****· u**

*θ* +**(ξ****· u)**^{2}
*2θ*^{2} *−* *u*^{2}

*2θ*] (3.6)
After choosing some suitable weighting and quadrature points to recover the
original moment integration, it will become the frequently used LB model. We
use D2Q9 model for example for presenting these procedures. First, we need
*to discretize the momentum space ξ properly. The integration is replaced by*
the summations as:

∫

**ψ(ξ)f**^{eq,C}* (x, ξ, t) =*∑

*a*

*W*_{a}**ψ(ξ**_{a}*)f*^{eq,C}**(x, ξ**_{a}*, t),* (3.7)
* where ψ(ξ) is the polynomial of ξ. Then, the above integral is evaluated by a*
Gaussian-type quadrature:

*I =*

∫

*exp(−ξ*^{2}

*2θ )ψ(ξ)dξ =*∑

*a*

*W*_{a}*exp(−ξ*^{2}

*2θ )ψ(ξ*

_{a}*),*(3.8) A cartesian coordinate system is used to recover the D2Q9 model. Therefore,

**we set ψ(ξ) = ξ***x*

^{m}*ξ*

^{n}

_{y}*, where ξ*

*x*

*and ξ*

*y*

**are the x and y components of ξ. Such***that, the integrand I will be*

*I = (√*

*θ)*^{(m+n+2)}*I*_{m}*I*_{n}*,* (3.9)

*where I**m* =∫_{+}_{∞}

*−∞* *e*^{−ς}*and ς = ξ**x**/√*

*2θ* *or ς = ξ**y**/√*

*2θ*. For D2Q9 model, the
*third-order Hermite formula for evaluating I**m* *is chosen, i.e. I**m* =∑3

*j=1**ω*_{j}*ς*_{j}* ^{m}*.
The three abscissas of the quadrature are:

*ς*_{1} =*−*√
*3/2*
*ς*_{2} = 0
*ς*_{3} =√

*3/2,* (3.10)

3.2. The Derivations of Conventional LB method 25

and the corresponding weight coecients are:

*ω*_{1} =*√*
*π/6*
*ω*2 = 2*√*

*π/3*
*ω*_{3} =*√*

*π/6.* (3.11)

After this integrating, the equilibrium distribution function of the D2Q9 model is:

*f*_{a}^{(0),C}*= ω*_{a}*n[1 + 3ξ*

_{a}

**· u***c*^{2} +**9(ξ**_{a}**· u)**^{2}

*2c*^{4} *−* **3u**^{2}

*2c*^{2}*],* (3.12)
*where c =√*

*3θ≡* ^{δx}_{δt}*. δx is the space between two adjacent lattices. Similarly,*
other models like D2Q7 and D3Q27 can be derived in a similar manner.

3.2.2.2 Hermite Polynomials Expansion

Hermite polynomials expanding of the distribution function for LB method was rstly presented in [17], then well developed in [47] and [16]. This method is inspired by [5], [18] and [48]. The brief descriptions are described below.

Following the approaches in [17][16][47], we adopt the Grad's moment ap- proach and seek solutions to (3.4) by expanding f(x, ξ, t) in terms of Hermite polynomials,

**f (x, ξ, t) = ω(ξ)**

∑*∞*
*n=0*

1

*n!***a**^{(n)}*(x, t)H*^{(n)}* (ξ)* (3.13)

**And the expansion coecients a**

*are given by*

^{(n)}**a**^{(n)}**(x, t) =**

∫

**f (x, ξ, t)**H^{(n)}* (ξ)dξ* (3.14)

**where ω(ξ) =**

_{(2π)}^{1}

*3/2*

*e*

^{−ξ}^{2}

^{/2}**is the weighting function, a**

^{(n)}*and H*

^{(n)}*are rank- n tensors and the product on the right-hand side denotes full contraction. Here and throughout the dissertation, the shorthand notations of Grad [5] for fully symmetric tensors have been adopted.*

**(ξ)**Some of the rst few tensor Hermite polynomials are given here,
*H*^{(0)}**(ξ) = 1**

*H*^{(1)}*i* **(ξ) = ξ**_{i}*H*^{(2)}*ij* **(ξ) = ξ**_{i}*ξ*_{j}*− δ**ij*

*H*^{(3)}_{ijk}**(ξ) = ξ**_{i}*ξ*_{j}*ξ*_{k}*− ξ**i**δ*_{jk}*− ξ**j**δ*_{ik}*− ξ**k**δ*_{ij}

The rst few expansion coecients can be easily identied with the familiar hydrodynamic variables:

**a**^{(0)} =

∫

**f dξ = n****a**^{(1)} =

∫

**f ξdξ = nu****a**^{(2)} =

∫

**f (ξ**^{2}**− δ)dξ = P + n(u**^{2}*− δ)*
**a**^{(3)} =

∫

**f (ξ**^{3}**− ξδ)dξ = Q + u(a**^{(2)}**− 2nu**^{2}) (3.15)
**Where the momentum ux tensor is dened as P = P***ij* *≡*∫

*f c*_{i}*c*_{j}* dc*, the

**heat ux tensor is dened as Q = g***ijk*

*≡*∫

*f c*_{i}*c*_{j}*c*_{k}**dc**, and We also have ε =

1

2**[a**^{(2)}_{ii}**−n(uu−3)]. The orthogonality of Hermite polynomials implies that the***leading moments of a distribution function up to the Nth order are preserved*
by truncations of the higher-order terms in its Hermite expansion. Thus, a
distribution function of the Boltzmann-BGK equation can be approximated by
*its projection onto a Hilbert space spanned by the rst N Hermite polynomials*
*without aecting the rst N moments. Here, up to Nth order, f*^{N}**(x, ξ, t)*** has exactly the same velocity moments as the original f(x, ξ, t) does. This*
guaranties that a LB method gas dynamic system can be constructed by a

nite set of macroscopic variables. As a partial sum of Hermite series with

*nite terms, the truncated distribution function f** ^{N}* can be completely and
uniquely determined by its values at a set of discrete abscissae in the velocity