数理解析研究所講究録 1413 巻 2005 年 193-200

193

Finite-Difference Lattice Boltzmann Methods for Binary Fluids
Aiguo Xu

Department of Physics, Yoshida-South Campus,
Kyoto University, SakyO-ku, Kyoto, 606-8501, Japan

In this proceeding we summarize our recent studies on two fluid lattice Boltzmann methods forbinary fluids. We first clearify Sirovich’s kinetic theory, then based on which three multispeed
discrete velocity models are formulated, which are for the Euler equations, isothermal Navier-Stokes
equations and the complete Navier-Stokes equations, respectively. Each formulated discrete velocity
model, together with an appropriate finite-difference scheme, composes afinite-difference lattice
Boltzmann method. The validity of the methods is verified by investigating (i) the Couette flowand (ii) the uniform relaxation process of the two components.
,,
PACS numbers:$47.11.+\mathrm{j}$$51.10.+\mathrm{y}$$05.20.\mathrm{D}\mathrm{d}$I.

INTRODUCTION

Lattice Boltzmann Method (LBM) is anumerical scheme to simulate kinetic systems. The
LBM recovers the hydrodynamic descriptions in the small Knudsen number limit. It has become
aviable and promising numerical scheme for simulating fluid flows. There are several options todiscretize the Boltzmann equation: (i) Standard LBM (SLBM)[I]; (ii) Finite-Difference LBM(LBM) [1-3]; (iii) Finite-Volume LBM $[1, 4]$ ;(iv) finite Element LBM[I, 5]; etc. These kindsof schemes are expected to be complementary in the LBM studies.
Even though various LBMs for multicomponent fluids 8-20] have been proposed and develbased on the single-fluid
oped , (i) most existing methods belong to the SLBM[6, 8-17,
theory[8-15, 17, 18, 21]; (ii) in Refs. $[6, 7]$ two SLBMs are proposed, but these two models are
not convenient (if not impossible) to simulate thermal and compressible systems, even isothermal and incompressible systems only if the two components have different particle masses. In
this study we develop two fluid FDLBMs for thermal and compressible binary fluids.
$\mathrm{a}\mathrm{n}\mathrm{d}/\mathrm{o}\mathrm{r}$II.

FORMULATION AND VERIFICATION OF THE FDLBMSThe formulation of aFDLBM consists of three steps: (i) select or design an appropriate
discrete velocity model (DVM), (ii) formulate the discrete local equilibrium distribution function, (iii) choose afinite-difference scheme. The continuous Boltzmann equation has infinite
velocities, so the rotational invariance is automatically satisfied. Recovering rotational invariant
macroscopic equations from adiscrete finite velocity microscopic dynamics imposes constraints
on the isotropy of DVM used. In our studies, the proposed FDLBMs are based on the two
DVMs described below.
DVMI$:\mathrm{v}_{0}=0$,$\mathrm{v}_{\mathrm{k}1}=\mathrm{v}_{\mathrm{k}}[\cos(\frac{\mathrm{i}\pi}{6}),$ $\sin(\frac{\mathrm{i}\pi}{6})],\mathrm{i}=1,2$, \cdots ,12,(1)

where indicates the -th group of particle velocities and indicates the direction of the particle
speed. It is easy find that (i) its odd rank tensors are zero, and (ii) its initial four even rank
tensors satisfy
$i$

$k$

$k$

$\sum_{0=1}^{12}v_{ki\alpha}v_{ki\beta}=6v_{k}^{2}\delta_{\alpha\beta}$,$\sum_{i=1}^{12}v_{ki}\alpha kv\dot{\iota}\beta v_{ki\gamma}v_{ki\delta}$$\sum_{i=1}^{12}vk\dot{\iota}\alpha vki\beta v_{k:}v\cdot v_{ki\mu}v_{k\nu}\gamma b\delta|.=\not\supset^{v_{k}\Delta_{\alpha\beta\gamma\delta\mu\nu}}16$$= \frac{3}{2}v_{k}^{4}\Delta_{\alpha\beta\gamma\delta}$(2),$\sum_{\dot{l}=1}^{12}v_{k}|\alpha vk_{\dot{l}}\beta vki\gamma vk_{\dot{l}}\delta v_{k\dot{\cdot}\mu}v_{k:\nu}v_{k\lambda}|.v_{ki\pi}=\frac{1}{32}v^{8}k\Delta_{\alpha\beta\gamma\delta\mu\nu\lambda\pi}$where , ,$\alpha$$\beta$$\cdots$

indicate

$x$

or

$y$

,

,

component and$\Delta_{\alpha\beta\gamma\delta}=\delta_{\alpha\beta}\delta_{\gamma\delta}+\delta_{\alpha\gamma}\delta_{\beta\delta}+\delta_{\alpha\delta}\delta_{\beta\gamma}$

,

(3)

194

,

(4)

$\Delta_{\alpha\beta\gamma\delta\mu\nu}=\delta_{\alpha\beta}\Delta_{\gamma\delta\mu\nu}+\delta_{\alpha\gamma}\Delta_{\beta\delta\mu\nu}+\delta_{\alpha\delta}\Delta_{\beta\gamma\mu\nu}+\delta_{\alpha\mu}\Delta_{\beta\gamma\delta\nu}+\delta_{\alpha\nu}\Delta_{\beta\gamma\delta\mu}$$\Delta_{\alpha\beta\gamma\delta\mu\nu\lambda\pi}=\delta_{\alpha\beta}\Delta_{\gamma\delta\mu\nu\lambda\pi}+\delta_{\alpha\gamma}\Delta_{\beta \mathit{5}\mu\nu\lambda\pi}+\delta_{\alpha\delta}\Delta_{\beta\gamma\mu\nu\lambda\pi}+\delta_{\alpha\mu}\Delta_{\beta\gamma\delta\nu\lambda\pi}$.

(5)$+\delta_{\alpha\nu}\Delta_{\beta\gamma\delta\mu\lambda\pi}+\delta_{a\lambda}\Delta_{\beta\gamma\delta\mu\nu\pi}+\delta_{\alpha\pi}\Delta_{\beta\gamma\delta\mu\nu\lambda}$

It is clear that this DVM is isotropic up to, at least, its 9th rank tensor.
DVM2 :$\mathrm{v}_{0}=0$,$\mathrm{v}_{\mathrm{k}1}=\mathrm{v}_{\mathrm{k}}[\mathrm{c}\mathrm{o}\mathrm{e}$$( \frac{\mathrm{i}\pi}{4})$,$\sin(\frac{\mathrm{i}\pi}{4})],\mathrm{i}=1,2$,$\cdots,8$.

(6)

Similarly, (i) its odd rank tensors are zero, and (ii) its initial three even rank tensors satisfy$\sum_{i=1}^{12}v_{k\iota\alpha}v_{k_{l}\beta}=4v_{k}^{2}\delta_{\alpha\beta}$,$\sum_{i=1}^{12}v_{ki\alpha}v_{ki\beta}v_{ki\gamma}v_{k\delta}|.=v_{k}^{4}\Delta_{\alpha\beta\gamma\delta}$$\sum_{\dot{l}=1}^{12}v_{k\dot{\iota}\alpha}v_{kt\beta}v_{ki\gamma}v_{ki\delta}v_{ki\mu k_{t}\nu}v$$= \frac{1}{6}v_{k}^{6}\Delta_{\alpha\beta\gamma\delta\mu\nu}$

,

(7)

.

DVM 2 is isotropic up to its 7th rank tensor.
We consider a binary mixture with two components, and , where the masses and temperatures of the two components are not significantly different. The interparticle collisions can bedivided into two kinds: collisions within the same species (self-collision) and collisions among
different species (cross-collision) [22]. Based on the DVM (1), the 2-dimensional BGK [23] kinetic
equation for species reads,
$A$

$B$

$A$

.

(8)

$\partial_{t}f_{ki}^{A}+\mathrm{v}_{ki}^{A}\cdot\frac{\partial}{\partial \mathrm{r}}f_{ki}^{A}-\mathrm{a}^{A}\cdot\frac{(\mathrm{v}_{ki}^{A}-\mathrm{u}^{A})}{\Theta^{A}}f_{k}^{A(0)}|.=J_{k\iota}^{AA}+J_{k}^{AB}|$

where$J_{k\dot{l}}^{AA}=-[f_{ki}^{A}-f_{ki}^{A(0)}]/\tau^{AA}$$f_{ki}^{A(0)}= \frac{n^{A}}{2\pi\Theta^{A}}\exp[-\frac{(\mathrm{v}_{kj}^{A}-\mathrm{u}^{A})^{2}}{2\Theta^{A}}]$,$\Theta^{A}=k_{B}T^{A}/m^{A}$, $J_{ki}^{AB}=-[f_{ki}^{A}-f_{ki}^{AB(0)}]/\tau^{AB}$

(9)

(10)

$f_{ki}^{AB(0)}= \frac{n^{A}}{2\pi\Theta^{AB}}\exp[-\frac{(\mathrm{v}_{k\dot{\mathrm{a}}}^{A}-\mathrm{u}^{AB})^{2}}{2\Theta^{AB}}]$, $\Theta^{AB}=k_{B}T^{AB}/m^{A}$

(11)

,
,
and $f^{AB(0)}$ are the corresponding Maxwellian distribution functions.
arethe local density, hydrodynamic velocity and temperature of species A.
,
are thehydrodynamic velocity and temperature of the mixture after equilibration process.
is theacceleration of species due to the effective external field.
For species , we have$f^{A(0)}$$n^{A}$$\mathrm{u}^{AB}$$\mathrm{u}^{A}$$T^{A}$$T^{AB}$$\mathrm{a}^{A}$

$A$

$A$

$f_{ki}^{A}$$n^{A}= \sum_{k_{\dot{1}}}$,$n^{A} \mathrm{u}^{A}=\sum_{ki}\mathrm{v}_{ki}^{A}f_{ki}^{A}$,$P^{A}(e_{\mathrm{i}\mathrm{n}\mathrm{t}}^{A}=n^{A}k_{B}T^{A})= \sum_{ki}\frac{1}{2}m^{A}(\mathrm{v}_{k}^{A}|.

whereis the local pressure (internal energy). For species
For the mixture, we have$P^{A}(e_{1\mathrm{n}\mathrm{t}}^{A})$,$\mathrm{u}^{AB}=(\rho^{A}\mathrm{u}^{A}+\rho^{B}\mathrm{u}^{B})/\rho$$B$-\mathrm{u}^{A})^{2}f_{ki}^{A}$(12), we have similar relations.

$nk_{B}T^{AB}= \sum_{kt}\frac{1}{2}[(\mathrm{v}_{k\mathrm{s}}^{A}-\mathrm{u}^{AB})^{2}m^{A}f_{ki}^{A}+(\mathrm{v}_{kv}^{B}-\mathrm{u}^{BA})^{2}m^{B}f_{ki}^{B}]$(13)
where $\rho^{A}=n^{A}m^{A}$ , $n=n^{A}+n^{B}$ and $\rho=\rho^{A}+\rho^{B}$ . Three sets of hydrodynamic quantities(for the two components , and for the mixture) are involved, but only two sets of them are
$A$

$B$

195

independent. So this is a tw0-fluid model. Without lossing generality, we focus on hydrodynamics of the two individual species. By expanding the local equilibrium distribution function$f^{AB(0)}$ aroundto the first order in flow velocity and temperature, the BGK model (8-11)$f^{A(0)}$becomes(14)$\partial_{t}f_{ki}^{A}+\mathrm{v}_{ki}^{A}\cdot\frac{\partial}{\partial \mathrm{r}}f_{ki}^{A}-\mathrm{a}^{A}\cdot\frac{(\mathrm{v}_{k\mathrm{n}}^{A}-\mathrm{u}^{A})}{\Theta^{A}}f_{kt}^{A(0)}=Q_{ki}^{AA}+Q_{ki}^{AB}$(15)$Q_{k\iota}^{AA}=-( \frac{1}{\tau^{AA}}+\frac{1}{\tau^{AB}})[f_{ki}^{A}-f_{ki}^{A(0)}]$$Q_{k}^{AB}|.=- \frac{f_{ki}^{A(0)}}{\rho^{A}\Theta^{A}}\{\mu_{D}^{A}(\mathrm{v}_{ki}^{A}-\mathrm{u}^{A})\cdot(\mathrm{u}^{A}-\mathrm{u}^{B})$$+ \mu_{T}^{A}[\frac{(\mathrm{v}_{kt}^{A}-\mathrm{u}^{A})^{2}}{2\Theta^{A}}-1](T^{A}-T^{B})-M^{A}[\frac{(\mathrm{v}_{k_{1}}^{A}\cdot-\mathrm{u}^{A})^{2}}{2\Theta^{A}}-1](\mathrm{u}^{A}-\mathrm{u}^{B})^{2}\}(16)$, $\mu_{T}^{A}=k_{B}n^{A}n^{B}/(\tau^{AB}n)$ , $M^{A}=n^{A}\rho^{A}\rho^{B}/(2\tau^{AB}n\rho)$ .
possesses. The continuous Maxwellian
Now, we go to the second step: formulate
an infinite sequence of moment properties. The Chapman-Enskog analysis[24] shows that,to follow the initial eight ones is sufficient to describe the samerequiring the discrete
Navier-Stokes equations,where$\mu_{D}^{A}=\rho^{A}\rho^{B}/(\tau^{AB}\rho)$$f^{A(0)}$$f_{ki}^{A(0)}$$f_{ki}^{A(0)}$

,

(17)

$\frac{\partial\rho^{A}}{\partial t}+\frac{\partial}{\partial r_{\alpha}}(\rho^{A}u_{\alpha}^{A})=0$$\frac{\partial}{\partial t}(\rho^{A}u_{\alpha}^{A})+\frac{\partial}{\partial \mathrm{r}_{\beta}}(\rho^{A}u_{\alpha}^{A}u_{\beta}^{A})+\frac{\partial P^{A}}{\partial r_{\alpha}}-\rho^{A}a_{\alpha}^{A}-\frac{\partial}{\partial r_{\beta}}[\eta^{A}(\frac{\partial u^{A}}{\partial r_{\beta}}+\frac{\partial u_{\beta}^{A}}{\partial r_{a}}-\frac{\partial u}{\partial r}\mathrm{L}\delta_{\alpha\beta})\gamma A]$$+ \frac{\rho^{A}\rho^{B}}{\tau^{AB}\rho}(u_{\alpha}^{A}-u_{\alpha}^{B})=0$

(18)

,

$- \partial \mathrm{e}_{\frac{A}{t}+\frac{\partial}{\partial r_{\alpha}}}\not\supset[(e^{A}+P^{A})u_{\alpha}^{A}]-\rho^{A}\mathrm{a}^{A}\cdot \mathrm{u}^{A}-\frac{\partial}{Tr_{\alpha}^{-}}[k^{A}\frac{\partial(k_{B}T^{A})}{\partial r_{\alpha}}+\eta^{A}u_{\beta}^{A}(^{\partial}*_{\beta}^{u^{A}}+\frac{\partial \mathrm{u}_{\beta}^{A}}{\partial r_{\alpha}}-\frac{\partial u_{\gamma}^{A}}{\partial r_{\gamma}}\delta_{\alpha\beta})]$$+_{\tau}^{\epsilon^{A}}*_{\rho}^{B}[(u^{A})^{2}-\mathrm{u}^{A}\cdot \mathrm{u}^{B}]+_{\overline{\tau}n}n_{\mathrm{I}^{n}arrow k_{B}(T^{A}-T^{B})-n^{A}\frac{\rho^{A}\rho^{B}}{2\tau^{AB}n\rho}(\mathrm{u}^{A}-\mathrm{u}^{B})^{2}=0}^{AB}$

,

(19)

where
, $\eta^{A}=P^{A}\tau^{AA}\tau^{AB}/(\tau^{AA}+\tau^{AB})$ , $k^{A}=2n^{A}\Theta^{A}\tau^{AA}\tau^{AB}/(\tau^{AA}+\tau^{AB})$ .
(20)$P^{A}=n^{A}k_{B}T^{A}$ , and neglecting
(17),
Eq.
using
quantity.
By
is a small
Recall that,
, Eq. (18) shows that the diffusion velocity,
the second and higher order terms in.
andis related to the gradients of, andare referred to Eq. (12) with . replaced by
The first three requirements on
the remaining five are$e^{A}=e_{\mathrm{i}\mathrm{n}\mathrm{t}}^{A}+ \frac{1}{2}\rho^{A}(u^{A})^{2}$$\mathrm{u}^{A}(\mathrm{u}^{B})$$u_{\alpha}^{B}-u_{\alpha}^{A}$$\mathrm{u}^{A}$$n^{A}$$T^{A}$$f_{ki}^{A(0)}$$f_{k}^{A}|$$\sum_{ki}m^{A}v_{ki\alpha}^{A}v_{kv\beta}^{A}f_{ki}^{A(0)}=P^{A}\delta_{\alpha\beta}+\rho^{A}u_{\alpha}^{A}u_{\beta}^{A}$$\sum_{k_{1}}m^{A}v_{ki\alpha}^{A}v_{ki\beta}^{A}v_{ki\gamma}^{A}f_{ki}^{A(0)}=P^{A}(u_{\gamma}^{A}\delta_{\alpha\beta}+u_{\alpha}^{A}\delta_{\beta\gamma}+u_{\beta}^{A}\delta_{\gamma\alpha})+\rho^{A}u_{\alpha}^{A}u_{\beta}^{A}u_{\gamma}^{A}$$f_{k\dot{l}}^{A(0)}$

(21)

(22)

196

(23)

$\sum_{ki}\frac{1}{2}m^{A}(v_{ki}^{A})^{2}v_{ki\alpha}^{A}f_{ki}^{A(0)}=2n^{A}k_{B}T^{A}u_{\alpha}^{A}+\frac{1}{2}\rho^{A}(u^{A})^{2}u_{\alpha}^{A}$$\sum_{kj}\frac{1}{2}m^{A}(v_{ki}^{A})^{2}v_{ki\alpha}^{A}v_{ki\beta}^{A}f_{ki}^{A(0)}=2P^{A}\Theta^{A}\delta_{\alpha\beta}+\frac{1}{2}P^{A}(u^{A})^{2}\delta_{\alpha\beta}$(24)$+3P^{A}u_{\alpha}^{A}u_{\beta}^{A}+ \frac{1}{2}\rho^{A}(u^{A})^{2}u_{\alpha}^{A}u_{\beta}^{A}$(25)$\sum_{k_{1}}\frac{1}{2}m^{A}(v_{ki}^{A})^{4}v_{ki\alpha}^{A}f_{ki}^{A(0)}=[12P^{A}\Theta^{A}+6P^{A}(u^{A})^{2}+\frac{1}{2}\rho^{A}(u^{A})^{4}]u_{\alpha}^{A}$

The requirement equation (25) contains the fifth order of the flow velocity
in polynomial up to the fifth order of :sufficient to expand$\mathrm{u}^{A}$. So it is$\mathrm{u}^{A}$$f_{k_{\dot{l}}}^{A(0)}$$f_{k\dot{\iota}}^{A(0)}=n^{A}F_{k}^{A} \{[1-\frac{(u^{A})^{2}}{2\Theta^{A}}+\frac{(u^{A})^{4}}{8(\Theta^{A})^{2}}]+\frac{v_{k_{2}\xi}^{A}u_{\xi}^{A}}{\Theta^{A}}[1-\frac{(u^{A})^{2}}{2\Theta^{A}}+\frac{(u^{A})^{4}}{8(\Theta^{A})^{2}}]$$+ \frac{v_{ki\xi}^{A}v_{ki\pi}^{A}u_{\xi}^{A}u_{\pi}^{A}}{2(\Theta^{A})^{2}}[1-\frac{(u^{A})^{2}}{2\Theta^{A}}]+\frac{v_{k\not\in}^{A}v_{\mathrm{k}t\pi}^{A}v_{k\iota\eta}^{A}u_{\xi}^{A}u_{\pi}^{A}u_{\eta}^{A}}{6(9^{A})^{3}}[1-\frac{(u^{A})^{2}}{2\Theta^{A}}]$$+ \frac{v_{ki\xi}^{A}v_{k\iota\pi}^{A}v_{ki\eta}^{A}v_{ki\lambda}^{A}u_{\zeta}^{A}u_{\pi}^{A}u_{\eta}^{A}u_{\lambda}^{A}}{24(\Theta^{A})^{4}}+\frac{v_{ki\epsilon^{v_{ki\pi}v_{ki\eta}v_{ki\lambda}v_{ki\delta}u_{\xi}u_{\pi}u_{\eta}u_{\lambda}u_{\delta}}}^{AAAAAAAAAA}}{120(\Theta^{A})^{5}}\}$(26)$+\cdots$

where$F_{k}^{A}= \frac{1}{2\pi 8^{A}}\exp[-\frac{(v_{k}^{A})^{2}}{2\Theta^{A}}]$

.

(27)

The truncated equilibrium distribution function $f_{ki}^{A(0)}(26)$ contains the fifth rank tensor of theand the requirement (22) contains its third rank tensor. Thus, a DVM being
particle velocity
isotropic up to its 8th rank tensors is enough to recover the physical isotropy of the continuous
Boltzmann equations to the Navier-Stokes level. So DVM (1) is an appropriate choice. To, one first needs calculate the factor.
is determined by thecalculate the discreteand the isotropic properties of the DVM (1). We finally obtaineight requirements on$\mathrm{v}^{A}$$f_{ki}^{A(0)}$$F_{k}^{A}$$F_{k}^{A}$$f_{ki}^{A(0)}$$\sum_{k\dot{\iota}}F_{k}^{A}=1$,$\sum_{k}F_{k}^{A}(v_{k}^{A})^{2}=\frac{\mathrm{e}^{A}}{6}$$\sum_{k}F_{k}^{A}(v_{k}^{A})^{6}=4(\Theta^{A})^{3}$,,$\sum_{k}F_{k}^{A}(v_{k}^{A})^{4}=\frac{2}{3}(\mathrm{e}^{A})^{2}$$\sum_{k}F_{k}^{A}(v_{k}^{A})^{8}=32(\Theta^{A})^{4}$,$\sum_{k}F_{k}^{A}(v_{k}^{A})^{10}=320(\Theta^{A})^{5}$(28), and other five nonzero ones, $v_{k}^{A}(k=1,2,3,4,5)$ are chosen,
Once a zero speed,$(k=0,1,2,3,4,5)$ will be fixed.
We come to the third step: finite-difference implementation of the discrete kinetic method.
There are more than one choices $[2, 18]$ available. One possibility is shown below,$v_{0}^{A}=0$$F_{k}^{A}$$f_{ki}^{A,(n+1)}=f_{ki}^{A,(n)}+[ \mathrm{a}^{A}\cdot\frac{(\mathrm{v}_{ki}^{A}-\mathrm{u}^{A})}{9^{A}}f_{ki}^{A(0)}+Q_{k\dot{\iota}}^{AA,(n)}+Q_{ki}^{AB,(n)}-\mathrm{v}_{ki}^{A}\cdot\frac{\partial f_{k\dot{l}}^{A,(n)}}{\partial \mathrm{r}}]\Delta t$

where the second superscripts , $n+1$ indicate the consecutive two iteration steps,
step; the spatial derivatives are calculated$n$$\Delta t$

,

(29)

the time$\mathrm{a}\epsilon$$(3f_{b,I}^{A,(n)}-4f_{ki,I-1}^{A,(n)}+f_{k\iota,t-2}^{A,(n)})/(2\Delta\alpha)$if$v_{k\iota\alpha}^{A}\geq 0$$(3f_{k_{\dot{l}},I}^{A,(n)}-4f_{ki,I+1}^{A,(n)}+f_{k\dot{\mathrm{r}},t+2}^{A,(n)})/(-2\Delta\alpha)$if$v_{k\mathrm{i}\alpha}^{A}<0$$\frac{\partial f_{ki}^{A,(n)}}{\partial\alpha}=\{$

(30)

137

where $=x$ , , the third subscripts $I-2$ , $I-1$ , , $I+1$ , $I+2$ indicate consecutive mesh nodes
in the direction.
If the kinetic numerical scheme is required to recover the hydrodynamics only up to theisothermal Navier-Stokes level, Eqs. (17)-(18) or the Euler level, Eqs.(17)-(19) with $\eta^{A}=k^{A}=$, following the same procedures, it is easy to find that DVM 2 is enough. For the isothermal
Navier-Stokes equation, Eq. (28) is replaced by$\alpha$$I$$y$$\alpha$$0$$\sum_{ki}F_{k}^{A}=1$,$\sum_{k}F_{k}^{A}(v_{k}^{A})^{2}=\frac{\mathrm{e}^{A}}{4}$,$\sum_{k}F_{k}^{A}(v_{k}^{A})^{4}=(\mathrm{e}^{A})^{2}$(31)$\sum_{k}F_{k}^{A}(v_{k}^{A})^{6}=6(\mathrm{e}^{A})^{3}$

For the complete Euler equation, Eq. (28) is replaced by$\sum_{ki}F_{k}^{A}=1$,$\sum_{k}F_{k}^{A}(v_{k}^{A})^{2}=\frac{\mathrm{e}^{A}}{4}$$\sum_{k}F_{k}^{A}(v_{k}^{A})^{6}=6(\Theta^{A})^{3}$.

,$\sum_{k}F_{k}^{A}(v_{k}^{A})^{4}=(\mathrm{e}^{A})^{2}$(32)$\sum_{k}F_{k}^{A}(v_{k}^{A})^{8}=48(\Theta^{A})^{4}$

In summary, to recover the tw0-dimensional complete Navier-Stokes equations, a 2model is sufficient to recovermodel is needed; adimensional 61 velocity
the tw0-dimensional Euler equations; recovering the tw0-dimensional isothermal Navier-Stokes
model. In principle, a DVM with lower isotropy can beequations can resort on a simpler
replaced by one with higher isotropy. But in practical simulations, one generally needs choose$\mathrm{D}2\mathrm{V}33$$(\mathrm{D}2\mathrm{V}61)$$\mathrm{D}2\mathrm{V}25$

the simplest one.
The validity of the formulated the FDLBMs is verified through two test examples. (The
Boltzmann constant $k_{B}=1.$ ) The first one is the isothermal and incompressible Couette flow
with a single component. In this case, $A=B$ . The initial state of the fluid is static. The
distance between the two walls is $D$ . At time $t=0$ they start to move at velocities , $-U$ ,work for such a
respectively. It is clear that all the three models
system. The horizontal velocity profiles of species or along a vertical line agree with thefollowing analytical solution,$U$$(\mathrm{D}2\mathrm{V}25, \mathrm{D}2\mathrm{V}33, \mathrm{D}2\mathrm{V}61)$

$A$

u

$B$

$= \gamma y-\sum_{j}(-1)^{j+1}\frac{\gamma D}{j\pi}\exp(-\frac{4j^{2}\pi^{2}\eta}{\rho D^{2}}t)\sin(\frac{2j\pi}{D}y)$

,

(33)

where $\gamma=2U/D$ is the imposed the shear rate, is an integer, the two walk locate at $y=\pm D/2$ .
(For example, see Fig. 1.)
The second one is the uniform relaxation process, which is an ideal process to indicate theequilibration behavior of the mixture. By neglecting the force terms and terms in spatial
derivatives, the Navier-Stokes equations (17)-(19) give$j$$\frac{\partial}{\partial t}\rho^{A}=0$

(34)

,

$\frac{\partial}{\partial t}(\mathrm{u}^{B}-\mathrm{u}^{A})=-\frac{1}{\rho}(\frac{\rho^{A}}{\tau^{BA}}+\frac{\rho^{B}}{\tau^{AB}})(\mathrm{u}^{B}-\mathrm{u}^{A})$,$\frac{\partial(T^{B}-T^{A})}{\partial t}=-\frac{1}{n}(\frac{n^{A}}{\tau^{BA}}+\frac{n^{B}}{\tau^{AB}})(T^{B}-T^{A})+\frac{\rho^{A}\rho^{B}}{2k_{B}n\rho}(\frac{1}{\tau^{AB}}-\frac{1}{\tau^{BA}})(\mathrm{u}^{B}-\mathrm{u}^{A})^{2}$

(35)

(36)

(The Euler equations play the same role as the Navier-Stokes equations in this case. Both theworks. ) The flow velocities of the two components equilibrate exponentiallyand) The equilibration of flow velocities also affects thatwith time. (For example, see Fig.
velocitydifference is zero, the temperatures equilibratetheflowtemperatures.
Whenof the) The simulation results agree well withexponentially with time. (For example, see Fig.
Eqs. (35) and (36)$\mathrm{D}2\mathrm{V}61$$\mathrm{D}2\mathrm{V}33$$2(\mathrm{a}).$$2(\mathrm{b}).$198$\supset\geq$$\mathrm{y}/\mathrm{D}$FIG. 1: Horizontal velocity profiles along a vertical line for the two species, and , at time $t=8$ .
The symbols are for simulation results. The solid line corresponds to the theoretical result, Eq. (33).
Parameters used in the tw0-fluid FDLBM are $m^{A}=m^{B}=1$ , $T=1$ , $n^{A}=n^{B}=1$ , $\gamma=0.001$ ,$\tau^{AA}=\tau^{BB}=\tau^{AB}=\tau^{BA}=0.2$. Parameters used
Eq. (33) are $\eta=\eta^{A}=0.1$ , $\rho=\rho^{A}=1$ .
$A$

$B$

$1\mathrm{n}|$$\mathrm{t}$$\mathrm{t}$FIG. 2: Uniform relaxation processes, (a) Equilibration of velocities; (b) Equilibration of temperatures. The symbols are for simulation results. The solid lines possess the theoretical slopes. Common
parameters for the simulations in (a) and (b) are $n^{A}=10$ , $n^{B}=1$ , $m^{A}=1$ , $m^{B}=10,7^{-AA}=\tau^{BB}=1$ ,$\tau^{AB}=10$ , $\tau^{BA}=1$ .
the initial conditions are $u_{x}^{A(0)}=-u_{x}^{B(0)}=-0.3$, $u_{y}^{A(0)}=u_{y}^{B(0)}=0$ , and$T^{A(0)}=1.3$ , $T^{B(0)}=0.7$ . The slope of the solid line in (b) is -11/20, which is consistent with Eq.
(35). In (b) the initial conditions are, and $T^{A(0)}=1.3$ , $T^{B(0)}=0.7$ . The slope ofthe solid line in (b) is -10.1/11, which is consistent with the first term of right-hand side of Eq.(36).
The second superscript “(0)” denotes the corresponding initial value. This figure shows an example
where the particle masses of the two species are significantly different.
$\ln(\mathrm{a})$$\mathrm{u}^{A(0)}=\mathrm{u}^{B(0)}=0$

III.

CONCLUSIONS AND REMARKSThe Chapman-Enskog analysis shows what properties the discrete Maxwellian distributionfunctionshould follow. Those requirements tell the lowest order of the flow velocity. The highest rank of tensors of the particle velocity
in the Taylor expansion ofin the requirements on the truncateddetermines the needed isotropy of the DVM. The
incorporation of the force tems makes no additional requirement on the isotropy of the DVM.
The present approach works for binary neutral fluid mixtures. One possibility to introduc$f_{ki}^{A(0)}$$\mathrm{u}^{A}$$\mathrm{v}^{A}$$f_{k\dot{v}}^{A(0)}$$f_{k\dot{\iota}}^{A(0)}$199interfacial tension is to modify the pressure tensors[14], which is implemented by changing
the force terms[3]. The specific force terms or pressure tensors depend on the system underconsideration, which are out of the scope of this Letter, but can be resolved under the same tw0. For binary fluids with disparate-mass components, saydimensional 61-velocity model$m^{A}\ll m^{B}$ , only if the total masses and temperatures of the two species are not significantly
different, Sirovich’s kinetic theory works, so do the corresponding FDLBMs. (See Fig. 2 forthe temperatures of the two components are greatly
an example.) When the masses
different, the tw0-fluid kinetic theory should be modified. In those cases, the Navier-Stokes
equations and the FDLBMs are not symmetric about the two components, but the FDLBMs
model. The formulation procedure is straightforward.
can still be resolved under the
A more detailed description is referred to $[25, 26]$ . We finally emphasize that the numericalerrors ffom the finite-difference schemes result in artificial viscosities in the simulation. The
comparison of various finite-difference schemes, discussion on numerical accuracy and stability
are referred to Refe. [2, 3, 18].
$(\mathrm{D}2\mathrm{V}61)$$\mathrm{a}\mathrm{n}\mathrm{d}/\mathrm{o}\mathrm{r}$$\mathrm{D}2\mathrm{V}61$

Acknowledgments

Aiguo Xu acknowledges Prof. G. Gonnella for guiding him into the LBM field and thanks
Profe. H. Hayakawa, V. Sofonea, S. Succi for helpful discussions. This work is partially supported by Grant-in-Aids for Scientific Research (Grand No. 15540393) and for the 21-th Century
COE “Center for Diversity and Universality in Physics” ffom the Minstry of Education, Culture
and Sports, Science and Technology (MEXT) of Japan.

[1] Succi S., The Lattice Boltzmann Equation, Oxford University Press, New York, 2001; Benzi R.,
Succi S. Vergassola M., Phys. Rep. 222, 3 (1992); Higuera F., Succi. S. Benzi R., Europhys. Lett.
9, (1989) 345.
, (1997)R2l; Seta T. Takahashi R. , J. Stat.
[2] Cao N., Chen S.,Jin S. Martinez D., Phys. Rev.
Phys. 107, (2002) 557; Sofonea V. Sekerka R. F., J. Comp. Phys. 184, (2003) 422; Watari M.
, (2004), (2003) 36306; Kataoka T. Tsutahara M., Phys. Rev.
Tsutahara M., Phys. Rev.
, (2004) R35701.
56702; Kataoka T. Tsutahara M., Phys. Rev.
, (2004)046702.
[3] Sofonea V., Lamura A., Gonnella G., Cristea A., Phys. Rev.
[4] Nannelli F., and Succi S., J. Stat. Phys. 68 , (1992) 401; Xi H., Peng G. Chou S. H., Phys. Rev.
, (2003) 16701; See also references
60, (1999) 3380; Ubertini S., Bella G., Succi S., Phys. Rev.
in [1].
, (2004) 065701(R).
[5] Li Y., LeBoeufe E. Basu P.K., Phys. Rev.
, (2002) R35301; ibid, 67, (2003) 36302[6] Luo L. S. Girimaji S. S., Phys. Rev.
[7] Malevanets A. and Yeomans J., Faraday Discuss. 112, (1999) 237.
, (2003) R35302.
[8] Guo Z. Zhao T. S., Phys. Rev.
, (2003) 56105; Physica A 331, (2004) 10;[9] Xu Aiguo, Gonnella G. Lamura A., Phys. Rev.
Physcia A 344, (2004) 750; c0nd-mat/0404205; Xu Aiguo, Commun. Theor. Phys. 39, (2003) 729.
[10] Gunstensen A. K. , Rothman D. H., Zaleski S., and Zanetti G., Phys. Rev. A 43, (1991) 4320.
, (1993) 4247.
[11] Flekkoy E. G., Phys. Rev.
[12] Grunau D., Chen S., Eggert K., Phys. Fluids A 5, (1993) 2557.
, (1996)3614.
[13] Shan X. Doolen G., J. Stat. Phys. 81, (1995)379; Phys. Rev.
[14] Orlandini E., Osborn W. R., Yeomans J. M., Europhys. Ltt. 32, (1995)463; Osborn W. R.,
Orlandini E., Swift M. R., Yeomans J. M., Banavar J. , Phys. Rev. Lett., 75, (1995) 4031; Swift
M. R., Orlandini E., Osborn W. R., Yeomans J. M., Phys. Rev. , 54, (1996) 5041; Gonnella, (1998) 480;
G., Orlandini E., Yeomans J. M., Phys. Rev. Lett. 78, (1997)1695, Phys. Rev.
45,
(1999)
314.
Lett.
Europhys.
M.,
G.,
A.,
YeomansJ.
Gonnella
Lamura[15] Kendon V. M., Desplat J. C., Bladon P., Cates M. E., Phys. Rev. Lett., 83, (1999) 576; Kendon
V. M., Cates M. E., Pagonabarrage I., Desplat J. , Bladon P., J. Fluid Mech. 440, (2001) 147.
[16] Yu H., Luo L. S., Girimaji S. S., Int. J. Comp. Eng. Sci. 3, (2002) 73.
$\mathrm{E}55$$\mathrm{E}69$$\mathrm{E}67$$\mathrm{E}69$$\mathrm{E}70$$\mathrm{E}$$\mathrm{E}68$$\mathrm{E}69$$\mathrm{E}66$$[]$$\mathrm{E}68$$\mathrm{E}67$$\mathrm{E}47$$\mathrm{E}54$$\mathrm{R}$$\mathrm{E}$$\mathrm{E}58$$\mathrm{C}$200[17] Facin P. C., Philippi P. , dos Santos L. O. E. in ICCS 2003, LNCS 2657 edited by Sloot P. M.
A. et al., 2003, pp.1007-1014.
[18] Cristea A. Sofonea V., Central European J. Phys. 2,(2004)382.
, (1993)1815; ibid, 49, (1994)2941.
[19] Shm X. Chen H., Phys. Rev.
, (2001)21503; Love P.
[20] Love P. J. Coveney P.V., Phys. Rev.
J. B., Coveney P.
V., ibid, 64, (2001)61302; Chin J. Coveney P. V., ibid, 66, (2002) 16303; Gonzalez-Segredo ,Nekovee M., Coveney P.V., ibid, 67, (2003)46304; Gonzalez-Segredo N. Coveney P. V., ibid, 69,(2004)61501.
[21] Sofonea V. Sekerka R. F., Physica A 299, (2001)494.
[22] Sirovich L., Phys. Fluids 5, (1962) 908; ibid, 9, (1966)2323.
[23] Bhatnagar P. L., Gross E. P., Krook M., Phys. Rev., 94, (1954)511.
[24] Chapmann S. and Cowling T. G., The Mathematical Theory of Non-Unifom Gases, 3rd ed.,Cambridge University Press, Cambridge 1970.
. 0rg/abs/c0nd-mat/04060lB].
[25] Xu Aiguo, Phys. Rev. (in press) [[26] Xu Aiguo, Europhys. Lett, (in press) [ $http//amiv$. 0rg/abs/c0nd-mat/04l04Sl].
$\mathrm{C}$$\mathrm{E}47$$\mathrm{E}64$$\mathrm{J}.,\mathrm{M}\mathrm{a}\mathrm{i}\mathrm{l}\mathrm{l}\mathrm{e}\mathrm{t}$$\mathrm{N}$$\mathrm{E}$$h\mathrm{f}\mathrm{f}\mathrm{l}p.\cdot//amiv$