Green-Kubo relations for lattice gas cellular automata

Volume 138, number 8

PHYSICS LETTERS A

10 July 1989

GREEN-KUBO RELATIONS FOR LATTICE GAS CELLULAR AUTOMATA M.H. ERNST and J.W. DUFTY’ Institute for Theoretical Physics, Princetonplein 5, P.O. Box 80006, 3508 TA Utrecht, The Netherlands Received 17 February 1989; revised manuscript received 21 April 1989; accepted for publication 2 May 1989 Communicated by A.R. Bishop

Green—Kubo relations are derived for linear transport coefficients (viscosities, diffusion) in lattice gases using the correlation function description of transport in fluids. The only ingredients are the local microscopic conservation laws without any further specification of the microdynamical laws.

The use of lattice gas cellular automata as microscopic models for macroscopic fluid dynamics has been illustrated in many recent applications [1]. The correspondencebetween these dynamical models and those for real fluids extends beyond the macroscopic level. Consequently, it is important to develop more explicitly the relationship between hydrodynamic equations and the underlying microscopic dynamics for cellular automata (CA). Our goal here is to present an elementary denvation of Green—Kubo time correlation function formulae for linear transport coefficients in fluids from a unifying point of view, followingmethods used for continuous fluids. The explicit form of the local microscopic conservation laws (see eq. (2)) is the only required ingredient, without any further specification of the detailed microdynamics. Here we give an outline and several applications; details are given elsewhere [2]. For specific models scattered results exist for transport coefficients in terms of Green—Kubo formulae. Hardy et al. [3] introduced the very first CAlattice gas, defined on a square lattice (referred to as HPP-model). These authors obtained Green—Kubo formulae for the viscosity coefficient in their model using the Landau—Lifshitz method of fluctuating hydrodynamics. This requires the explicit construction of Navier—Stokes type hydrodynamic equations and

the addition of a fluctuating Langevin force. The same method has been followed by Frisch et a!. [1] to obtain a Green—Kubo relation for the shear viscosity in the so-called FHP-models, defined on a hexagonal lattice. For the same FHP-models Rivet [1] derived a Green—Kubo formula by a different method, also standardly used in continuous fluids. Starting from the local equilibrium solution of the Liouville equation he applies the Chapman—Enskog method. Our derivation is more direct and quite general. The diffusion coefficient ofa tagged particle in CA-fluids or CA-Lorentz gases as a time sum of the velocity autocorrelation function is a well-known result, also contained in the present derivation. Pomeau [4] and Katz et al. [5] also derived Green—Kubo formulae for respectively the heat conductivity and nonlinear conductivity, that are not contained in our method. The reason is that there the local conservation laws do not possess the simple form (see eq. (2)) required here as a starting point. The microstate oflattice gas models is specified by the set of occupation numbers n(c, r, t), with values 0 or 1 (“Fermi exclusion rule”), for the velocity state c at lattice site r at the integer-valued time t. The local conserved densities Aa(r, t) = ~ aan(c, r; t)

(1)

C

Permanent address: Department of Physics, University of Florida, Gainesville, FL 32611, USA.

are determined by the collisional invariants aa (c)

0375-9601/89/$ 03.50 © Elsevier Science Publishers B.V. (North-Holland Physics Publishing Division)

391

Volume 138, number 8

PHYSICS LETTERS A

(e.g. mass aa(c) = 1, momentum aa(c)=c) and obey the local conservation laws: ~ a~(c)[n(c, r+c; t+l)—n(c, r; t)]=o. (2) No further specification of the microscopic dynamics is required, except that it allows an equilibrium state, uniform over all states {c, r}. Consequently, both stochastic and deterministic dynamics are included in the analysis. The equilibrium correlation functions of these densities, Gae(r’ 1; rt)

(r’,

=
t’ )M~(r,t)>

(3)

‘

with &4=A— determine the linear response of the average value of Aa to small deviations from equilibrium. It can be shown that they obey the unearized hydrodynamic equations for large space and time scales when such equations exist. We follow a method due to Zwanzig [6] with minor modifications for the discrete nature ofspace and time, and obtain for large time 1 and small wave number k the hydrodynamic equations for the FouriertransformsGai3(k,t)=of(3) in matrix form (the inner product is defined as =V~<öF*~G>): (a,+ikQ+k2A)G(k, t)=0, (4) where Q and A describe linearized Euler and Navier—Stokes hydrodynamics. The transport coefficientsofAap if they exist expressions, are given by the small-s limit the Green—Kubo —

