# Finite-Difference Lattice Boltzmann Methods for Binary Fluids

###### Abstract

We investigate two-fluid BGK kinetic methods for binary fluids. The developed theory works for asymmetric as well as symmetric systems. For symmetric systems it recovers Sirovich’s theory and is summarized in models A and B. For asymmetric systems it contributes models C, D and E which are especially useful when the total masses and/or local temperatures of the two components are greatly different. The kinetic models are discretized based on an octagonal discrete velocity model. The discrete-velocity kinetic models and the continuous ones are required to describe the same hydrodynamic equations. The combination of a discrete-velocity kinetic model and an appropriate finite-difference scheme composes a finite-difference lattice Boltzmann method. The validity of the formulated methods is verified by investigating (i) uniform relaxation processes, (ii) isothermal Couette flow, and (iii) diffusion behavior.

PACS numbers: 47.11.+j, 51.10.+y, 05.20.Dd

## I Introduction

Gas kinetic theory plays a fundamental role in understanding many complex processes. To make solutions possible, many of the kinetic models for gases are based on the linearized Boltzmann equation, especially based on the BGK approximationPR94511 . Even thus, only in very limited cases are analytic solutions available. Basically speaking, there are two options to simulate Boltzmann equation systems. First, one can design procedures based on the fundamental properties of rarefied gas alone, like free flow, the mean free path, and collision frequency. Such a scheme does not need an a priori relationship with the Boltzmann equation, but the scheme itself will reflect many ideas and/or concepts used in the derivation of Boltzmann equation. In the best case, such a simulation will produce results being consistent with the solution of Boltzmann equation. The second option is to start from the Boltzmann equation and design numerical schemes as accurate as possibleCarlo . The discrete Boltzmann equation approach or lattice Boltzmann method (LBM) has been becoming a viable and promising scheme for simulating fluid flowsSucci .

LBMs for single-component fluids have been well studied, while for binary mixtures still need more clarificationSuccip1 . For binary fluids, although various LBMs have been proposed PRE6635301 ; IJCES373 ; Yeomans2fm ; PRE6835302 ; PRE6756105 ; PRA434320 ; PRE474247 ; PFA52557 ; JSP81379 ; EPL32463 ; PRL83576 ; ICCS2003 ; CEJP2382 ; VictorIJMPC ; VictorPRE ; ShanChen ; Coveney ; XuEPL , most of them PRE6835302 ; PRE6756105 ; PRA434320 ; PRE474247 ; PFA52557 ; JSP81379 ; EPL32463 ; PRL83576 ; ICCS2003 ; CEJP2382 ; VictorIJMPC ; VictorPRE are based on the single-fluid theoryPhysA299494 . For systems with different component properties, a two-fluid theory is necessary. Sirovich’s two-fluid kinetic theoryPF5908 works for (approximately) symmetric systems where the two components have (approximately) the same total masses and local temperatures. A LBM based on Sirovich’s theory and for the complete two-dimensional Navier-Stokes equations(NSE) is given in XuEPL . This LBM is based on a two-dimenaional model with sixty-one discrete velocities (D2V61). Many compressible fluids can be well described by the Euler equationsDalton . In fluid mechamics of low-speed flow, the temperature remains nearly constant and consequently the isothermal NSE description is extensively usedDalton . From the Chapman-Enskog procedureChap the Euler equation is a lower-order approximation compared with the NSE. The isothermal NSE is a simplified case of the complete NSE. For the above two kinds of systems, using the LBM for complete NSE system is not neccessary and computationally inefficient. In this study we generalize Sirovich’s theory so that it works also for asymmetric systems where the total masses and/or local temperatures of the two components are greatly different, then formulate LBMs for the two kinds of systems. The LBMs formulated here require simpler discrete velocity models(DVMs). For the Euler-equation system a DVM with thirty-three discrete velocities (D2V33) is enough. For the isothermal NSE system, a D2V25 is sufficient.

