Abstract.: A general lattice Boltzmann method for simulation of fluids with tailored transport coefficients is presented. It is based on the recently introduced quasi-equilibrium kinetic models, and a general lattice Boltzmann implementation is developed. Lattice Boltzmann models for isothermal binary mixtures with a given Schmidt number, and for a weakly compressible flow with a given Prandtl number are derived and validate

Ansumali, S.

Arcidiacono, S.

Chikatamarla, S. S.

Gorban, A. N.

Karlin, I. V.

Prasianakis, N. I.

RERO DOC Digital Library

Eur. Phys. J. B 56, 135–139 (2007)DOI: 10.1140/epjb/e2007-00100-1 THE EUROPEANPHYSICAL JOURNAL BQuasi-equilibrium lattice Boltzmann methodS. Ansumali1, S. Arcidiacono2, S.S. Chikatamarla3, N.I. Prasianakis3, A.N. Gorban4, and I.V. Karlin3,a1 CBE, Nanyang Technological University, 639798 Singapore, Singapore2 Paul Scherrer Institute, Combustion Research Unit, 5232 Villigen PSI, Switzerland3 Aerothermochemistry and Combustion Systems Lab, ETH Zurich, 8092 Zurich, Switzerland4 Department of Mathematics, University of Leicester, Leicester LE1 7RH, UKReceived 29 September 2006 / Received in ﬁnal form 2 February 2007Published online 12 April 2007 – c© EDP Sciences, Societa` Italiana di Fisica, Springer-Verlag 2007Abstract. A general lattice Boltzmann method for simulation of ﬂuids with tailored transport coeﬃcientsis presented. It is based on the recently introduced quasi-equilibrium kinetic models, and a general latticeBoltzmann implementation is developed. Lattice Boltzmann models for isothermal binary mixtures witha given Schmidt number, and for a weakly compressible ﬂow with a given Prandtl number are derived andvalidated.PACS. 47.11.-j Computational methods in ﬂuid dynamics – 05.20.Dd Kinetic theory1 IntroductionThe lattice Boltzmann (LB) method is a powerful ap-proach to hydrodynamics, with applications ranging fromlarge Reynolds number ﬂows to ﬂows at a micron scale,porous media and multi-phase ﬂows [1]. However, in spiteof its rapid development over a decade, the method is stilllacking a systematic and ﬂexible construction of numeri-cal schemes for situations beyond a simple ﬂuid, such asmixture models with a given Schmidt number, or thermalmodels with a given Prandtl number.In this paper, we introduce a general method of con-structing the lattice Boltzmann models. Our approach isbased on recently introduced quasi-equilibrium (QE) ki-netic models [2,3]. The structure of the paper is as fol-lows: ﬁrst, for the sake of completeness, we remind theconstruction of continuous time-space QE models. Afterthat, we derive a lattice Boltzmann discretization schemeof the QE models. The resulting quasi-equilibrium latticeBoltzmann method (QELBM) is illustrated with two ex-amples of particular interest. The ﬁrst example is the ther-mal model with a prescribed Prandtl number. The secondexample is the isothermal binary mixture model with aprescribed Schmidt number. Basic steps of constructionof these models are presented and numerical validation isprovided.a e-mail: karlin@lav.mavt.ethz.ch2 Quasi-equilibrium kinetic models2.1 GeneralWe address the class of lattice Boltzmann modelsequipped with the Boltzmann entropy function of theform,H =n∑i=1fi ln(fiWi), (1)where fi ≥ 0 are populations of the discrete velocities vi,i = 1, . . . , n, and Wi > 0 are corresponding weights. Awide class of relevant entropy functions (1) pertinent tosimulation of hydrodynamics was described in [4–8]. Inthe examples below we shall use the H-function for theisothermal model [4] and the recently introduced weaklycompressible thermal model [6] in two dimensions, forwhich the set of nine velocities vi and correspondingweights Wi arevx = {0, 1, 0,−1, 0,+1,−1,−1,+1}vy = {0, 0, 1, 0,−1,+1,+1,−1,−1}W = (1/36) {16, 4, 4, 4, 4, 1, 1, 1, 1} . (2)We begin with a generic construction of the quasi-equilibrium kinetic models [2,3] speciﬁed for the discretevelocity case. Let M = {M1, . . . ,MkM } be a set of lo-cally conserved ﬁelds, and N = {N1, . . . , NkN} be a setof quasi-conserved slow ﬁelds. Functions Mm and Nl areassumed linear functions (moments) of the populations,136 The European Physical Journal BMm =∑i Mmifi and Nl =∑i Nlifi. The choice of Mand N depends on the particular problem but we al-ways assume that the density is included into the list M ,ρ = M1 =∑ni=1 fi. The equilibrium population vectorf eq(M) is deﬁned as the minimizer of the H-function (1)under ﬁxed M . The quasi-equilibrium f∗(M,N) is deﬁnedas the minimizer of H under ﬁxed M and N . By construc-tion, functions f∗ and f eq satisfy consistency relationsM(f∗(M,N)) = M,N(f∗(M,N)) = N,M(f eq(M)) = M. (3)The quasi-equilibrium kinetic model reads∂tfi + viα∂αfi = − 1τ1[fi − f∗i (M(f), N(f))]− 1τ2[f∗i (M(f), N(f))− f eqi (M(f))] . (4)Denoting the right hand side (the collision integral) as Qi,it is easy to see that consistency condition implies localconservation laws, M(Q) = 0. Note that the ﬁrst partof the collision integral which describes relaxation to thequasi-equilibrium, also conserves the N -ﬁelds, N(f−f∗) =0. Moreover, it is straightforward to prove the H-theorem.For that, it suﬃces to rewriteQi = − 1τ2[fi − f eqi ]−(τ2 − τ1τ1τ2)[fi − f∗i ]. (5)The entropy production σ =∑ni=1(∂H/∂fi)Qi becomesσ = − 1τ2n∑i=1ln(fif eqi)(fi − f eqi )−(τ2 − τ1τ1τ2) n∑i=1ln(fif∗i)(fi − f∗i ), (6)and is non-positive semi-deﬁnite provided relaxation timessatisfy the hierarchy,τ1 ≤ τ2. (7)Thus, in the QE models, relaxation to the equilibrium issplit in two steps. In the ﬁrst step, the distribution func-tions relaxes to the quasi-equilibrium with the faster re-laxation time τ1. In the second step, the quasi-equilibriumrelaxes to the equilibrium with the slower relaxation timesτ2. If τ1 = τ2, the intermediate relaxation step to thequasi-equilibrium disappears from (4), and it reduces tothe Bhatnagar-Gross-Krook (BGK) model.2.2 Triangle entropy methodFor a practical implementation, explicit form of the func-tions f eqi (M) and f∗i (M,N) are required. While for mostof the cases, the equilibrium f eqi can be found explicitlyeither in a closed form or in a form of expansion, explicitconstruction of the quasi-equilibrium is case-dependent.Here we suggest a simple way to ﬁnd quasi-equilibria inexplicit form by perturbation around the equilibrium [2,9].Let us assume that the equilibrium f eq(M) has been foundexplicitly, and that near the equilibriumf∗ = f eq(M) + δf∗(M,N). (8)Expanding the entropy function H (1) at the equilibriumto quadratic terms, we obtainH(f eq + δf) = H(f eq) +n∑i=1δfi(ln(f eqiWi)+ 1)+12n∑i=1δf2if eqi+ O(δf3). (9)The quadratic expansion of the entropy function (9) main-tains convexity, and we ﬁnd its minimizer δf∗ subject tothe linear constraintsn∑i=1Mkiδf∗i = 0,n∑i=1Nliδf∗i = Nl −N eql , (10)where N eq = N(f eq(M)) are the values of the non-conserved ﬁelds at equilibrium. Solution to this minimiza-tion problem has the form,f∗ = f eqi (M) + δf∗i= f eqi (M)(1 +kM∑s=1λsMsi +kN∑s=1χsNsi), (11)where Lagrange multipliers λs and χs are found upon sub-stituting function (11) into the constrains (10), and solv-ing the (kM + kN )× (kM + kN ) linear algebraic problem.Introducing a (kM + kN )-dimensional vector of Lagrangemultipliers, (λ, χ) = (λ1, . . . , λkM , χ1, . . . , χkN ), and ma-trices(AMM )kl =n∑i=1Mkifeqi Mli, k, l = 1, . . . , kM(AMN )kl =n∑i=1Mkifeqi Nli, k = 1, . . . , kM , l = 1, . . . , kN(ANN )kl =n∑i=1Nkifeqi Nli, k, l = 1, . . . , kN , (12)we ﬁnd the solution in the matrix form,(λχ)=(AMM AMNATMN ANN)−1 (0N −N eq), (13)where T denotes transposition, and N − N eq is thekN -dimensional vector (N1 −N eq1 , . . . , NkN −N eqkN ). Notethat the solution depends linearly on the deviation of thenon-conserved ﬁelds N −N eq and non-linearly on the lo-cally conserved M . For this quasi-linear quasi-equilibrium,S. Ansumali et al.: Quasi-equilibrium lattice Boltzmann method 137the entropy production becomesσ = − 1τ2n∑i=1(fi − f eqi )2f eqi−(τ2 − τ1τ1τ2) n∑i=1(fi − f∗i )2f eqi+ O(δf3). (14)Thus, with the use of the triangle entropy method, thekinetic model satisﬁes the entropy production inequality(both of the two quadratic forms in (14) are non-positivesemi-deﬁnite) once the populations remain close to thelocal equilibrium. This is suﬃcient for most of the appli-cations.3 Quasi-equilibrium lattice BoltzmannmethodWe shall now derive a second-order time discretization forthe generic kinetic equation (4). Since the derivation onlyuses the consistency condition (3), it is equally applica-ble to exact quasi-equilibria and to those obtained by thetriangle entropy method (8). Following [10], kinetic equa-tions (4) are integrated in time from t to t+δt along char-acteristics, and the time integral of the right hand side isevaluated by trapezoidal rule. Introducing a mapfi → gi = fi − δt2 Qi(f), (15)the result is written asgi(x + ciδt, t + δt) = gi(x, t) − ω1[gi(x, t)− f∗i (x, t)]−ω1τ1τ2[f∗i (x, t)− f eqi (x, t)], (16)whereω1 =2δt2τ1 + δt,f∗(x, t) = f∗(M(f), N(f)),f eqi (x, t) = feqi (M(f)). (17)Note that the right hand side in (16) still contains termswhich depend on populations f rather than on the func-tion g (15). In order to obtain a closed-form equationsolely in terms of functions g, we ﬁrst remark that, takinglocally conserved moments of the map (15) we have (sameas in the BGK case) M(f) = M(g), thus f eq(f) = f eq(g)in the second term in (16). For the quasi-conserved ﬁeldsN , the situation is slightly diﬀerent. Evaluating the mo-ments N of the map (15), we obtainN(g) = N(f) +δt2τ2(N(f)−N eq(f)). (18)Inverting this relation, and substituting it into (16), aftersimple transformations we obtain:gi(x + ciδt, t + δt) = (1− ω1)gi + ω1(τ1τ2)f eqi (M)+ω1(τ2 − τ1τ2)f∗i (M,N′), (19)where M = M(g), while N ′ is evaluated according to theruleN ′ =(1− δt2τ2 + δt)N(g) +δt2τ2 + δtN eq(g). (20)Equation (19) is the basic second-order time stepping algo-rithm for the quasi-equilibrium lattice Boltzmann models,and is the main result of this paper. It is important to notethe shift in the dependence of the quasi-equilibrium pop-ulation f∗, it does not depend just on N(g) but ratheron a convex linear combination between N(g) and theequilibrium value N eq(g). Discretization in space dependson the problem at hand (in the simplest case, the latticeBoltzmann discretization is readily applicable). We shallnow proceed with speciﬁc examples of QELBM.4 Examples4.1 One-component fluid with a given Prandtl numberNavier-Stokes equations for a one-component compress-ible ﬂuid are characterized by the Prandtl number, Pr =(Cpµ)/κ, where Cp is the speciﬁc heat at constant pressure(Cp = (D + 2)/2 is speciﬁc heat of ideal gas below), µ isviscosity coeﬃcient, and κ is thermal conductivity. Whena single relaxation time kinetic equation is used (for ex-ample, the BGK model), this results in a ﬁxed Prandtlnumber. As we shall see it below, physical consistencyof the QE models, the relaxation time hierarchy (7), im-plies that two models with two diﬀerent quasi-equilibriaare required in order to cover the entire range of Prandtlnumber. For the sake of concreteness, we consider a sim-ple, weakly compressible ﬂow model on the standard two-dimensional D2Q9 lattice [6,13]. The locally conservedﬁelds M are density ρ, momentum density j, and pres-sure density p,n∑i=1{1, vi, v2i }fi = {ρ, j, 2p + ρ−1j2}. (21)To second order in the momentum, the equilibriumreads [6,13]:f eqi (ρ, j, p) = ρ(1− pρ)2 ( pρ2(1− pρ))v2i [1 +viαjαp+jαjβ2p2(viαviβ −4(pρ)2 + v2i (1 − 3 pρ )2(1− pρ )δαβ)]. (22)This model operates in a small temperature window T =p/ρ, around the reference temperature T0 = 1/3 (seeRefs. [6,13] for an estimate of the temperature window).We remark that, due to a low symmetry of the velocity set,the lattice BGK (LBGK) model with the equilibrium (22)gives PrLBGK = 4 [13] (not Pr = 1, as in the continuouskinetic theory or thermal lattice Boltzmann models witha higher symmetry [5,8]). In the ﬁrst case (Pr < PrLBGK),138 The European Physical Journal Bthe speciﬁed slow variables N are the components of theheat ﬂux q,qα =n∑i=1(viα − uα) (vi − u)2 fi, (23)where u = j/ρ is the mean velocity. In the opposite case(Pr > PrLBGK), the slow variables N are the componentsof the stress tensor ΘΘαβ =n∑i=1[(viα − uα) (viβ − uβ)− 2Dδαβ(vi − u)2]fi.(24)Explicit expressions for the pertinent quasi-equilibria,f∗i (ρ, j, p, q) (case (23), Pr ≤ PrLBGK) and f∗i (ρ, j, p,Θ)(case (24), Pr ≥ PrLBGK) are easily found using theexplicit formulas of the triangle entropy method, equa-tions (11) and (13), and are not displayed here. Using theChapman-Enskog method, we derive the hydrodynamicequations (Navier-Stokes-Fourier) for the density, momen-tum and temperature from the continuous QE models (4)with the corresponding quasi-equilibria f∗i (ρ, j, p, q) andf∗i (ρ, j, p,Θ) (note that the explicit form of f∗i is not re-quired to perform this analytical computation, the consis-tency condition 3 suﬃces). The viscosity µ and the ther-mal conductivity κ coeﬃcients thus obtained, imply thefollowing Prandtl number:µ = ρT0τ1, κ =12ρT0τ2 ⇒ Pr = 4τ1τ2≤ PrLBGK (23),(25)µ = ρT0τ2, κ =12ρT0τ1 ⇒ Pr = 4τ2τ1≥ PrLBGK (24).(26)A comment to these formulas is in order: the hierarchyof relaxation times (τ1 ≤ τ2) implies a restriction onthe range of admissible Prandtl number when a speciﬁedquasi-equilibrium is used. In the case where the heat ﬂuxis considered as the slow variable, the thermal conductiv-ity (relaxation rate of the heat ﬂux) is proportional to theslow relaxation time, κ ∼ τ2, while µ ∼ τ1. In the oppositecase when the stress tensor is chosen as a slow variable,the dependence is inverted, µ ∼ τ2 and κ ∼ τ1.Implementation of the present QE models is a straight-forward application of equation (19). In Figure 1, simula-tion of Couette ﬂow between parallel plates at diﬀerenttemperatures is compared with the analytical solution.Diﬀusive boundary condition was used [14]. Agreementbetween simulation and analytical solution is excellent.4.2 Binary mixture at a given Schmidt numberOur second example is the isothermal binary mixture ofideal ﬂuids, A and B, with particles masses mA and mB.The locally conserved variables M are the densities ofthe components, ρA,B, and the momentum of the mix-ture j = jA + jB, where jA,B are the momenta of the0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 100.511.522.5xT NormFig. 1. Steady state temperature proﬁle in Couette ﬂow be-tween parallel plates moving with a relative velocity U , at thetemperature diﬀerence δT . Lower curve: QELBM (23) withPr = 0.710 (air); middle curve: LBGK model with Pr = 4; up-per curve: QELBM (24) with Pr = 8. Symbol: simulation; line:analytical solution. The Eckert number, Ec = (U2)/(CpδT )is 1.5.components. The eﬃciency of diﬀusion mixing versus vis-cous dissipation of momentum in various ﬂuids is charac-terized by the Schmidt number, Sc = µ/(ρDAB), whereµ is the viscosity coeﬃcient, DAB is the binary diﬀu-sion coeﬃcient, and ρ = ρA + ρB. Let us further intro-duce the molar fractions XA,B = nA,B/(nA + nB) wherenA,B = ρA,B/mA,B is the number density of components,and the reference Schmidt number, S˜c = mAB/(ρXAXB),with mAB = (ρAρB)/(ρA + ρB) the reduced mass density.We shall now consider a pair of quasi-equilibrium modelswith two diﬀerent quasi-equilibria which cover the entirerange of Schmidt number, Sc ≶ S˜c. For concreteness, weuse the standard isothermal D2Q9 model. For the single-component case, the equilibrium f eqi (ρ, j) to second orderin momentum is obtained by setting p/ρ = T0 = 1/3 inequation (22). The equilibrium of the mixture is then de-scribed by the populations f eqA,Bi = feqi (ρA,B, j).The case Sc < S˜c has been already considered re-cently [11,12], and is mentioned here for the sake of com-pleteness. In this case we consider the momenta of thecomponents, jA,B as the quasi-conserved variables N ,and the corresponding quasi-equilibrium is immediatelyread oﬀ the equilibrium of the single-component ﬂuid;to second order in momentum it is f∗A,Bi(ρA,B, jA,B) =f eqi (ρA,B, jA,B). The continuous QE model recovers equa-tions of the isothermal hydrodynamics that is, theadvection-diﬀusion equations for the densities ρA,B andthe Navier-Stokes equation for the momentum j withµ = τ1p0, DAB = τ2p0 XAXBmAB , where p0 = T0(nA + nB)is the pressure at the reference temperature T0. Thismodel recovers Sc = (τ1/τ2)S˜c, and according to theS. Ansumali et al.: Quasi-equilibrium lattice Boltzmann method 139Fig. 2. Binary diﬀusion for the case Sc > S˜c and mB/mA =500 at diﬀerent time steps t. Symbol: simulation; line: analyt-ical solution. Triangles: t = 500; circles: t = 3000; squares:t = 6000; diamonds: t = 9000.hierarchy of the relaxation times (7), pertains to the ﬂuidswith Sc ≤ S˜c.In the opposite case, Sc > S˜c, we choose the stresstensors of the components, ΘA,B, as the quasi-conservedﬁelds N (see Eq. (24)). Corresponding quasi-equilibria,f∗A,Bi(ρA,B, j,ΘA,B), are constructed using the triangleentropy method in the same way as the case Pr > PrLBGK.The lattice Boltzmann implementation is based on (19)and an interpolation step, as explained in [11,12]. In Fig-ure 2, we present a simulation of diﬀusion of two ﬂuidswith a high mass ratio mB = 500mA in a setup where ini-tially a half-space is ﬁlled with the 90–10% mixture andthe other half-space — with the 10–90% mixture. Agree-ment with the analytical solution is excellent.5 ConclusionIn conclusion, we suggested a systematic, physically trans-parent and realizable approach to constructing latticeBoltzmann models for hydrodynamic systems. All themodels considered herein require only the choice of thequasi-equilibria appropriate to the physical context of theproblem. Following the same pattern, it is straightforwardto construct kinetic models for bulk viscosity and chemicalreactions.References1. S. Succi, The Lattice Boltzmann Equation For FluidDynamics and Beyond (Oxford University Press, Oxford,2001)2. A.N. Gorban, I.V. Karlin, Physica A 206, 401 (1994)3. C.D. Levermore, J. Stat. Phys. 83, 1021 (1996)4. I.V. Karlin, A. Ferrante, H.C. O¨ttinger, Europhys. Lett.47, 182 (1999)5. S. Ansumai, I.V. Karlin, H.C. O¨ttinger, Europhys. Lett.63, 798 (2003)6. S. Ansumali, I.V. Karlin, Phys. Rev. Lett. 95, 260605(2005)7. S.S. Chikatamarla, S. Ansumali, I.V. Karlin, Phys. Rev.Lett. 97, 010201 (2006)8. S.S. Chikatamarla, I.V. Karlin, Phys. Rev. Lett. 97,190601 (2006)9. A.N. Gorban, I.V. Karlin, Physica A 360, 325 (2006)10. X. He, S. Chen, G.D. Doolen, J. Comput. Phys. 146, 282(1998)11. S. Arcidiacono, S. Ansumali, I.V. Karlin, J. Mantzaras,K.B. Boulouchos, Mathematics and Computers inSimulation 72, 79 (2006)12. S. Arcidiacono, J. Mantzaras, S. Ansumali, I.V. Karlin,C. Frouzakis, K.B. Boulouchos, Phys. Rev. E 74, 056707(2006)13. N.I. Prasianakis, S.S. Chikatamarla, I.V. Karlin, S.Ansumali, K.B. Boulouchos, Mathematics and Computersin Simulation 72, 179 (2006)14. S. Ansumali, I.V. Karlin, Phys. Rev. E 66, 026311 (2002)

Quasi-equilibrium lattice Boltzmann method

http://doc.rero.ch/record/319401/files/10051_2006_Article_100.pdf

Quasi-equilibrium lattice Boltzmann method

Abstract

Similar works

Full text

Available Versions

RERO DOC Digital Library