—

t (x~)

Aa 6(S) =

(5)

e_s

where the asterisk indicates that the term t 0 only counts half. The projected current is defined as Qapap(c)16n(c, r;

of velocity states per site. The inner product is defined as (jig) = ~~f(c)g(c). We give an outline of the derivation. Fourier transformation allows us to write eq. (2) formally as a balance equation: Aa(k, 1+ 1) Aa(k, t) = jkBc~(k, t+ 1) where the current is (with c,=i~c) Ba

,

(8)

( k~1) = ~cjaa( c) [(e ~kC_ 1) /ikc] ~(c, k, 1,) C

By taking the inner product

(9) eq. (8)> we


obtain an equation for G~(k, 1) in terms of = because of stationarity. By applying eq. (8) once more we eliminate A~to obtain a current—current correlation function, yielding for the Laplace transform of <~~(t)>

[es_ 1 +L(k, s)] =

.

(10)

For small k the matrix L(k, s) is 2[ L(k,s)=ik]’, (11) where =>~ e~’. The s-dependent terms reduce in the small-k limit to the r.h.s. of eq. (5) minus the term at t = 0. Subsequent k-expansion of ik yields ikk 2kO; the first term equals ikQ. 0+ ~kNext we consider the long-time or small-s limit of eq. (10) where es_ 1 is replaced by s+ ~s2... =s— ~k2Q2 i.e. the second time derivatiye of slow variables is eliminated using the Euler equations [1]. Combination of both 0(k2) terms yields A”=~k_o’—~Q2, (12)

(6)

which represents the contributions of the t=0 term in eq. (5). Note that it is half the correlation func-

1= c and ~n = n . The Euler and susceptibility matrices are defined respectively by

tion of subtracted currents at 1 = 0. The quantity II” is in the CA-context usually referred to as the propagation part [1] of the matrix of transport coefficients, since it originates from the Euler part of the hydrodynamic equations. A more detailed derivation is given elsewhere [2]. The Green—Kubo relations (5) are the main results of this Letter. Next some comments and applications will be presented. The primary differences

Ja(t) =

~ [ciaa(c)

10 July 1989

t)

,

with c —

Qa# =K( Zap = K(

c/aa I a~)[(ala)—’] a I a) afi = K( aa I ap)

~p,

(7)

where the fluctuation K= < [~n(c, 0)2> = x (1 ) on account of the Fermi exclusion rule. The average density p = b < n> where b is the number —

392

—

Volume 138, number 8

PHYSICS LETTERS A

between eq. (5) and the equations for continuous fluids are due to the discrete sum over time instead of an integration, and the fact that only one half of the 1=0 term occurs. The existence of Navier—Stokes hydrodynamics requires that A(s) remains finite as s—0. For typical fluid models A(0) is finite for dimensionality d> 2. In d= 2 and d= 1 typical matrix elements of A(s) show small-s divergences proportional to ln s and s 1/3 respectively, leading to the well-known long-time tails [7]. Important applications are the viscosities of one component fluids. Here the collisional invariants are associated with number and momentum conservation, —

10 July 1989

instance i= or ~ K— 2i~/d,or they may be expressed in invanants. In an isotropic fluid there are only two viscosities (0=0). On certain lattices (such as the hexagonal lattice, used in the FHP-models [1]), fourth rank tensors have the same symmetry as in an isotropic fluid, and so do the corresponding hydrodynamic equations. In such systems the shear viscosity ~ and bulk viscosity C~ expressed in terms of invariants, read: ~(s)=p[(d— l)(d+2)X, V]’ X~

*

e_st

t=0

a0(c): a~=1

,

a1=~C~ê=c1,

2xtV) —‘

*