This paper is arranged in the following way: In section II we review and develop the two-fluid BGK kinetic theory. Sirovich’s original treatments are clarified and summarized in models A and B. For asymmetric systems three kinetic models (C, D and E) are derived. The hydrodynamics and diffusion behavior of the model systems are discussed. In section III the kinetic models are discretized based on a multispeed discrete velocity model. Then, possible FD schemes are given and the corresponding numerical viscosities and diffusivities are analyzed. Numerical tests are shown in Section IV. Section V concludes the present paper.

## Ii Two-fluid BGK kinetic theory

In a binary system with two components, and , roughly speaking, the approach to equilibrium can be divided into two processes. One is referred to as Maxwellization (i.e., each species equilibrates within itself so that the local distribution function approaches to its local Maxwellian). The other is the equilibration of species (i.e., the differences in hydrodynamic velocities and local temperatures of the two components eventually vanish). Correspondingly, the interparticle collisions fall into two categories: self-collisions (collisions within the same species) and cross-collisions (collisions between different species)PF5908 ; PRE6635301 .

### ii.1 General description

For a two-dimensional binary gas system the BGK kinetic equations readPRE6635301 ,

(1) |

(2) |

where

(3) |

(4) |

(5) |

(6) |

( ) and () are the distribution function and particle velocity of the component (); and are the local Maxwellians which work as references for the self- and cross-collisions; , , are the local number density, hydrodynamic velocity and temperature of the species ; , are the local hydrodynamic velocity and temperature of the mixture after equilibration process; is the acceleration of the species due to the effective external field. For species , we have

(7) |

(8) |

(9) |

(10) |

(11) |

where and () are the local mass density and internal mean kinetic energy (hydrostatic pressure) of species , is the Boltzmann constant. For species , we have similar relations.

For the mixture, we have

(12) |

(13) |

(14) |

(15) |

where , , , , , are the total number density, total mass density, barycentric velocity, mean temperature, total internal energy, and total hydrostatic pressure, respectively. It is easy to find the following relations,

(16) |

(17) |

Here three sets of hydrodynamic quantities [(, , ), (, , ) and (, , )] are involved. If assume that the two components are in local equilibrium, implying that , , can be replaced by and , , can be replaced by in the definitions of , and , we arrive at the one-fluid theory and Eq. (17) recovers Dalton’s lawDalton , where and are the temperature and the velocity of the system in the complete equilibrium. It is clear that the one-fluid theory is conditionally valid. If the differences among , , and/or among , , are not small, the above replacements result in large errors. Since each set of the hydrodynamic quantities can be described by the other two sets, in such cases, a two-fluid theory is preferable. Without loss of generality, we require the description to be dependent on (, , ) and (, , )foot1 .

A key point to complete the two-fluid kinetic description is how to calculate the local Maxwellian (). Within Sirovich’s original treatments, it is Taylor expanded around () to the first order of flow velocity and temperaturePF5908 . This treatment is reasonable when the hydrodynamic properties of the two components are nearly symmetric, i.e., , . To make a general theory working also for asymmetric systems where the hydrodynamic properties of the two components are greatly different, we introduce the reference distribution function in a general way and do the Taylor expansion around it.

For , we choose the reference distribution function as

(18) |

where the second superscript “r” means “reference” and the corresponding quantities are the reference hydrodynamic quantities which take values in the following way,

(19) |

(20) |

(21) |

Let us make the solutions more explicit. Firstly for , from Eq.(13) we have

(22) |

Then for , from Eq. (16), when we have

(23) | |||||

when we have

(24) | |||||

Considering together (20)-(24) gives

(25) |

(26) |

In the case of and , gets back to .

