3.3 The Derivations of Single Relaxation Time Semiclassical LB
3.3.1 Summary of Single Relaxation Time Semiclassical LB
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:
∂tf + ξ· ∇f = RBGK(f ) =−1
τ[f− feq] (3.1) Rewrite the Boltzmann BGK Equation (3.1) in the form:
Dtf + 1 τf = 1
τfeq (3.2)
where Dt ≡ ∂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
et′/τfeq(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(δt2), also with τ∗ = τ /δt, the standard LB equation is derived as:
f (x + ξδt, ξ, t + δt)− f(x, ξ, t) = RBGK(f ) =−1
τ∗[f (x, ξ, t)− feq(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 discretizfollow-ing the equilibrium distribution function feq(x, ξ, t) on velocity space .
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:
feq,C = n
(2πθ)D/2exp[−c2
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: 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:
∫
ψ(ξ)feq,C(x, ξ, t) =∑
a
Waψ(ξa)feq,C(x, ξa, t), (3.7) where ψ(ξ) is the polynomial of ξ. Then, the above integral is evaluated by a Gaussian-type quadrature: A cartesian coordinate system is used to recover the D2Q9 model. Therefore, we set ψ(ξ) = ξxmξny, where ξx and ξy are the x and y components of ξ. Such that, the integrand I will be
I = (√ third-order Hermite formula for evaluating Im is chosen, i.e. Im =∑3
j=1ωjςjm. 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 coe cients are:
ω1 =√
After this integrating, the equilibrium distribution function of the D2Q9 model is:
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, And the expansion coe cients a(n) are given by
a(n)(x, t) =
∫
f (x, ξ, t)H(n)(ξ)dξ (3.14) where ω(ξ) = (2π)13/2e−ξ2/2is the weighting function, a(n)and H(n)(ξ)are rank-n terank-nsors arank-nd the product orank-n the right-harank-nd side derank-notes full corank-ntractiorank-n. 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 coe cients can be easily identi ed with the familiar hydrodynamic variables:
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, fN(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 fN can be completely and uniquely determined by its values at a set of discrete abscissae in the velocity
3.2. The Derivations of Conventional LB method 27
space. This is possible because with f truncated to order N, the integrand on the right-hand side of (3.14) can be expressed as:
fN(x, ξ, t)H(n)(ξ) = ω(ξ)q(x, ξ, t) (3.16)
where q(x, ξ, t) is a polynomial in ξ of a degree no greater than 2N. Using the Gauss-Hermite quadrature, a(n) can be precisely calculated as a weighted sum of functional values of q(x, ξ, t):
a(n)(x, t) = a Gauss-Hermite quadrature of degree ≥ 2N. Thus, fN is completely de-termined by the set of discrete functional values, fN(x, ξa, t); a = 1, ..., l, and therefore its rst N velocity moments, and vise versa. The set of discrete dis-tribution functions fN(x, ξa, t)now serve as a new set of fundamental variables (in physical space) for de ning the uid system in place of the conventional hydrodynamic variables.
Next, we expand the equilibrium distribution feq,C in the Hermite polyno-mial basis to the same order as fN, i.e., feq,C(x, ξ, t) ≈ f(0),C,N(x, ξ, t), and
These coe cients a(n)0 can be evaluated exactly and we have
a(0)0 = n, a(1)0 = nu, a(2)0 = n[uu + (θ− 1)δ],
a(3)0 = n[uuu + (θ− 1)δu] (3.20) where n,u and θ are in non-dimensional form hereinafter.
Denote fa(0),C ≡ wafeq,C(ξa)/ω(ξa) and for N = 3, we get the explicit Hermite expansion of the MB distribution at the discrete velocity ξa as:
fa(0),C,3 = wan{1 + ξa· u + 1
2[(u· ξa)2− u2+ (θ− 1)(ξa2− D)]
+ξ· u
6 [(u· ξa)2 − 3u2+ 3(θ− 1)(ξ2a− D − 2)]}. (3.21) where D = δii. We could check if θ = 1 which means isothermal here. Then (3.21) up to second order will be as the same as (3.6).
3.3 The Derivations of Single Relaxation Time Semiclassical LB method
The left hand side of UUB-BGK equation is as the same as the conventional Boltzmann-BGK equation (3.1). The only dierence of this two equation is the equilibrium distribution function on the right hand side, that means, the classical equilibrium distribution function feq,C (3.5) is replaced by the quan-tum equilibrium distribution function feq,Q (3.22) in UUB-BGK equation.
The quantum equilibrium distribution function which includes BE equilib-rium distribution function and FD equilibequilib-rium distribution function has been derived as (2.12) and rewritten here as:
feq,Q={z−1exp[(ξ− u)2
2θ ]− η}−1, (3.22)
3.3. The Derivations of Single Relaxation Time Semiclassical LB
method 29
First, The macroscopic variables is de ned as
Φ(x, t) =
∫ m3dξ
h3 ϕf. (3.23)
For examples, the number density
n(r, t) = generalized FD/BE function which has been de ned in Ch2 (2.19) together with the calculation approximations is listed as
gν(z)≡ 1 Notice that the series is valid only when z < 1. In deriving semiclassical LB method, the particle velocity dimension d and space dimension D do not need to be the same, such that we could have the more generalized results.
Following the same procedures applied to classical gases in [16] and already brie y described above, de ning ζ = mhddξ, and dζ = mddξ/hd, such that the coe cients a(n)0 which is de ned as (3.14) of the semiclassical expansion can be calculated as
Notice that the coe cients a(N )0 here is not yet non-dimensional, and δ has the same dimension as uu, ˆδ represents delta without dimension. To
nondi-mensionalize those coe cients we choose some reference variables as:
and notice the non-dimensional coe cients are:
ˆ
Combine (3.27) and (3.32), the non-dimensional coe cients are derived as below:
3.3. The Derivations of Single Relaxation Time Semiclassical LB Finally, we take oˆon every variables and all the non-dimensional coe -cients are listed as:
Put (3.37) into (3.13) and notice the variable ξ has been replaced by ζ now, we can get the expansion of the quantum equilibrium distribution function.
fa(0),Q,3 = ωan{1 + ζa· u + 1
To validate this result, we note that in classical limit z << 1, the general-ized FD function gν(z)will approximately equal z no matter what the order ν is. That means, gd/2+1g (z)
d/2(z) = 1in (3.38), such that, the semiclassical equilibrium
distribution function (3.38) in classical limit will become Comparing this result with (3.21), when d equals D, (3.21) and (3.39) are the same in the classical limit.
3.3.1 Summary of Single Relaxation Time Semiclassical LB Method
Once we have obtained fN and feq,Q,N at the discrete velocity abscissae ζa, we are ready to set up the whole system for solving (3.1) in the physical con guration space. We use (3.4) for solving (3.1) in the following simulations.
For clearly (3.4) is listed below.
fa(x + ζa, t + 1)− fa(x, t) =−1
τ[fa− fa(0),Q]. (3.40) where fa(0),Q is given by (3.38) and τ could be calculated from (2.20). Ap-plying Gauss-Hermite quadrature to the moment integration, we have the macroscopic quantities, the number density, number density ux, and energy density. And the macroscopic variables become:
n(x, t) = In summary, (3.40) and (3.41) form a closed set of dierential equations gov-erning the set of variables fa(x, t)in the physical con guration space. All the
3.3. The Derivations of Single Relaxation Time Semiclassical LB
method 33
macroscopic variables and their uxes can be calculated directly from their corresponding moment summations.
We should note that there appears the fugacity z in this semiclassical LB method that not shown in classical LB method, such that, we need some addi-tional procedures to calculate the fugacity z. The detail calculating procedures are shown below. First, the density can be calculated as
∑
a
fa = n (3.42)
and the velocity is calculated as
∑
a
faζa= nua (3.43)
Then, recall the two formulas:
n = θd/2gd/2 (3.44)
∑
a
faζa2 = ndθgd/2+1
gd/2 + nu2 = 2E, (3.45) where E is the total energy. Combining (3.44) and (3.45), we can get
2E− d( n gd/2
)2+dd gd/2+1− nu2 = 0 (3.46)
By solving(3.46) with Newton method, we can get the fugacity z. then we can get temperature θ from (3.44). Since the derivation is based on expanding distribution function f onto Hermite polynomials as [16], the lattices and corresponding weights listed in [16] also could be used. The D2Q9 lattice model for two dimensional ow and D3Q19 lattice model for three dimensional
ow are listed in Fig. 3.3, Table 3.2 and Fig. 3.4, Table 3.3. The subscript f s denotes a fully symmetric set of points.
The present semiclassical LB method is derived from expanding distribu-tion funcdistribu-tion onto Hermite polynomials. The accuracy of semiclassical LB method should be carefully discussed. According to [16], the second order N = 2 kinetic theory of hydrodynamics is necessary for representing n, u, θ and the momentum ux tensor P . On the other hand, to describe the dynam-ics of the internal energy of a uid system, the third-order N = 3 Hermite terms must be retained. The most important principle we should notice is that the higher order N we expand the equilibrium distribution function f, the higher quadrature degree n of lattice structure we should use. The rela-tion between them is n ≥ 2N. Such that, an accurate Navier-Stokes level of description for an isothermal momentum equation requires a quadrature of a degree of precision greater than 6. And we should include the Hermite expan-sion up to N ≥ 4, and the quadrature degree up to n ≥ 8 to ensure accuracy up to the Burnett order for isothermal uids. See the details in [16]. Since the semiclassical LB method is also derived from expanding the equilibrium distribution function onto Hermite polynomials, the principle of n ≥ 2N may be followed. But, we can not forget in deriving both classical LB method and semiclassical LB method, there are still higher order terms have been drop o, and the total in uences of dierent order method are worth studying. A two dimensional cylinder ow is considered here for studying this issue of accuracy.
we consider a uniform two-dimensional viscous ow over a circular cylinder in a BE quantum gas to illustrate the present semiclassical LB method in dierent order on D2Q9 lattices. N = 2 and N = 3 expansion equations are both set for this case. The computation domain is (−1, 1) × (−1, 1) and set by 201 × 201 lattices, and the cylinder is set at the center of the com-putation domain with the radius D = 0.1. Uniform Cartesian grid system is used. The free stream velocity is u∞= 0.1, free stream temperature θ∞ = 0.5