e~’
C(s) =p(d a,=~,c=c, (1=1,2,..., d—1) (13) where ~i, ~2 I } is a set of d orthonormal unit vectors. The susceptibility matrix Xap=ôapXa with x~ = bK and Xi = X~= bKc~,where the 2.sound yeThe Eulocity is defined through = (1 /bd) ~C c elements ler matrix has only twoco~non-vanishing Q,~=1 and Q~,,= c~.The only non-vanishing elements ofthe transport matrix A 11, A1, and A,,. are given in terms of the viscosity tensor

{1c,

e~’,

~1zjrs(5)

1(0)J,1(t)>

,

(14)

Xi

with

(17)

t=0

,

where

summation convention is used and ô,,J~/dis a traceless tensor. U— U~ Several observations areC in here. The formula for the bulk viscosity is aorder new result. It shows directly that the bulk viscosity vanishes identically in single speed models, because the trace J= 0 since J° —J

Ic 2 = dc~.A Green—Kubo formula for the shear viscosity in the special case of the FHP-model has been derived by Rivet. If all particles have the same speed (FHP-I), Rivet’sresultfor~agreeswith(17). However, for the models FHP-II and FHP-III, that inelude a rest particle, Rivet’s expression namely, eq. —

~ (c,c1 —ô,3c~).3n(c,r; 1)

(15)

(17) with J°replaced by J does not represent the shear viscosity, but equals i~~+~c.

through the relations pA,,= pA,1= ~7~ii and pA,,. where subscripts land t are defined as in

The general formula (17) for the viscosities also apply to a one-dimensional fluid type cellular automaton, with five velocities per site (c=0, ±1, ±2;

general form of the fourth rank viscosity tensor depends in general on thelattice symmetries of the derlying lattice. On a cubic (as used in unthe HPP-model [3]) this tensor has in general three independent viscosity coefficients

b= considered invariants by Qian et al. [8]. The model has 5), twoas2/5=2. collisional 1, a,=c) and As there are only(an= longitudinal cornc~=~c ponents, there is no shear viscosity and the bulk or longitudinal viscosity is given by (12) with Xi= ~(S—p) This viscosity is expected to diverge strongly in one dimension, K(s) —~s’13[7]. The resulting sound wave damping constant is expected to be divergent at large wave length, as is actually observed in computer simulations [81. Next we consider lattices with cubic symmetries. The very first CA-fluid considered, the two-dimensional HPP-rnodel, has this symmetry. Also Green—

J~(t) =

,

~,

~

—

‘.

~1,jrs~1[ôirôjs +ôjsöjr — (2/d)5,~,ö,~) + KÔUÔ,.S + 0ö~?,

(16)

where ö is an invariant tensor under cubic symmetry transformation with elements ö~ = 1 if ir~J=r=s and vanishes otherwise. There are many ways to express the viscosities in tensor elements; for (4)

393

Volume 138, number 8

PHYSICS LETTERS A

Kubo formulae have been given for this model [3]. If such models admit only zero velocities (rest partides) and velocities that are multiples of the 2dnearest neighbor lattice vectors, then Jr,, in (10) vanishes and ~ on account of (11). If one considers in addition a single speed model (where JcI2=dc~o),then J,1=J,°~ is traceless and eq. (11) implies K= Old. With these constraints (that all apply to the HPP-model) the viscosity tensor takes the simple form

10 July 1989

fusion equation with diffusion coefficient given by the small-s limit of D(s)=(bd)’ ~

x

*

e~’

~ c’ .c
1)>

r.c,~

—

= 0 [ ó,~— (lId) öijtärs]

,

(18)

where 0 is a positive viscosity coefficient. Combining this with eqs. (14) and (4) yields the hydrodynamic equations for the HPP-model, as quoted in ref. [3] apart from an inconsistent sign in front of their symbol 0. From the relation ~ and (14) one obtains a possible version of the Green— Kubo formula, 0(s) =p(2x 1 V) —‘

x 1=0 ~ * e ~‘< [~(0) —Jv%,(0)] [J~~(t) —J~~(t)]>. (19) This result resembles that of ref. [3] apart from a 2x missing factor ( 1) —‘ and the missing relative weight ~ for the 1=0 term. The relation 0= (d— i ) —‘ > in combination with (14) represents the viscosity as a cubic invariant. Other important applications concern diffusion: for instance, diffusion in binary mixtures of different particles or diffusion in colored mixtures of otherwise identical particles; self-diffusion of a tagged particle in a CA-fluid or in a CA-Lorentz gas of immobile scatterers; collective diffusion in a one component lattice gas (non-fluid type CA), where the microscopic dynamics only conserves numbers, but not momentum. These applications will be reported elsewhere. Here we treat only the examples of tagged particle diffusion and diffusion in a colored mixture. To study tagged particle diffusion, let n*(c, r; t) denote the occupation number for a tagged particle in a one component fluid. The labelling of this partide is such that it does not change its dynamics from that of the other particles. Then, for long times and large space scales the probability density to find the tagged particle at a particular location obeys a dif-

394

=d

~ 0 * e ~‘

(20)

I=

where the last line identifies the velocity autocorrelation function. This relationship of the self-diffusion coefficient to the velocity autocorrelation function is the expected one, in analogy to the continuous fluid result. It applies for both fluid models such as FHP and HPP and lattice Lorentz models. Tagged particle diffusion in lattice gas cellular automata can also be studied through the mean square 2>-+2dDt as t—~cc. In the case of<[r(t)—r(0)] diffusion in a color mixture we condisplacement, sider first a one component lattice gas and label a fraction, x, of the particles. As above the labelling is done such that there is no change in the dynamics of these particles. To distinguish the red and white partides we introduce a color variable, a(c, r; t), that takes the value 1 for red and 0 for white for every occupied state. The relevant conserved quantity is now the variable >.. n (c, r; 1) a(c, r; t). Its average is px, where p is the average number density with 0.~op
D~(s)= (dx,. V) —‘

* 1=0

e

“



J,.(t)= ~ cn(c, r; t) [a(c, r;

t)

—x]

(21)

.

r,,

The susceptibility here is given by x,.—px (1 x). This is also a new result! In concluding this Letter we emphasize again the unifying features and general validity of the Green— Kubo formula (5) for CA-fluid models in which the microdynamics obeys the following requirements: (i) all collisions have zero impact parameter (point partides), (ii) all collisions obey certain conservation —

Volume 138 number 8

PHYSICS LETTERSA

laws (e.g. number, momentum); (iii) there exists a spatially uniform equilibrium state. Restriction (i) yields the microdynamics: n (c, r+ c, I + 1) = n (c, r, 1) + 4 (c, {n}) where A represents the collision term. The collision rules may be deterministic or stochastic; the models may be single speed or multiple speed, property (ii) then guarantees the validity of the local conservation law eq. (2). Several CA-fluids such as lattice Lorentz gases and the FHP-models have extra global invariants leading to staggered momentum densities that also satisfy eq. (2). Therefore the diffusivities of the corresponding soft modes can also be expressed in Green—Kubo formulae In any sensible CA-fluid a uniform equilibrium state must exist. For the validity of the present derivation it is irrelevant whether the model violates or obeys detailed balance (with possibly different rate constants for the forward and backward collisions). The research of J.W.D. was supported by the Pieter Langerhuizen Lambertuszoon Fonds, and the

10 July 1989

hospitality of the Instituut voor Theoretische Fysica is gratefully acknowledged. Also, M.H.E. is grateful to the organizers for their invitation to participate in the CECAM workshop on Cellular automata and Brownian motion (Orsay, August 1988), where some preliminary results were obtained.

References [11U. Frisch, D. d’Humières, B. Hasslacher, P. Lallemand, Pomeau and J.P. Rivet, Complex Syst. 1 (1987) 649;

Y.

J.P. Rivet, Complex Syst. 1 (1987) 839. [2] Dufty0. anddeM.H. Ernst, Phys. Chem.,Phys. to be Rev. published. [3] J.W. J. Hardy, Pazzis andJ.Y. Pomeau, A 13 (1976) 1949. [4] Y. Pomeau, J. Phys. A (1984) L4 15. [5] S. Katz, J. Lebowitz and H. Spohn, J. Stat. Phys. 34 (1984)

[61R. Zwanzig, Phys. Re~40 (1964) 2427. [71 Y. Pomeau and P. Resibois, Phys. Rep. 19

(1975) 69. [8] Y.H. Qian, D. d’Humières and P. Lallemand, in: 17th ICTAM, August 1988); Y.H. Qian, private communication (1988).

395