Both of and are local quantities. Their values are functions of position and time. It is possible for such a phenomenon, but , to occur, where and are two different positions in the system. While in a theory it is not convenient to use the reference state in such a way: and . Instead, we prefer to use one of the two possibilities, or , in the whole system, where is an arbitrary position in the system. For we have the same preference. To that aim, ,, and in the criteria (25) and (26) are replaced by their spacially averaged values, , , and , respectively. This treatment is reasonable from a statistical sense.

### ii.2 Kinetic models for symmetric systems

For systems with and , we can use , , i.e, Sirovich’s kinetic theory. In this case the equations for the two components are symmetric. The cross-collision term in (1) becomes

(27) | |||||

where

(28) |

If we concern the hydrodynamics only up to the NSE level, () in the force term can be replaced by (). The BGK model (1-6) can be rewritten as

(29) |

(30) |

where

(31) |

(32) | |||||

the expressions of and are obtained from Eqs. (31) and (32) via formal replacements of the superscripts and . In the isothermal case, , the expression of is simplified as

(33) |

For the convenience of description, the kinetic model with (29)-(32) is referred to as kinetic model A; the one with (29)-(31),(33) is referred to as kinetic model B.

### ii.3 Kinetic models for asymmetric systems

#### ii.3.1 Kinetic model C: for isothermal systems with

For such a system, , and

(34) |

(35) | |||||

Thus, within kinetic model C,

(36) |

(37) |

(38) |

(39) |

where

(40) |

#### ii.3.2 Kinetic model D: for systems with and

The references are and

(41) |

Since

(42) | |||||

within kinetic model D

(43) |

(44) |

(45) | |||||

(46) | |||||

where

(47) |

#### ii.3.3 Kinetic model E: for systems with and

In this case, the reference velocity and reference temperature for both and are and , respectively.

(48) |

(49) |

Within the kinetic model E

(50) |

(51) |

(52) | |||||

(53) | |||||

### ii.4 Hydrodynamics and diffusion

#### ii.4.1 Hydrodynamics

A connection between a kinetic model and corresponding hydrodynamics is the Chapman-Enskog analysisChap . All above kinetic models contribute to (i) the same continuity equation at the Euler and the NSE levels,

(54) | |||||

(55) |

(ii) the same Euler momentum equations,

(56) |

(57) |

and (iii) the same NSE momentum equation for component ,

(58) |

where

(59) |

describes the momentum transferred from component to , and it is also the diffusion flux density which will be clear from a later equation (72);

(60) |

is the stress tensor,

(61) |

is the viscous stress tensor, and

(62) |

is the viscosity.

Moldes A and B contributes to symmetric hydrodynamics for the two components. The Euler energy equation of model A for component reads,

(63) |

where

(64) |

is the local total energy, and

(65) |

is the heat transfered from component to .

The NSE momentum equation for component from model C reads

(66) |

where the definition of is similar to that of , and

(67) |

is an additional stress tensor due to the asymmetry of densities of the two components.

The Euler energy equation of model D for component is the same as Eq. (63) and for component reads

(68) |

where the definition of is similar to that of and

The Euler energy equations from kinetic model E are as follows,

(69) |

(70) |

#### ii.4.2 Diffusion

From the continuity equations (54)-(55) we have

(71) |

and

(72) |

where is given in Eq. (59) and it is the amount of the component transported relative to the component by diffusion through unit area in unit time. For the incompressible fluids where is a constant, the continuity equation (72) is equivalent to the following diffusion-convection equation,

(73) |

where . The diffusion velocity is determined by the momentum equation. We can find a simple relation for it in the following case: We consider a binary system without external forces and where the flow velocities , are small and their derivatives can be regarded as higher-order small quantities. From the momentum Eq. (56) or (58), by neglecting the second and higher-order terms in and/or , then using the definition (11), we obtain

(74) | |||||

(75) |

If further assume the system to be isothermal, the density flux of component reads

(76) |

where

(77) |

is the diffusivity of component . Eq. (76) is Fick’s first lawWeb . From Eqs. (73) and (75) we have