- Email: [email protected]
Kinetic equation for liquids with a multistep potential of interaction. H-theorem
PHYV,,A ELSEVIER
Physica A 234 (1996) 89-107
Kinetic equation for liquids with a multistep potential of interaction. H-theorem I.P. Omelyan*, M.V. Tokarchuk Institute for Condensed Matter Physics, National Ukrainian Academy of Sciences, 1 Svientsitsky St., UA-290011 Lviv, Ukraine Received 12 December 1995
Abstract
A new kinetic equation for dense gases and liquids with the interparticle potential of interaction in a form of the multistep function is proposed. The equation is consistently derived on the basis of Zubarev's nonequilibrium statistical operator method using a modification for weakening correlations principle which takes into account local conservation laws. In particular cases this equation coincides with the well known kinetic equations for hard-sphere and square-well potentials. Application of the obtained kinetic equation to real systems is discussed. A proof of a global H-theorem is given.
PACS: 05.20.Dd, 51.10.+y, 82.20.Mj Keywords: Kinetic equation; Dense gases; Multistep potential; Boundary condition; H-theorem
1. Introduction
The construction of kinetic equations for dense gases and liquids is one of the main problems in the kinetic theory of classical systems. The Boltzmann equation gives the foundations for the description of nonequilibrium properties of a dilute gas, but at the same time we do not have an analogous kinetic equation for dense systems with a realistic interparticle potential of interaction. In order to give more rigorous theoretical assumptions than Boltzmann's, Bogoliubov [1] proposed a new approach to the problem of obtaining the kinetic equations. It allows to generalize the resuits to higher densities. In the framework of this method the Boltzmann equation appears as a zeroth approximation in the expansion in powers of density. But if one tries to investigate higher orders of the expansion in density, it appears that terms beginning with the second-order contribution are divergent [2]. The reason for the * Corresponding author. E-mail: [email protected]. 0378-4371/96/$15.00 Copyright (~ 1996 Elsevier Science B.V. All rights reserved PII S 0 3 7 8 - 4 3 7 1 ( 9 6 ) 0 0 2 8 6 - 5
90
I.P. Omelyan, M.V. Tokarchuk/ Physica A 234 (1996) 89-107
divergences lies in the fact that the expansion was based on dynamics of isolated groups of particles in the infinite space, regardless of the medium, i.e., the rest of the particles. Theories of Kirkwood [3,4], Born-Green [5], Cohen [6,7], Garcia-Colin and Flores [8], Sandri [9] and Frieman [10] lead to results which are similar to those obtained by the Bogoliubov's method. In order to obviate the divergences, the region of integration must be cut off in coordinate space. As a result of this truncation, an expansion of the transport coefficients in density contains nonanalytic logarithmic terms [11]. In a more rigorous approach it is necessary to perform a partial resummation of the BBGKY hierarchy [12, 13]. This problem is practically unsolvable for high densities. On the other hand, the solution to kinetic equations for dense gases is a very complicated problem. These equations contain the many-particle scattering operators and this fact does not allow to solve them in general. That is why the construction of kinetic equations for dense gases and liquids with model interparticle potentials of interaction is a problem of great importance. The first step in this direction was made in the semiphenomenological standard Enskog theory (SET) [14] of dense gases. Simple model of hard spheres, used in this theory, was very successful in explaining the dependence of the transport coefficients on density. In order to apply the results of SET to realistic systems, Enskog proposed to substitute hydrostatic pressure of the hard spheres system by the thermodynamic pressure of a real system. Taking this suggestion as the basic one, Hanley et al. [15] elaborated the modified Enskog theory (MET). In this theory the diameter of hard spheres is determined by the second virial coefficient of the equation of state for the real system. Thus, it depends on temperature and density. Using different equations of state, BH (Barker, Henderson) [16], WCA (Weeks, Chandler and Andersen) [17], one obtains the corresponding versions of MET theory. In 1961, Davis, Rise and Sengers [18] proposed the kinetic theory DRS. Within the framework of this theory an interparticle potential of interaction is chosen in the square-well form. An attractive part of the real potential is approximated here by a finite height wall. The case where the potential between the hard sphere and attractive walls is considered as smooth tail (instead of the constant potential in DRS) has also been investigated [19]. In the kinetic mean-field theory (KMFT) [20] a smooth attractive tail is considered along with the hard sphere potential. However, as is now well established, all these potential models do not describe real systems satisfactorily. In particular, as was shown in recent calculations [21], KMFT does not give a correct description of the temperature dependence for the transport coefficients. From the point of view of the statistical theory of nonequilibrium processes, SET and DRS contain two essential drawbacks. The first one is that their kinetic equations are not obtained in the frames of some consequential theoretical scheme and one does not know how to improve these theories. The second one is that the H-theorem is not proved. Nevertheless, not long ago the way of constructing of SET entropy functional was given [22]. In order to overcome these drawbacks Karkheck and Stell [23] have
I.P. Ornelyan, M.V. TokarchuklPhysica A 234 (1996) 89-107
91
obtained the kinetic equation for a revised version of the Enskog theory (RET), using the principle of maximum of the entropy functional. This equation was proposed for the first time by Van Beijeren and Ernst [24] on the basis of the diagrammatic method. Resibois [25,26] proved the H-theroem for it. In 1985, the revised version of DRS theory (RDRS) has been proposed [27]. The kinetic equation of RDRS satisfies the H-theorem as well. In the paper of Stell et al. [28], the authors finished construction of the kinetic variational theory (KVT) introducing a whole family of theories, namely, KVT I, KVT II and KVT III. The KVT III refers to the local energy constraint in maximizing entropy as distinguished from the global energy constraint (KVT II) or hard-core constraint (KVT I). According to this classification, the kinetic equation of RDRS [27] is obtained applying KVT III to the square-well potential and KMFT [20] is derived from KVT I as an application to the potential with a smooth tail. In Ref. [28], the KVT III version of KMFT was proposed. In this case, the quasi-equilibrium binary distribution function (QBDF) of hard spheres is to be substituted for the full QBDF, which takes into account the full potential of interaction. The main conclusion of the KVT I version of KMFT [20] is that the smooth part of the potential does not contribute explicitly to the transport coefficients. This fact is caused by the approximations of the theory. When the full QBDF is used, as in the KVT III version of KMFT, the transport coefficients are determined by the soft part of the potential. Zubarev and Morozov [29] using the Zubarev's nonequilibrium statistical operator (NSO) method [30] have obtained kinetic equations for dense gases by formulation of a new boundary condition for the BBGKY hierarchy. This condition takes into account the correlations connected with the local conservation laws of hydrodynamic quantities. Such an approach leads to results similar to KVT III. In particular, it was shown [31] that within this approach the kinetic equation for hard spheres coincides with the kinetic equation of RET. In the present paper, in order to choose a more realistic model, a kinetic equation for dense classical systems is proposed for the potential of interaction in the form of a multistep function. This function consists of the hard-sphere part and a system of repulsive and attractive walls. The derivation of the kinetic equation is performed on the basis of theoretical scheme of the kinetic theory construction by the NSO method with the help of a modification to the weakening correlations principle [29]. It is shown that the obtained kinetic equation satisfies a global H-theorem. For certain parameters of the model multistep potential, this equation coincides with the kinetic equations of RET, RDRS and KVT III version of KMFT theories. Thus, the kinetic theory proposed in this paper can be considered as a general one, which covers the results of previous theories for particular cases and, consequently, shows a connection of these theories with one another. We also discuss the region of application of the generalized kinetic equation and ways of possible improvement of the theory.
I.P. Omelyan, M. V. Tokarchuk/ Physica A 234 (1996) 89-107
92
2. General derivation of kinetic equations Let us consider a system of identical classical particles with mass # in volume V. We shall describe the system in the grand nonequilibrium canonical ensemble, so that the number of particles N is not fixed. A nonequilibrium state of the system is completely defined by the set N =- 1,2 .... of the N-particle distribution functions p(Xl,X2 ..... XN; t) =- p(x N, t) which satisfy the Liouville equation
-~+LN
p(xN;t)=O.
(1)
Here, xi =-- (ri, vi) denotes the coordinate ri and velocity vi of particle i, LN = ~iUl L(i)-~Ni< j L(i,j) is the Liouville operator, where L(i) = viO/Ori, L(i,j) = 1/# O~ij/Ori (c~/~3vi- O/Ovj) and ~ij - cb([ri- rj[) is the interparticle potential of interaction. The functions p(xN; t) determine the density of probability to find the system with N particles in the point x u of phase space at time t, so that ~ ucxz = l f p ( x N ;t)dFu = 1, where diN = dxN/N!, dx N =~ d X l " " dxN and dxi = dri dvi. We note that p(x N, t) are symmetric functions with respect to change xi ~ xj of the phase variables for arbitrary pair of particles. The formal solution for (1) is of the form
p(xN; t) = e x p ( - ( t
-
to)LN)p(xN; to).
(2)
However, the following problem appears. How to choose the initial time to and the initial value p(x N, to) of the distribution function? The NSO method gives a solution of this problem by introducing the boundary condition [30] lim
lim exp(toLN)[p(xN; to) - pq(x N, to)] = 0,
(3)
is the quasi-equilibrium statistical operator. According to the hypothesis of an abbreviated description [30, 32], in order to construct this operator, we restrict ourselves a priori to a certain set of the slowest physical quantities which are supposed to determine the nonequilibrium state. In the usual Bogoliubov's scheme of weakening correlations, the function p q w a s taken as a product of one-particle distribution functions fl (xl; t) = (hU(xl; (xt) N))t, where where
pq
I~N(xs" (xt.N. ~ = 1 N s ..... ~.1 Z
~_[ - - •(Xk --Xj,(t)) j~¢...~j~ k=l
(4)
is the s-particle operator of phase density and (...)t = ENC~=I f . . . p(xN; t)dFu denotes the nonequilibrium averaging. In the modified scheme [29], the density of interaction energy e(rl;t) = (e(rl;(xt)N)) t is considered as an additional independent parameter of the abbreviated description, where
~(rl;(X') N) = ~
dr2
dVl
dr2 ~1212N(x2;(xt)N).
(5)
I.P. Omelyan, M. V. Tokarchuk/ Physica A 234 (1996) 89-107 The set of quasi-equilibrium distribution functions be found, provided that the informational entropy
S(p) = - ~
pq ~ (pq(xl; t ) , pq(X2; t) .... )
/pq(xN;t) In pq(xN;t)dI'N
93 can
(6)
N=I is an extremum with the constraint that the mean values of the dynamic variables fiN and ~ are fixed [28], i.e.,
(?~N(xl ; (X! )U ))q = fl (Xl; t),
(7)
(e(P'l; (xt)N))q -~ e(rl; t)
(8)
and preserving the normalization condition ~-~N=]fpq(xN; t)dFN = l, where (...)q = ~'~N~=l f . . . pq(xN; t)dFN denotes the quasi-equilibrium averaging. Solution of the conditional extremum problem 6S(p)/Sp = 0 is of the form [32]
pq(xN;t)
- i=] 2(xi; t)
= z-l(t)exp
- Zi
(9)
where 3ij = ½(fl(ri;t) + fl(rj;t)) and the Lagrange multipliers 2(x;t),fl(r;t) and Z(t) are determined from the self-consistent Eqs. (7), (8) and the normalizing condition, respectively. The effect of the boundary condition (3) is equivalent to transition from the Liouville Eq. (1) to the equation with a source in the right-hand part, which destroys the invariance of this equation with respect to the time inversion [33]
-~ + LN p(xN;t) = --e(p(xN;t) -- pq(xN;t)),
(10)
where the limit e ~ +0 is taken after the thermodynamic limit V ~ w while calculating the average values. Let us integrate Eq. (10) over dFu-s at certain N. Then, taking into account that p is equal to zero on the boundary of the phase space and performing the summation over the number of particles N, we obtain the BBGKY hierarchy with the boundary condition
( ~ + Ls) f~(xS;t)-~-~ f L(i,s + l)f,+l(X~+';t)dx~+] i=1 = -e(f~(xS; where
f~(xS; t)
t) - fq,(XS; t)),
= (fiU(x~; (x')N)) t and
(11 )
fqs(xS; t) = {hUCx , s , s", ,(x!~N~ I ,!qt
are the s-particle nonequilibrium and quasi-equilibrium distribution functions, respectively. The quasiequilibrium distribution function can be presented in a more explicit form f q ~ (xS" ,t)
=
Gs(r';t) ~ I f](xi;t) i--I
(12)
94
LP. Omelyan, M.V. Tokarchuk / Physica A 234 (1996) 89-107
through the s-particle quasi-equilibrium function in coordinate space ~s(rS; t)
Gs(rS;t)- I-Ii=l ~l(ri;t)
(13)
S
where ~s(rS;t) = z - l ( t ) Z
/drN_,exp
N=s+ 1
i=s+ 1
, ~ ( x i ; t ) - Z flij~PiJ i
(14)
is the auxiliary function. It can be shown easily that the dependence of 2(xl;t) on velocity is determined by the one-particle function only: 2(Xl;t) = In ~l(rl; t ) - l n fl(xl;t); therefore, correlations in the space of velocities are not taken into account in the quasi-equilibrium statistical operator (9). However, the absence of these correlations in the infinitely distant past (see Eq. (3)) does not make essential limitations for the system, because, for an arbitrary moment of time they can appear as a result of interactions. Account of the configurational correlations leads to automatical satisfaction of the conservation laws for hydrodynamic variables (number density of particles, momentum and total energy). This is important for the construction of kinetic equations for dense gases and liquids, where the small parameter is absent. In the "pair collisions" approximation we put f3 = 0 for the second equation (s -- 2) of the hierarchy (11 ). Obviously, this is not in the spirit of the two-particle approximation as in the Boltzmann theory, because the essential part of many-particle correlations is taken into account implicitly by the factor G2(rl,r2;t). Then the second equation becomes simpler: ~-~ + L2 + e f2(xl,x2; t) =
eG2(rl, r2; t)fl (Xl; t)fl(x2; t)
(15)
and its formal solution is of the form 0
f2 ( X2 ; t) = e /
dz exp (( e + L2 )Z) G2(r l, r2; t + z)fl (xl ; t + z)fl (x2; t + z).
--oo
(16) According to Abel's theorem, this solution can be rewritten as follows:
f2(x2;t)=
lim
z---} - - ( x 3
exp(L2z)G2(rl,r2;t + z)fl(xl;t + z)fl(x2;t + z).
(17)
Substituting (17) into the first equation (s = 1) of the hierarchy (11), we obtain a general form of the kinetic equation
(-~ + L(1)) fl(Xl;t)= I(x,;t),
(18)
I.P. Omelyan, M. K TokarchuklPhysica A 234 (1996) 89-107
95
where I ( x l ; t ) = f d x 2 L ( 1 , 2 ) lim exp(L2z)G2(rl,r2;t + z) J × f l ( x l ; t + z ) f i ( x 2 ; t + z) "~----+ - -
O~
(19)
is the collision integral.
3. Multistep potential of interaction
The collision integral (19) is still sufficiently complicated to be written in an explicit form for an arbitrary potential of interaction. As is now known [14, 18], only in two particular cases, namely, for the hard-sphere and square-well potentials, this integral has been reduced to an analytical form. In order to extend those previous results, we consider the interparticle potential in the form of a multistep function ~Pij = ¢~s = lim ~ij~0(rij),
(20)
where
eo
/ 80,
(I)ij(rij) =
f,k,
t O,
rij < ¢70, Gk-I < rij < ak,
k = 1.... , N *
(21)
ao* < rij
and N* is the total number of walls except that of the hard sphere. The potential with a hard core contains the strong singularity (operator L ( i , j ) is badly defined). That is why we shall deal with the potential #ij,0(rij) and only in final expressions we shall put e0 ~ ~ . Let us separate for convenience systems of attractive and repulsive walls. Let n* be the number of repulsive walls at the distance ari - ai and of the height Ae~ = e~ - e~+l > 0, where e~ -= ei, i = 1. . . . . n*; and m* be the number of attractive walls with the corresponding parameters aaj =- aj+,,., Ae~ = e~+1 -e~ > 0, where e~ = ej+,,., j = 1..... m*; en.+l =e~, e m a . + l = 0 , N* : n* + m*. Thus, the parameters a0 (hard r * ,O'aj and A~ completely determine the geometry of sphere diameter), n*, O'ri,Aei,m the multistep potential (see Fig. 1, where a specific case n* = 1, m* = 3 is designed and the multistep function approximates some smooth real potential). The limit r ~ -cx~ in the collision integral (19) is a formal one in the mathematical sense. On the physical level of description this means that ITI ~>~0, where z0 > 0 is some characteristic interval of time. For different values of r0 we can investigate different stages of evolution for the system. In the Boltzmann theory [34] [z] ~zc, where Zc is the time of binary interactions. On the other hand, I~1 is far less than a characteristic scale of time Zm for hydrodynamic variables. The above situation is possible because of the fact that for dilute gases kinetic (z ,-~ %) and hydrodynamic ('C '~ "~m) stages of evolution are far from one another in time, i.e., Zc ,~'cf ~Tm, where zf is the characteristic time of free motion. The pattern is qualitatively different in dense gases and liquids with a realistic smooth potential of interaction. The kinetic and r
I.P. Omelyan, M. V. Tokarchuk / Physica A 234 (1996) 89-107
96
EO
E1
E2
li\l ~0
ffl
y..--_...~..-E4 t ...........
if2
E3
G3
~4
Fig. 1. The multistep potential at n* = 1 and m* = 3. The dashed curve corresponds to a smooth real potential.
hydrodynamic stages appear to be closely connected. Besides, such quantities as the length of free motion and time of interaction are not defined in the usual way, because all particles influence the dynamics of interaction for some chosen pair of particles. However, for special types of potentials, such classification of times remains valid even for high densities. For the multistep potential the region t2 of binary interactions consists of the following subregions [ t r k - Aro, ak ÷ Aro], k = 0 , 1 , . . . , N * , where Aro ~ +0 due to the singular nature of the potential under consideration. This potential has the finite range max{ak} = trN. > 0 of action. We can introduce the following set of specific time intervals: c0 -- Aro/9o ~ +0 is the time of pair interactions on the walls, Az = min{ak - a k - l } / 9 o > 0 is the time of motion between the two nearest neighbouring walls and T = aN*~90 > 0 is the time of motion of whole systems of walls for some pair of particles, where 90 is the averaged relative velocity of two particles. Because the potential contains the horizontal parts, where the force of the interparticle interaction is equal to zero, it is possible to introduce also time .t-f as an average time of free motion of particles in the system. Changing geometry of the potential and increasing density, we can make this time arbitrarily small in order to fulfil the relation Zo~f~Az
< T<
"Cm .
(22)
The kinetic stage of the evolution corresponds to small times of order z ~ zf. Therefore, as far as the relation (22) is satisfied, the formal limit z ~ -cx~ must be considered as z/z0 ~ - ~ or merely as I~1 ~>~0. For smooth nonsingular potentials the region of interaction has some nonzero size and z0 is finite. For singular potentials, the case
L P. Omelyan, M. K Tokarchuk / Physica A 234 (1996) 89-107
97
z0 -~ + 0 is possible and the limit Z/Zo ~ - o o remains valid even for z -~ - 0 if only z0 is a value of the higher order infinitesimal
lim { lim+oZ/Zo}--~-c~.
~----~-- 0
(23)
"C(
Thus, the formal limit z --~ - 0 can be applied for the collision integral (19) in the sense of (22) and (23). As is now well established, the limit z -* - 0 leads to the kinetic equations of the RET [23,31] and RDRS [27] theories• The limit r ~ - 0 allows one to ignore the time retardation in (19). By the construction, G2 depends on time only functionally through the local density n(r; t) = f fl(x; t ) d v and "reciprocal temperature" fl(r; t), i.e., G2(rl, r2; t) = G2(rl, r2 [n(r; t), fl(r; t ) ) .
(24)
Then lim G 2 ( r l , r 2 ; t + z ) = ~-,-o
G2(rhr2;t)+
ao2
limz If ~--,-o [ J J
}=
dr [ 6(72 an(r;t) L&(r;t) at O2(r
(25)
,r2;t)
because during the time z ,-~ %f <~1:m hydrodynamic parameters n(r; t) and fl(r; t) could not essentially change. Besides,
(
[
f l ( x l ; t + z) = e -L(1)r f l ( X l ; t ) + z L(1) + ~
f[(Xl;t)
}
,
(26)
where 1
p
T
(xl; t) = - ] dz' eL(l)*'f l (xl; t) TJ 0
(27)
are the average values of the one-particle function during the time interval z. But • -3 llm,~_0 f l (Xl; t) = fl (xl; t) and
~li~m-of~(xl;t + z) = ,~-01ime - t O ) * f l ( x l ; t )
(28)
because it follows from (22) that during the essential part of time z the motion of particles is free and during the very short interval z0 they interact. Noting that fl(x; t) is a continuous function of r, we obtain
af, lim e-L~l)~ f l ( X l ; t ) = f l ( x l ; t ) -- ~lim_ --0 "cvl ~Crl = f l ( x l ; t ) .
(29)
T~+ - - 0
As a result of (25), (28) and (29) the collision integral (19) transforms into the form l(Xl;t)
f d x 2 L(1,2) lim_o et2~ Gz(rt, r2; t ) f l (xl; t ) f l (x2;t).
(30)
We consider the dynamics of interaction for two particles in the centre of mass system. Then L(1,2) = - I/p* d a ~ f f a r l 2 0 / 0 9 , L2 = -g3/art2 - L(1,2), where rj2 =
98
L P. Omelyan, M. V. Tokarchuk l Physica A 234 (1996) 89-107
i'1 -- r2, 0 ~- V2 -- !~1, ~ -~ r12/r12 and #* = #/2 is the reduced mass. Integrating (30) over dr2 in the spherical system of coordinates with the centre coinciding with the centre of the first particle yields oo
l<,r,&i
I =
2 1 ^'+12 _~_~ eL2*G2(rl, r2, t )fl (Xl; t)fl (X2; t), d g r 1 2 7 f f Orl2 09 ~lim -*-0
o (31) where d6 is the element of the surface sphere of unit radius. According to the definition (21) of the multistep potential we have ( ' 0 -- ' ] ) ~ ( r 1 2 -- 0"0) -- Z AI3~(r12 - ffri) -t- Z i=1 j=l
Or12 --
Ae~b(ri2 -
O'aj)
(32)
and, therefore, it is quite sufficient to perform the integrations in (31) over the coordinate space only inside the region 12 of interaction. Taking into account that fl is continuous in an arbitrary point of space and G2 is continuous over r12 out of the region ~, one can show that the following equalities take place:
i dr12L(1,2) lim eL2* G2(rl,r2; t)fl(Xl; t)fl(x2; t) ~-*-0
f2
= -- [ drl2 9_~._0 l i m e L2~G2(rl, r2;t)fl (Xl ;t)fl (x2; t) d Q
= -S
or12 z---+-O
d D 9 ,-~lim-oeL2~Gz(rl, r2; t)fl(Xl; t)fl(x2; t),
(33)
D
where the theorem about the gradient has been used, so that D is the surface, spanned on the volume f2. In view of Eqs. (32) and (33), we can write the collision integral for the multistep potential in the form N*
I(xl; t) = Z I k ,
(34)
k=0
where
lk=
lim z--+--O dro-*+O
a
27 S d9
dtr trgeL'~ G2(rl,r2;t)fl(xl;t)fl(x2;t)
'ri2=ak-- Ar°
(35)
rl2=trk+Ar°
and I0 = Ins is the hard-core part of the collision integral, whereas If = Ii and 17 = Ij+n. are the parts caused by interaction on the ith repulsive and jth attractive walls, respectively.
99
I.P. Omelyan, M.V. TokarchuklPhysica A 234 (1996) 89-107
Owing to the condition (22) and the evidence of the limit r ~ - 0 , the collision integral (34) breaks down into several parts. Physically, this means that during the time interval of order of free motion zf, two particles can interact only at one of the walls of the multistep potential. Interaction of same pair of particles during this time interval in several subregions f2 is impossible. Moreover, if one considers a separate pair of particles, it appears that in the velocity space their interaction at one of the walls does not depend on the interaction at the remaining ones. It is caused by the fact that during the time of motion between neighbouring walls, a pair of particles interacts many times with the medium and velocity correlations weaken. So, during some interaction of the pair of particles at a potential wall, this pair does not "remember" its velocities after the previous "meeting". A large part of the configurational correlations is conserved and is taken into account by the factor G2. The property of additivity for collision integrals will break, if condition (22) is not fulfilled. In this case, one should consider complex interference phenomena. As we can see from Eq. (35), all parts of the collision integral have the same structure. In order to obtain analytical expressions for them, it is necessary to consider the action of the two-particle operator of evolution l i m ~ _ 0 e z2~ at the external r12 tTk ~-Aro and internal r12 = a k - Aro surfaces of the sphere of radius ak. The results for Ih~ and I] are known. For the hard sphere wall the collision integral is [23] :
I h s = O -2
//
d~[~910(~9)[F2(rl,V'l,rt + a 0 a , v2;t )
dg
-F2(Xl, rl - a+ 6, v2; t)],
(36)
where v~=vl+6(69),
l = 1,2
(37)
are the velocities of a pair of particles after reflection from the hard sphere wall, F2(xl,x2; t) = G2(rl,r2; t)fl(xl; t)fl(x2; t) ~ fqz(X2; t), o-± = limzr0~+o(a 4- Aro) and O denotes the Heaviside function. In the case of a finite attractive wall we have [27]
I; =
jfdo/d +0(-69
lagl{ --
O(tTo)[F2(¥1, Vatl,rl-{-'°'aj~ vtt" .2, t ) - F 2 ( x .
r~ - ~ ,
v2; t)]
a))[F2(rl,Val,rl +a+.~, Va2 Ii,., t) -- F2(Xl, rl -- O'aj~,V2; t)] tit
+O(--~g)O(c~g + o~j a )[F2(rl, vl,rl ' + O'aj~ - ^ 102, " t ) - - F 2 ( x l , r t - - a ~ j f f , v2;t)]}
(38) where ,,
10al -~- ~l -4-
+
l = 1,2
(39)
1.P. Omelyan, M.V. Tokarchuk/Physica A 234 (1996) 89-107
100
are velocities of particles after leaving the region of interaction, whereas '" : I)ai
1~i -4-
½e(eg + v/(eg)
2
-- (~aj)'
(40)
i = 1,2
are velocities of particles after entering the region of interaction and aaj =
~/..
Considering the dynamics of two particles at a finite repulsive wall, we have obtained
17 =a2ri / dg / d~l~gl{ O ( - ~ 9 ) -[- ^
If,
>[F2(rlVtrtl,rl q- (TritT, Vr2,t ) -- F2(Xl,r
I - o-~Tt~, v2; t)]
q-O((rg -- (~[)[F2(rl,vrl '" ,rl q- O'rT.(7, vr2, m. t /-~- F2(Xl, rz - a+t~, r2; t)]
+O(t~g)O(-~g + ¢~)[F2(rl, v'1, rl + a+ ~, v'2; t) - F2(Xl, rl - tr+tr, v2; t)]}. (41) Here I~ri = 13i d:
+ ~ri),
i = 1,2
(42)
are velocities of particles after leaving the region of interaction, where 0~ri VaAs~/[,I. This process does not require additional conditions in 9 and takes place for arbitrary 89 < 0. At the same time ----
'" lJri
=-
Vi 4-
}e(eg + v/(eg)
-
-
O~ri),
i = 1,2
(43)
are velocities of particles after entering the region of interaction. This process is possible when the relative velocity of two particles is sufficient to pass the potential wall. Otherwise, when the relative velocity is too low, velocities of particles change during interaction as after reflection (see Eq. (37)). Therefore, except for the terms Ihs (RET type) and I] (RDRS type), we have obtained new terms I r (41) which are caused by contribution of repulsive walls of finite heights into kinetics. Obviously, the repulsive wall of a finite height is a generalization of the infiniteheight hard sphere wall. Taking into account (13), (14) we have lim
_ = 0.
G2(rl,r2;t)
ZIg~--+OO
(44)
r12 =o'ri
This means that particles cannot penetrate the infinite-height wall. Then it is easy to verify that the term I r (40) at Ae[ --+ oo and ari =-=-ao reduces to the hard core type (36), i.e., /hs = Ae~---+oo lim /ir -,,-o0"
(45)
The collision integral (34) can be expressed in a more compact form
1(xi; t) = f dx2 "PjzG2(rz, re; t)fl(Xl;t)fl (xz; t) J
(46)
LP. Omelyan, M.V. Tokarchuk / Physica A 234 (1996) 89-107
101
in terms of the two-particle pseudo-Liouville operator n*
m*
(47)
E ?:+E Z ?; j = l p=b,c,d
i=1 p=b,c,d
which consists of the corresponding pseudo-Liouville operators
?hs =
~_ ^
^a
o-20 d(~ [(~g[O((~g){a(r,2 4- o-o o-)B ( 6) - 6(r12 -- O'g-a)}
(48)
for the hard sphere wall, ?rf = o-ri2
d6
lp^
[(~olOP(~g){a(rt2 +
o-r;
^P ^ 2p^ ° ' ) B r i ( ° ' ) - ~(¥12 - o-r; o-)}
(49)
for the ith repulsive wall and
J dd I@lO~.(@){5(r,2 +
^ = o-aj f 7;P 2
l , p ^ ^P
o-aj o-)Baj(o-)
^
2, p ^
- a(l*12 - o-aj (7)}
(50)
for the jth attractive wall. Here, o~,( @ ) = o ( - @
o~,(@)
)
o~.(@) = o(@)o(~
=
o ( @
-
~[),
Oabj(@) = O(@),
- @),
O~j( @ ) = 0 ( - - @ ) 0 ( 4 + @ )
O~aj(@ ) = O ( - - @ - ~ ) ,
(51)
are the corresponding Heaviside functions, l,b
Gri
1,c
+
= o-ri, -
o-ri ~ o-ri , I,d
¢7ri
+
~- o-ri,
2,b o-ri = O'r/, 2,c
o-ri
2,d
Gri
+
1,b o-aj = o-aj, 1,c
+
~ o-ri,
o-aj = o-aj,
+ = o-r/,
o-aj
1,d
~- o-aj'
2,b a~ 0"aj ~ o- , 2,c -O'aj ~ o-aj, 2,d
-
o-aj ~ Gaj
(52)
and /~P are operators of velocity displacement caused by the interaction of p-type at each of walls, namely, /~a(~)~(Vl,l~2)
= ~fi(~)~(Vl,V2)
= /~dj(~)~(Vl,P2)
~bri(~)~(Vl,V2
) = aa'a" 11 11 ' " (Vrl, Vr2 )
B^b / ~ ) . . ~ ( V ~ , V 2 )
hci(~)~(Vl,V2)
~ (Vrl11!, Vr2 lit ) , = ~/~
/~cj(~)~(I)l,
= ...~(Vl,l~2) , ~ t- ,
~ (I)a,,1, I~a2 ,, ) , = "~
!~2 ) ~__ ~----------------~r. (Val ,ttt Va2ttt ) ,
(53)
where ~ is an arbitrary function of velocities. Physically, the values o-ri 1,p,to-aj 1,p,) (52) correspond to the distances between the particles just after the interaction of p-type on the ith repulsive (jth attractive) wall of the potential, whereas the values o-ri2,p.( o - a2j, )p . correspond to the distances just before the collision.
1.P. Omelyan, M. 1I.. Tokarchuk/ Physica A 234 (1996) 89-107
102
4. H-theorem
Let us show that the kinetic equation obtained satisfies a global H-theorem. In order to do this, we introduce the operator 7~+ conjugated to 7~n by the transformation
S avisdx'~(xi'xllt)'12dk(xl'x';t)=i dVlidxldC(xi.x,;t)'£~(xl.x,;O. (54) where (p and ~b are arbitrary functions of xl,x2. By straightforward calculation we obtain n*
^
i=1
m*
p=b,c,d
j=l
^p+
(55)
p=b,c,d
where
"I~ = f~glO(~g)f(r12 - a+)(Ba(~) - 1),
^P+ "~ I~alOg (#a),~ ( r12 Tri aj aj
- O ' r2,p~ i aj/
(56)
(l~Pj(~)~.~:~p -aj
~a2r~P~ • aj/
(57,
Here ~ is the operator of coordinate displacement. It acts on the relative distance between particles as follows: ~ r il(2)'P ~ ( 1r2 ) = ~ ( O'lr~2)'p I ' aj \ aj /
(58)
where q¢ is an arbitrary function of rn. Our task is to prove the H-theorem, i.e. the property dS(t) _- t)S(pq)>~0. dt St
(59)
Using (6), (9) and (14) we obtain t3S(pq) _ S(1) q- S(2) q- S(3) '
(60)
St where S (l) ----
dxl In ~l(rl;
. t3e(rl; t) S (3) = fJ d r l fl(rl;t) -~ •
it)
S (2) : -
dXl l n f l ,
St
(61)
Let us consider each term separately. Taking into account the identities (14), (54) and the kinetic equation (18) with the collision integral (46), we have for S (1) in the
I.P. Omelyan, M. V. Tokarchuk / Physica A 234 (1996) 89-107
103
explicit form
S(1)=-lidxlidx2F2(xi'xz;t){'i2vlcOi~<) 2 +8"9[~$(r12- "°v) "
+ ~
(b(rl2 - a{) - fi(rl2 - ak-))
])
(62)
.
k=l For the second term we find
S(2)= - I-2S dx,S
,/ i
= -2
dXl
dx2F2(Xl,X2;t)'F~In(fl(xl;t)fl(xz;t))
{o
dxz F2(xl,x2; t)grg
n*,m*
+Z
Z
(69)6(r12- o~) In Ba(gr)fl(xl(Xl, ;t)fl(xz;flt)fl(x2; t) t)
OPri(ag)6(ri2-aZri:') In^p fl(xl;t)fl(x2;t)
i,j=l p=b,c,d
aj
ali
B ri ( ~ r ) f l ( x 1 ; t ) f l ( x 2 ; t)
aj
(63) In order to calculate the third term it is necessary to obtain the equation of motion for the density of interaction energy. Comparing the first equation (s -- 1) of the BBGKY hierarchy (11) with the kinetic equation, multiplying the second (s = 2) equation of the hierarchy on 1/2~lz and integrating over dr1 dx2, we find, taking into account (5), (8), that
Oe(rl;t) cq---~ + ~i f
dr:
f dx2vl~-Vlr2tx:;x2;t ~,2.. )
= :'Id"Sd"O""'"'":":'S : dr: I dx2 F2(x:,x2; t)/~+CPl2 " .
2
i=1 m*
-~9(-6g)6(r12-
~r?)] + ~
As~[O(-~9 -
°(j)6(r12 - o"~.)
j=l
-- ~9(aO)a(r12 -- ~y~)] } "
(64)
104
LP. Omelyan, M.V. TokarchuklPhysica A 234 (1996) 89-107
Then the third term is of the form
S `3) ~---~
fdx2F2(x~,x2;t)
,221)10fl(rt; t)&_____~ +/~t218.gl
dXl
A~r[o(8.g - ~xr)~(rl2 - fir+) - O(-8.g)~(r12 - O-r~)]
×
L i=1 + ~
A 4 [ O ( - 8 . g - ~ ) a ( r , 2 - ~2j) -
0(8.~)6(r~2 -
~-)1
•
(65)
j=l Now we add the three terms (62), (63) and (65). Taking into account the inequality > 0 (the equality takes place only at h = z) and, since, the one-particle distribution function is positively defined, after some algebra one rewrites (60) in the form
ln(h/z) > 1 -z/h;h,z
n* m* d~(tt) ~Jhs ~- E Jri-Jr- E Jaj, i=1 j=l where
l/ /
Jhs = ~
dXl
(66)
dx2 (~(r12 - o'~){-8.g + 18.g [0(8"9)(1 - Ba(8.))}f2(xl,X2; t), (67)
Jri =
,/ dxl / dx2
"~
O'r; ) - ($(r12 - ar2[b))+
-
18.ol
o~(~a)
p=b,c,d
2p ^P ^ P / X~(rl2 - ~ri' )(1 - Bri(~r)exp(-flnEri)) f
Jaj = ~1 /
dxl
f
dx2
{-~g(5(r12
2p X ~(r12 - era) )(1 -
_
F2(Xl,X2;t ) ,
(68)
~)2 b ) _ 5(r,2 - (~'f))+
[Sg[ E oP(~g) p=b,c,d
^p . p } F2(xl,x2;t) Baj(cr)exp(-flnEaj))
(69)
,I and E~'~/= - ~ , ~ c = a~r, e g = 0, E~j = -A~), E:j = a ~ , E~ = 0. The term Jhs = 0 just as in the proof of H-theorem in the RET theory [25]. Considering the terms Jri,Jaj, one should use the fact that G2(r~,r2;t) is a discontinuous
105
I.P. Omelyan, M.V. TokarchukIPhysica A 234 (1996) 89-107
function. It is caused by the discontinuity of the initial potential ~12. It follows from (9), (12)-(14) that
/
drl2 F 2 ( x b x 2 ; t)O(r12 - O'+) r = exp(/312Ac r ) drl2 F2(Xl ,x2; t )6(rl2
- ~r+.) =
i
exp(-fllzAey)
drl2F2(Xl,X2; t ) ~ ( r l 2 - o-ri ) ,
fdr,2F2(Xl,X2; t)6(r12 - o-~). (70)
Then, one can show easily that each of the terms (68), (69) is equal to zero, i.e., Jri = 0, i = 1,n*; Jaj = 0, j = 1,m*. Then from (66) it follows that dS(t)/dt>,O and the H-theorem is proved. Let us determine the one-particle distribution function jit°)(xl;t) for which dS(t)/dt = 0. According to (66)-(69) such function must fulfil the following equations: ~(~g)6(r12 -- o-O~)[1 - - / T a ( r ~ ) ] f ( ° ) ( X l ;
t)f~°)(x2; t) = O,
OPri.(69)6 (rl2 - a2~ip) [exp (flizEPr}) -B~iy(~) 1 f(°)(xl;t)f~°)(xz;t) = O . ay
ayl
(71) By the straightforward calculation we can verify that the conditions (71) are fulfilled simultaneously at
l;°,(x,;,) _-,(,,.,)' k,('(')7 21t /#
exp
{-9(,)(o, - ,,))2}
,
where u = fvfl(x;t)dv is the hydrodynamic velocity. In the equilibrium spatially homogeneous case n,/7 = const., u = 0 and f(0) tends to the Maxwell distribution. In the presence of external field W(r) the kinetic equation transforms into the form (0
+L(1)
lOW(rl) ~ )
#
Orl
cGvi
fl(xl;t)
---- S dx2 l'12F2(Xl,X2;t),
(73)
i.e., only flow term is changed, while the collision integral remains the same (in the limit z ~ - 0 the external field does not change velocities of particles during collision). Putting fl - f(0) at n =- n(r), /3 = const., u = 0 and integrating (73) over drl, we obtain the integral equation for local density in an equilibrium state
± ~grl
[ln n(rl ) +/3 W(rl )] =
x
6(r12 - G~-) + Z ( 1
fdr2~n(,2
)Gz(rl, r2 )
-exp(-/312A~[)6(rl2
-- o-ri+)
i=1
+ Z j=l
where G 2 ( r l , r 2 )
(1 - exp(/312A~) ~(r12 -- O-+)
}
,
is the equilibrium distribution function.
(74)
106
I..P. Omelyan, M. II.. Tokarchuk/Physica A 234 (1996) 89-107
5. Conclusions We consider some special cases of choosing for the parameters of the model multistep potential (21). Choosing n* = m* = 0 corresponds to the hard-sphere potential. Here I[ =- 0 and I~ - 0, i.e., I -= Ihs and the kinetic equation is the equation of the RET theory [23]. In the case when n* = 0, m* = 1 and Ae~ ~ 0, the initial multistep potential transforms into the square-well one and I = Ihs + IF coincides with the collision integral of the RDRS theory [27]. It is interesting also to consider the case, when the initial potential transforms into a some smooth potential @t(r12) at r > ~0, i.e., LI~[ ~ O,
Ae r ~ O ,
i = 1..... n*; n*---*oo,
A~a ~ O,
Aej --~ 0,
j = 1,...,m*; m* --* c~,
(75)
AF.j--~ Cg@____L",2=~., t
(76)
at the condition -
-
-
-
--.)
O(IJt
0r12
r,2=ari
'
A¢~a
¢3r12
'
where Aa r = Aai, A ~ a = AM+,* and Aak = ak - ak-l. In such situation it can be shown that the collision integral transforms into the form I = Ihs + IMF, where IMF =
O~t ~F2(xl,x2; t)
dx2 t~rl
t~vt
(77)
coincides with the mean-field term of the KMFT theory [28]. The last case, however, does not fulfil the condition (22) for geometry of the interparticle potential ( d z becomes too low). As a result, the limit T --+ - 0 leads to the term IMF without collisions. At z ~ - 0 direct contribution of collisions into kinetics is provided only for singular potentials, when the distances between neighbouring walls are sufficiently big, and for high densities. On the other hand, besides the error due to the fact that Eq. (22) is not fulfilled for smooth potentials, there is also a model error caused by the fact that the multistep potential differs from real (as a rule smooth) potentials. As one can see, these two errors have opposite tendency. Thus, in order to use the obtained kinetic equation for real systems, it is necessary to find some compromise solution in the problem of geometry for the multistep potential. This problem can be solved as follows. Let the transport coefficients be as already obtained on the basis of the kinetic equation. Then for an optimal representation of a real potential, parameters of the multistep function can be accepted as those that give the best fit to the transport coefficients in a whole region of temperatures and densities. This will be the subject of a separate study. The condition (22) is not fulfilled when density decreases. In the limit of dilute gases the time of flee motion becomes too big: zf ~ T ~ d z and this fact does not give the possibility of the kinetic equation under consideration to transform into the Boltzmann one. Therefore, the kinetic theory proposed in this paper can be applied exclusively for dense gases and fluids.
I.P. Omelyan, M.V. Tokarchuk/Physica A 234 (1996) 89-107
107
Acknowledgements The authors warmly acknowledge stimulating discussions with Prof. D.N. Zubarev and Prof. V.G. Morozov concerning this work. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34]
N.N. Bogoliubov, in: Studies in Statistical Mechanics I (North-Holland, Amsterdam, 1962). J.R. Dorfman and E.G.D. Cohen, J. Math. Phys. 8 (1967) 282. J.G. Kirkwood, J. Chem. Phys. 14 (1946) 180. J.G. Kirkwood, J. Chem. Phys. 15 (1947) 72. M. Born and H.S. Green, Proc. Roy. Soc. A 188 (1946) 10. E.G.D. Cohen, Physica 28 (1962) 1045. E.G.D. Cohen, J. Math. Phys. 4 (1963) 183. L.S. Garcia-Colin and A. Flores, J. Math. Phys. 7 (1966) 254. G. Sandri, Ann. Phys. (N.Y.) 24 (1963) 332, 380. E.A. Frieman, J. Math. Phys. 4 (1963) 410. K. Kawasaki and J. Oppenheim, Phys. Rev. 139 (1965) 1763. J.R. Dorfman, in: Lectures in Theoretical Physics (Gordon and Breach, New York, 1967). E.G.D. Cohen, Physica A 194 (1993) 229. S. Chapman and T.G. Cowling, The Mathematical Theory of Non-Uniform Gases (Cambridge University Press, Cambridge, 1952). H.G.M. Hanley, R.D. McCarty and E.G. Cohen, Physica 60 (1972) 322. J.A. Barker and D. Henderson, J. Chem. Phys. 47 (1967) 4714. J.D. Weeks, D. Chandler and H.C. Anderson, J. Chem. Phys. 54 (1971) 5237. H.T. Davis, S.A. Rice and J.V. Sengers, J. Chem. Phys. 35 (1961) 2210. R.C. Castillo, E. Martina, M. Lopez, J. Karkheck and G. Stell, Phys. Rev. A 39 (1989) 3106. J. Karkheck and G. Stell, J. Chem. Phys. 75 (1981) 1475. R.C. Castillo and J.V. Orozco, Int. J. Thermophys. 11 (1990) 1025. M. Grmela, R. Rosen and L.S. Garcia-Colin, J. Chem. Phys. 75 (1981) 5474. J. Karkheck and G. Stell, Phys. Rev. A 25 (1982) 3302. H. Van Beijeren and M.H. Ernst, Physica 68 (1973) 437. P. Resibois, J. Star. Phys. 19 (1978) 593. P. Resibois, Phys. Rev. Lett. 40 (1978) 1409. J. Karkheck, H. Van Beijeren, I. de Schepper and G. Stell, Phys. Rev. A 32 (1985) 2517; see also H. Van Beijeren, J. Karkheck and J.V. Sengers, Phys. Rev. A 37 (1988) 2247. G. Stell, J. Karkheck and H. Van Beijeren, J. Chem. Phys. 79 (1983) 3166. D.N. Zubarev and V.G. Morozov, Teor, Mat. Fiz. 60 (1984) 270 [in Russian]. D.N. Zubarev, Nonequilibrium Statistical Thermodyamics (Nauka, Moscow, 1971) [English transl.: (Plenum Press, New York, 1974)]. D.N. Zubarev, V.G. Morozov, I.P. Omelyan and M.V. Tokarchuk, Teor. Matem. Fiz. 87 (1991) 113 [in Russian]. D.N. Zubarev, in: Itogi Nauki i Tekhniki, Sovr. Prob. Mat., Vol. 15 (VINITI, Moscow, 1980) [in Russian]. D.N. Zubarev and M.Y. Novikov, Teor. Mat. Fiz. 13 (1972) 406 [in Russian]. L. Boltzmann, in: Lectures on Gas Theory (University of California, Berkeley, 1964).
Physica A 234 (1996) 89-107
Kinetic equation for liquids with a multistep potential of interaction. H-theorem I.P. Omelyan*, M.V. Tokarchuk Institute for Condensed Matter Physics, National Ukrainian Academy of Sciences, 1 Svientsitsky St., UA-290011 Lviv, Ukraine Received 12 December 1995
Abstract
A new kinetic equation for dense gases and liquids with the interparticle potential of interaction in a form of the multistep function is proposed. The equation is consistently derived on the basis of Zubarev's nonequilibrium statistical operator method using a modification for weakening correlations principle which takes into account local conservation laws. In particular cases this equation coincides with the well known kinetic equations for hard-sphere and square-well potentials. Application of the obtained kinetic equation to real systems is discussed. A proof of a global H-theorem is given.
PACS: 05.20.Dd, 51.10.+y, 82.20.Mj Keywords: Kinetic equation; Dense gases; Multistep potential; Boundary condition; H-theorem
1. Introduction
The construction of kinetic equations for dense gases and liquids is one of the main problems in the kinetic theory of classical systems. The Boltzmann equation gives the foundations for the description of nonequilibrium properties of a dilute gas, but at the same time we do not have an analogous kinetic equation for dense systems with a realistic interparticle potential of interaction. In order to give more rigorous theoretical assumptions than Boltzmann's, Bogoliubov [1] proposed a new approach to the problem of obtaining the kinetic equations. It allows to generalize the resuits to higher densities. In the framework of this method the Boltzmann equation appears as a zeroth approximation in the expansion in powers of density. But if one tries to investigate higher orders of the expansion in density, it appears that terms beginning with the second-order contribution are divergent [2]. The reason for the * Corresponding author. E-mail: [email protected]. 0378-4371/96/$15.00 Copyright (~ 1996 Elsevier Science B.V. All rights reserved PII S 0 3 7 8 - 4 3 7 1 ( 9 6 ) 0 0 2 8 6 - 5
90
I.P. Omelyan, M.V. Tokarchuk/ Physica A 234 (1996) 89-107
divergences lies in the fact that the expansion was based on dynamics of isolated groups of particles in the infinite space, regardless of the medium, i.e., the rest of the particles. Theories of Kirkwood [3,4], Born-Green [5], Cohen [6,7], Garcia-Colin and Flores [8], Sandri [9] and Frieman [10] lead to results which are similar to those obtained by the Bogoliubov's method. In order to obviate the divergences, the region of integration must be cut off in coordinate space. As a result of this truncation, an expansion of the transport coefficients in density contains nonanalytic logarithmic terms [11]. In a more rigorous approach it is necessary to perform a partial resummation of the BBGKY hierarchy [12, 13]. This problem is practically unsolvable for high densities. On the other hand, the solution to kinetic equations for dense gases is a very complicated problem. These equations contain the many-particle scattering operators and this fact does not allow to solve them in general. That is why the construction of kinetic equations for dense gases and liquids with model interparticle potentials of interaction is a problem of great importance. The first step in this direction was made in the semiphenomenological standard Enskog theory (SET) [14] of dense gases. Simple model of hard spheres, used in this theory, was very successful in explaining the dependence of the transport coefficients on density. In order to apply the results of SET to realistic systems, Enskog proposed to substitute hydrostatic pressure of the hard spheres system by the thermodynamic pressure of a real system. Taking this suggestion as the basic one, Hanley et al. [15] elaborated the modified Enskog theory (MET). In this theory the diameter of hard spheres is determined by the second virial coefficient of the equation of state for the real system. Thus, it depends on temperature and density. Using different equations of state, BH (Barker, Henderson) [16], WCA (Weeks, Chandler and Andersen) [17], one obtains the corresponding versions of MET theory. In 1961, Davis, Rise and Sengers [18] proposed the kinetic theory DRS. Within the framework of this theory an interparticle potential of interaction is chosen in the square-well form. An attractive part of the real potential is approximated here by a finite height wall. The case where the potential between the hard sphere and attractive walls is considered as smooth tail (instead of the constant potential in DRS) has also been investigated [19]. In the kinetic mean-field theory (KMFT) [20] a smooth attractive tail is considered along with the hard sphere potential. However, as is now well established, all these potential models do not describe real systems satisfactorily. In particular, as was shown in recent calculations [21], KMFT does not give a correct description of the temperature dependence for the transport coefficients. From the point of view of the statistical theory of nonequilibrium processes, SET and DRS contain two essential drawbacks. The first one is that their kinetic equations are not obtained in the frames of some consequential theoretical scheme and one does not know how to improve these theories. The second one is that the H-theorem is not proved. Nevertheless, not long ago the way of constructing of SET entropy functional was given [22]. In order to overcome these drawbacks Karkheck and Stell [23] have
I.P. Ornelyan, M.V. TokarchuklPhysica A 234 (1996) 89-107
91
obtained the kinetic equation for a revised version of the Enskog theory (RET), using the principle of maximum of the entropy functional. This equation was proposed for the first time by Van Beijeren and Ernst [24] on the basis of the diagrammatic method. Resibois [25,26] proved the H-theroem for it. In 1985, the revised version of DRS theory (RDRS) has been proposed [27]. The kinetic equation of RDRS satisfies the H-theorem as well. In the paper of Stell et al. [28], the authors finished construction of the kinetic variational theory (KVT) introducing a whole family of theories, namely, KVT I, KVT II and KVT III. The KVT III refers to the local energy constraint in maximizing entropy as distinguished from the global energy constraint (KVT II) or hard-core constraint (KVT I). According to this classification, the kinetic equation of RDRS [27] is obtained applying KVT III to the square-well potential and KMFT [20] is derived from KVT I as an application to the potential with a smooth tail. In Ref. [28], the KVT III version of KMFT was proposed. In this case, the quasi-equilibrium binary distribution function (QBDF) of hard spheres is to be substituted for the full QBDF, which takes into account the full potential of interaction. The main conclusion of the KVT I version of KMFT [20] is that the smooth part of the potential does not contribute explicitly to the transport coefficients. This fact is caused by the approximations of the theory. When the full QBDF is used, as in the KVT III version of KMFT, the transport coefficients are determined by the soft part of the potential. Zubarev and Morozov [29] using the Zubarev's nonequilibrium statistical operator (NSO) method [30] have obtained kinetic equations for dense gases by formulation of a new boundary condition for the BBGKY hierarchy. This condition takes into account the correlations connected with the local conservation laws of hydrodynamic quantities. Such an approach leads to results similar to KVT III. In particular, it was shown [31] that within this approach the kinetic equation for hard spheres coincides with the kinetic equation of RET. In the present paper, in order to choose a more realistic model, a kinetic equation for dense classical systems is proposed for the potential of interaction in the form of a multistep function. This function consists of the hard-sphere part and a system of repulsive and attractive walls. The derivation of the kinetic equation is performed on the basis of theoretical scheme of the kinetic theory construction by the NSO method with the help of a modification to the weakening correlations principle [29]. It is shown that the obtained kinetic equation satisfies a global H-theorem. For certain parameters of the model multistep potential, this equation coincides with the kinetic equations of RET, RDRS and KVT III version of KMFT theories. Thus, the kinetic theory proposed in this paper can be considered as a general one, which covers the results of previous theories for particular cases and, consequently, shows a connection of these theories with one another. We also discuss the region of application of the generalized kinetic equation and ways of possible improvement of the theory.
I.P. Omelyan, M. V. Tokarchuk/ Physica A 234 (1996) 89-107
92
2. General derivation of kinetic equations Let us consider a system of identical classical particles with mass # in volume V. We shall describe the system in the grand nonequilibrium canonical ensemble, so that the number of particles N is not fixed. A nonequilibrium state of the system is completely defined by the set N =- 1,2 .... of the N-particle distribution functions p(Xl,X2 ..... XN; t) =- p(x N, t) which satisfy the Liouville equation
-~+LN
p(xN;t)=O.
(1)
Here, xi =-- (ri, vi) denotes the coordinate ri and velocity vi of particle i, LN = ~iUl L(i)-~Ni< j L(i,j) is the Liouville operator, where L(i) = viO/Ori, L(i,j) = 1/# O~ij/Ori (c~/~3vi- O/Ovj) and ~ij - cb([ri- rj[) is the interparticle potential of interaction. The functions p(xN; t) determine the density of probability to find the system with N particles in the point x u of phase space at time t, so that ~ ucxz = l f p ( x N ;t)dFu = 1, where diN = dxN/N!, dx N =~ d X l " " dxN and dxi = dri dvi. We note that p(x N, t) are symmetric functions with respect to change xi ~ xj of the phase variables for arbitrary pair of particles. The formal solution for (1) is of the form
p(xN; t) = e x p ( - ( t
-
to)LN)p(xN; to).
(2)
However, the following problem appears. How to choose the initial time to and the initial value p(x N, to) of the distribution function? The NSO method gives a solution of this problem by introducing the boundary condition [30] lim
lim exp(toLN)[p(xN; to) - pq(x N, to)] = 0,
(3)
is the quasi-equilibrium statistical operator. According to the hypothesis of an abbreviated description [30, 32], in order to construct this operator, we restrict ourselves a priori to a certain set of the slowest physical quantities which are supposed to determine the nonequilibrium state. In the usual Bogoliubov's scheme of weakening correlations, the function p q w a s taken as a product of one-particle distribution functions fl (xl; t) = (hU(xl; (xt) N))t, where where
pq
I~N(xs" (xt.N. ~ = 1 N s ..... ~.1 Z
~_[ - - •(Xk --Xj,(t)) j~¢...~j~ k=l
(4)
is the s-particle operator of phase density and (...)t = ENC~=I f . . . p(xN; t)dFu denotes the nonequilibrium averaging. In the modified scheme [29], the density of interaction energy e(rl;t) = (e(rl;(xt)N)) t is considered as an additional independent parameter of the abbreviated description, where
~(rl;(X') N) = ~
dr2
dVl
dr2 ~1212N(x2;(xt)N).
(5)
I.P. Omelyan, M. V. Tokarchuk/ Physica A 234 (1996) 89-107 The set of quasi-equilibrium distribution functions be found, provided that the informational entropy
S(p) = - ~
pq ~ (pq(xl; t ) , pq(X2; t) .... )
/pq(xN;t) In pq(xN;t)dI'N
93 can
(6)
N=I is an extremum with the constraint that the mean values of the dynamic variables fiN and ~ are fixed [28], i.e.,
(?~N(xl ; (X! )U ))q = fl (Xl; t),
(7)
(e(P'l; (xt)N))q -~ e(rl; t)
(8)
and preserving the normalization condition ~-~N=]fpq(xN; t)dFN = l, where (...)q = ~'~N~=l f . . . pq(xN; t)dFN denotes the quasi-equilibrium averaging. Solution of the conditional extremum problem 6S(p)/Sp = 0 is of the form [32]
pq(xN;t)
- i=] 2(xi; t)
= z-l(t)exp
- Zi
(9)
where 3ij = ½(fl(ri;t) + fl(rj;t)) and the Lagrange multipliers 2(x;t),fl(r;t) and Z(t) are determined from the self-consistent Eqs. (7), (8) and the normalizing condition, respectively. The effect of the boundary condition (3) is equivalent to transition from the Liouville Eq. (1) to the equation with a source in the right-hand part, which destroys the invariance of this equation with respect to the time inversion [33]
-~ + LN p(xN;t) = --e(p(xN;t) -- pq(xN;t)),
(10)
where the limit e ~ +0 is taken after the thermodynamic limit V ~ w while calculating the average values. Let us integrate Eq. (10) over dFu-s at certain N. Then, taking into account that p is equal to zero on the boundary of the phase space and performing the summation over the number of particles N, we obtain the BBGKY hierarchy with the boundary condition
( ~ + Ls) f~(xS;t)-~-~ f L(i,s + l)f,+l(X~+';t)dx~+] i=1 = -e(f~(xS; where
f~(xS; t)
t) - fq,(XS; t)),
= (fiU(x~; (x')N)) t and
(11 )
fqs(xS; t) = {hUCx , s , s", ,(x!~N~ I ,!qt
are the s-particle nonequilibrium and quasi-equilibrium distribution functions, respectively. The quasiequilibrium distribution function can be presented in a more explicit form f q ~ (xS" ,t)
=
Gs(r';t) ~ I f](xi;t) i--I
(12)
94
LP. Omelyan, M.V. Tokarchuk / Physica A 234 (1996) 89-107
through the s-particle quasi-equilibrium function in coordinate space ~s(rS; t)
Gs(rS;t)- I-Ii=l ~l(ri;t)
(13)
S
where ~s(rS;t) = z - l ( t ) Z
/drN_,exp
N=s+ 1
i=s+ 1
, ~ ( x i ; t ) - Z flij~PiJ i
(14)
is the auxiliary function. It can be shown easily that the dependence of 2(xl;t) on velocity is determined by the one-particle function only: 2(Xl;t) = In ~l(rl; t ) - l n fl(xl;t); therefore, correlations in the space of velocities are not taken into account in the quasi-equilibrium statistical operator (9). However, the absence of these correlations in the infinitely distant past (see Eq. (3)) does not make essential limitations for the system, because, for an arbitrary moment of time they can appear as a result of interactions. Account of the configurational correlations leads to automatical satisfaction of the conservation laws for hydrodynamic variables (number density of particles, momentum and total energy). This is important for the construction of kinetic equations for dense gases and liquids, where the small parameter is absent. In the "pair collisions" approximation we put f3 = 0 for the second equation (s -- 2) of the hierarchy (11 ). Obviously, this is not in the spirit of the two-particle approximation as in the Boltzmann theory, because the essential part of many-particle correlations is taken into account implicitly by the factor G2(rl,r2;t). Then the second equation becomes simpler: ~-~ + L2 + e f2(xl,x2; t) =
eG2(rl, r2; t)fl (Xl; t)fl(x2; t)
(15)
and its formal solution is of the form 0
f2 ( X2 ; t) = e /
dz exp (( e + L2 )Z) G2(r l, r2; t + z)fl (xl ; t + z)fl (x2; t + z).
--oo
(16) According to Abel's theorem, this solution can be rewritten as follows:
f2(x2;t)=
lim
z---} - - ( x 3
exp(L2z)G2(rl,r2;t + z)fl(xl;t + z)fl(x2;t + z).
(17)
Substituting (17) into the first equation (s = 1) of the hierarchy (11), we obtain a general form of the kinetic equation
(-~ + L(1)) fl(Xl;t)= I(x,;t),
(18)
I.P. Omelyan, M. K TokarchuklPhysica A 234 (1996) 89-107
95
where I ( x l ; t ) = f d x 2 L ( 1 , 2 ) lim exp(L2z)G2(rl,r2;t + z) J × f l ( x l ; t + z ) f i ( x 2 ; t + z) "~----+ - -
O~
(19)
is the collision integral.
3. Multistep potential of interaction
The collision integral (19) is still sufficiently complicated to be written in an explicit form for an arbitrary potential of interaction. As is now known [14, 18], only in two particular cases, namely, for the hard-sphere and square-well potentials, this integral has been reduced to an analytical form. In order to extend those previous results, we consider the interparticle potential in the form of a multistep function ~Pij = ¢~s = lim ~ij~0(rij),
(20)
where
eo
/ 80,
(I)ij(rij) =
f,k,
t O,
rij < ¢70, Gk-I < rij < ak,
k = 1.... , N *
(21)
ao* < rij
and N* is the total number of walls except that of the hard sphere. The potential with a hard core contains the strong singularity (operator L ( i , j ) is badly defined). That is why we shall deal with the potential #ij,0(rij) and only in final expressions we shall put e0 ~ ~ . Let us separate for convenience systems of attractive and repulsive walls. Let n* be the number of repulsive walls at the distance ari - ai and of the height Ae~ = e~ - e~+l > 0, where e~ -= ei, i = 1. . . . . n*; and m* be the number of attractive walls with the corresponding parameters aaj =- aj+,,., Ae~ = e~+1 -e~ > 0, where e~ = ej+,,., j = 1..... m*; en.+l =e~, e m a . + l = 0 , N* : n* + m*. Thus, the parameters a0 (hard r * ,O'aj and A~ completely determine the geometry of sphere diameter), n*, O'ri,Aei,m the multistep potential (see Fig. 1, where a specific case n* = 1, m* = 3 is designed and the multistep function approximates some smooth real potential). The limit r ~ -cx~ in the collision integral (19) is a formal one in the mathematical sense. On the physical level of description this means that ITI ~>~0, where z0 > 0 is some characteristic interval of time. For different values of r0 we can investigate different stages of evolution for the system. In the Boltzmann theory [34] [z] ~zc, where Zc is the time of binary interactions. On the other hand, I~1 is far less than a characteristic scale of time Zm for hydrodynamic variables. The above situation is possible because of the fact that for dilute gases kinetic (z ,-~ %) and hydrodynamic ('C '~ "~m) stages of evolution are far from one another in time, i.e., Zc ,~'cf ~Tm, where zf is the characteristic time of free motion. The pattern is qualitatively different in dense gases and liquids with a realistic smooth potential of interaction. The kinetic and r
I.P. Omelyan, M. V. Tokarchuk / Physica A 234 (1996) 89-107
96
EO
E1
E2
li\l ~0
ffl
y..--_...~..-E4 t ...........
if2
E3
G3
~4
Fig. 1. The multistep potential at n* = 1 and m* = 3. The dashed curve corresponds to a smooth real potential.
hydrodynamic stages appear to be closely connected. Besides, such quantities as the length of free motion and time of interaction are not defined in the usual way, because all particles influence the dynamics of interaction for some chosen pair of particles. However, for special types of potentials, such classification of times remains valid even for high densities. For the multistep potential the region t2 of binary interactions consists of the following subregions [ t r k - Aro, ak ÷ Aro], k = 0 , 1 , . . . , N * , where Aro ~ +0 due to the singular nature of the potential under consideration. This potential has the finite range max{ak} = trN. > 0 of action. We can introduce the following set of specific time intervals: c0 -- Aro/9o ~ +0 is the time of pair interactions on the walls, Az = min{ak - a k - l } / 9 o > 0 is the time of motion between the two nearest neighbouring walls and T = aN*~90 > 0 is the time of motion of whole systems of walls for some pair of particles, where 90 is the averaged relative velocity of two particles. Because the potential contains the horizontal parts, where the force of the interparticle interaction is equal to zero, it is possible to introduce also time .t-f as an average time of free motion of particles in the system. Changing geometry of the potential and increasing density, we can make this time arbitrarily small in order to fulfil the relation Zo~f~Az
< T<
"Cm .
(22)
The kinetic stage of the evolution corresponds to small times of order z ~ zf. Therefore, as far as the relation (22) is satisfied, the formal limit z ~ -cx~ must be considered as z/z0 ~ - ~ or merely as I~1 ~>~0. For smooth nonsingular potentials the region of interaction has some nonzero size and z0 is finite. For singular potentials, the case
L P. Omelyan, M. K Tokarchuk / Physica A 234 (1996) 89-107
97
z0 -~ + 0 is possible and the limit Z/Zo ~ - o o remains valid even for z -~ - 0 if only z0 is a value of the higher order infinitesimal
lim { lim+oZ/Zo}--~-c~.
~----~-- 0
(23)
"C(
Thus, the formal limit z --~ - 0 can be applied for the collision integral (19) in the sense of (22) and (23). As is now well established, the limit z -* - 0 leads to the kinetic equations of the RET [23,31] and RDRS [27] theories• The limit r ~ - 0 allows one to ignore the time retardation in (19). By the construction, G2 depends on time only functionally through the local density n(r; t) = f fl(x; t ) d v and "reciprocal temperature" fl(r; t), i.e., G2(rl, r2; t) = G2(rl, r2 [n(r; t), fl(r; t ) ) .
(24)
Then lim G 2 ( r l , r 2 ; t + z ) = ~-,-o
G2(rhr2;t)+
ao2
limz If ~--,-o [ J J
}=
dr [ 6(72 an(r;t) L&(r;t) at O2(r
(25)
,r2;t)
because during the time z ,-~ %f <~1:m hydrodynamic parameters n(r; t) and fl(r; t) could not essentially change. Besides,
(
[
f l ( x l ; t + z) = e -L(1)r f l ( X l ; t ) + z L(1) + ~
f[(Xl;t)
}
,
(26)
where 1
p
T
(xl; t) = - ] dz' eL(l)*'f l (xl; t) TJ 0
(27)
are the average values of the one-particle function during the time interval z. But • -3 llm,~_0 f l (Xl; t) = fl (xl; t) and
~li~m-of~(xl;t + z) = ,~-01ime - t O ) * f l ( x l ; t )
(28)
because it follows from (22) that during the essential part of time z the motion of particles is free and during the very short interval z0 they interact. Noting that fl(x; t) is a continuous function of r, we obtain
af, lim e-L~l)~ f l ( X l ; t ) = f l ( x l ; t ) -- ~lim_ --0 "cvl ~Crl = f l ( x l ; t ) .
(29)
T~+ - - 0
As a result of (25), (28) and (29) the collision integral (19) transforms into the form l(Xl;t)
f d x 2 L(1,2) lim_o et2~ Gz(rt, r2; t ) f l (xl; t ) f l (x2;t).
(30)
We consider the dynamics of interaction for two particles in the centre of mass system. Then L(1,2) = - I/p* d a ~ f f a r l 2 0 / 0 9 , L2 = -g3/art2 - L(1,2), where rj2 =
98
L P. Omelyan, M. V. Tokarchuk l Physica A 234 (1996) 89-107
i'1 -- r2, 0 ~- V2 -- !~1, ~ -~ r12/r12 and #* = #/2 is the reduced mass. Integrating (30) over dr2 in the spherical system of coordinates with the centre coinciding with the centre of the first particle yields oo
l<,r,&i
I =
2 1 ^'+12 _~_~ eL2*G2(rl, r2, t )fl (Xl; t)fl (X2; t), d g r 1 2 7 f f Orl2 09 ~lim -*-0
o (31) where d6 is the element of the surface sphere of unit radius. According to the definition (21) of the multistep potential we have ( ' 0 -- ' ] ) ~ ( r 1 2 -- 0"0) -- Z AI3~(r12 - ffri) -t- Z i=1 j=l
Or12 --
Ae~b(ri2 -
O'aj)
(32)
and, therefore, it is quite sufficient to perform the integrations in (31) over the coordinate space only inside the region 12 of interaction. Taking into account that fl is continuous in an arbitrary point of space and G2 is continuous over r12 out of the region ~, one can show that the following equalities take place:
i dr12L(1,2) lim eL2* G2(rl,r2; t)fl(Xl; t)fl(x2; t) ~-*-0
f2
= -- [ drl2 9_~._0 l i m e L2~G2(rl, r2;t)fl (Xl ;t)fl (x2; t) d Q
= -S
or12 z---+-O
d D 9 ,-~lim-oeL2~Gz(rl, r2; t)fl(Xl; t)fl(x2; t),
(33)
D
where the theorem about the gradient has been used, so that D is the surface, spanned on the volume f2. In view of Eqs. (32) and (33), we can write the collision integral for the multistep potential in the form N*
I(xl; t) = Z I k ,
(34)
k=0
where
lk=
lim z--+--O dro-*+O
a
27 S d9
dtr trgeL'~ G2(rl,r2;t)fl(xl;t)fl(x2;t)
'ri2=ak-- Ar°
(35)
rl2=trk+Ar°
and I0 = Ins is the hard-core part of the collision integral, whereas If = Ii and 17 = Ij+n. are the parts caused by interaction on the ith repulsive and jth attractive walls, respectively.
99
I.P. Omelyan, M.V. TokarchuklPhysica A 234 (1996) 89-107
Owing to the condition (22) and the evidence of the limit r ~ - 0 , the collision integral (34) breaks down into several parts. Physically, this means that during the time interval of order of free motion zf, two particles can interact only at one of the walls of the multistep potential. Interaction of same pair of particles during this time interval in several subregions f2 is impossible. Moreover, if one considers a separate pair of particles, it appears that in the velocity space their interaction at one of the walls does not depend on the interaction at the remaining ones. It is caused by the fact that during the time of motion between neighbouring walls, a pair of particles interacts many times with the medium and velocity correlations weaken. So, during some interaction of the pair of particles at a potential wall, this pair does not "remember" its velocities after the previous "meeting". A large part of the configurational correlations is conserved and is taken into account by the factor G2. The property of additivity for collision integrals will break, if condition (22) is not fulfilled. In this case, one should consider complex interference phenomena. As we can see from Eq. (35), all parts of the collision integral have the same structure. In order to obtain analytical expressions for them, it is necessary to consider the action of the two-particle operator of evolution l i m ~ _ 0 e z2~ at the external r12 tTk ~-Aro and internal r12 = a k - Aro surfaces of the sphere of radius ak. The results for Ih~ and I] are known. For the hard sphere wall the collision integral is [23] :
I h s = O -2
//
d~[~910(~9)[F2(rl,V'l,rt + a 0 a , v2;t )
dg
-F2(Xl, rl - a+ 6, v2; t)],
(36)
where v~=vl+6(69),
l = 1,2
(37)
are the velocities of a pair of particles after reflection from the hard sphere wall, F2(xl,x2; t) = G2(rl,r2; t)fl(xl; t)fl(x2; t) ~ fqz(X2; t), o-± = limzr0~+o(a 4- Aro) and O denotes the Heaviside function. In the case of a finite attractive wall we have [27]
I; =
jfdo/d +0(-69
lagl{ --
O(tTo)[F2(¥1, Vatl,rl-{-'°'aj~ vtt" .2, t ) - F 2 ( x .
r~ - ~ ,
v2; t)]
a))[F2(rl,Val,rl +a+.~, Va2 Ii,., t) -- F2(Xl, rl -- O'aj~,V2; t)] tit
+O(--~g)O(c~g + o~j a )[F2(rl, vl,rl ' + O'aj~ - ^ 102, " t ) - - F 2 ( x l , r t - - a ~ j f f , v2;t)]}
(38) where ,,
10al -~- ~l -4-
+
l = 1,2
(39)
1.P. Omelyan, M.V. Tokarchuk/Physica A 234 (1996) 89-107
100
are velocities of particles after leaving the region of interaction, whereas '" : I)ai
1~i -4-
½e(eg + v/(eg)
2
-- (~aj)'
(40)
i = 1,2
are velocities of particles after entering the region of interaction and aaj =
~/..
Considering the dynamics of two particles at a finite repulsive wall, we have obtained
17 =a2ri / dg / d~l~gl{ O ( - ~ 9 ) -[- ^
If,
>[F2(rlVtrtl,rl q- (TritT, Vr2,t ) -- F2(Xl,r
I - o-~Tt~, v2; t)]
q-O((rg -- (~[)[F2(rl,vrl '" ,rl q- O'rT.(7, vr2, m. t /-~- F2(Xl, rz - a+t~, r2; t)]
+O(t~g)O(-~g + ¢~)[F2(rl, v'1, rl + a+ ~, v'2; t) - F2(Xl, rl - tr+tr, v2; t)]}. (41) Here I~ri = 13i d:
+ ~ri),
i = 1,2
(42)
are velocities of particles after leaving the region of interaction, where 0~ri VaAs~/[,I. This process does not require additional conditions in 9 and takes place for arbitrary 89 < 0. At the same time ----
'" lJri
=-
Vi 4-
}e(eg + v/(eg)
-
-
O~ri),
i = 1,2
(43)
are velocities of particles after entering the region of interaction. This process is possible when the relative velocity of two particles is sufficient to pass the potential wall. Otherwise, when the relative velocity is too low, velocities of particles change during interaction as after reflection (see Eq. (37)). Therefore, except for the terms Ihs (RET type) and I] (RDRS type), we have obtained new terms I r (41) which are caused by contribution of repulsive walls of finite heights into kinetics. Obviously, the repulsive wall of a finite height is a generalization of the infiniteheight hard sphere wall. Taking into account (13), (14) we have lim
_ = 0.
G2(rl,r2;t)
ZIg~--+OO
(44)
r12 =o'ri
This means that particles cannot penetrate the infinite-height wall. Then it is easy to verify that the term I r (40) at Ae[ --+ oo and ari =-=-ao reduces to the hard core type (36), i.e., /hs = Ae~---+oo lim /ir -,,-o0"
(45)
The collision integral (34) can be expressed in a more compact form
1(xi; t) = f dx2 "PjzG2(rz, re; t)fl(Xl;t)fl (xz; t) J
(46)
LP. Omelyan, M.V. Tokarchuk / Physica A 234 (1996) 89-107
101
in terms of the two-particle pseudo-Liouville operator n*
m*
(47)
E ?:+E Z ?; j = l p=b,c,d
i=1 p=b,c,d
which consists of the corresponding pseudo-Liouville operators
?hs =
~_ ^
^a
o-20 d(~ [(~g[O((~g){a(r,2 4- o-o o-)B ( 6) - 6(r12 -- O'g-a)}
(48)
for the hard sphere wall, ?rf = o-ri2
d6
lp^
[(~olOP(~g){a(rt2 +
o-r;
^P ^ 2p^ ° ' ) B r i ( ° ' ) - ~(¥12 - o-r; o-)}
(49)
for the ith repulsive wall and
J dd I@lO~.(@){5(r,2 +
^ = o-aj f 7;P 2
l , p ^ ^P
o-aj o-)Baj(o-)
^
2, p ^
- a(l*12 - o-aj (7)}
(50)
for the jth attractive wall. Here, o~,( @ ) = o ( - @
o~,(@)
)
o~.(@) = o(@)o(~
=
o ( @
-
~[),
Oabj(@) = O(@),
- @),
O~j( @ ) = 0 ( - - @ ) 0 ( 4 + @ )
O~aj(@ ) = O ( - - @ - ~ ) ,
(51)
are the corresponding Heaviside functions, l,b
Gri
1,c
+
= o-ri, -
o-ri ~ o-ri , I,d
¢7ri
+
~- o-ri,
2,b o-ri = O'r/, 2,c
o-ri
2,d
Gri
+
1,b o-aj = o-aj, 1,c
+
~ o-ri,
o-aj = o-aj,
+ = o-r/,
o-aj
1,d
~- o-aj'
2,b a~ 0"aj ~ o- , 2,c -O'aj ~ o-aj, 2,d
-
o-aj ~ Gaj
(52)
and /~P are operators of velocity displacement caused by the interaction of p-type at each of walls, namely, /~a(~)~(Vl,l~2)
= ~fi(~)~(Vl,V2)
= /~dj(~)~(Vl,P2)
~bri(~)~(Vl,V2
) = aa'a" 11 11 ' " (Vrl, Vr2 )
B^b / ~ ) . . ~ ( V ~ , V 2 )
hci(~)~(Vl,V2)
~ (Vrl11!, Vr2 lit ) , = ~/~
/~cj(~)~(I)l,
= ...~(Vl,l~2) , ~ t- ,
~ (I)a,,1, I~a2 ,, ) , = "~
!~2 ) ~__ ~----------------~r. (Val ,ttt Va2ttt ) ,
(53)
where ~ is an arbitrary function of velocities. Physically, the values o-ri 1,p,to-aj 1,p,) (52) correspond to the distances between the particles just after the interaction of p-type on the ith repulsive (jth attractive) wall of the potential, whereas the values o-ri2,p.( o - a2j, )p . correspond to the distances just before the collision.
1.P. Omelyan, M. 1I.. Tokarchuk/ Physica A 234 (1996) 89-107
102
4. H-theorem
Let us show that the kinetic equation obtained satisfies a global H-theorem. In order to do this, we introduce the operator 7~+ conjugated to 7~n by the transformation
S avisdx'~(xi'xllt)'12dk(xl'x';t)=i dVlidxldC(xi.x,;t)'£~(xl.x,;O. (54) where (p and ~b are arbitrary functions of xl,x2. By straightforward calculation we obtain n*
^
i=1
m*
p=b,c,d
j=l
^p+
(55)
p=b,c,d
where
"I~ = f~glO(~g)f(r12 - a+)(Ba(~) - 1),
^P+ "~ I~alOg (#a),~ ( r12 Tri aj aj
- O ' r2,p~ i aj/
(56)
(l~Pj(~)~.~:~p -aj
~a2r~P~ • aj/
(57,
Here ~ is the operator of coordinate displacement. It acts on the relative distance between particles as follows: ~ r il(2)'P ~ ( 1r2 ) = ~ ( O'lr~2)'p I ' aj \ aj /
(58)
where q¢ is an arbitrary function of rn. Our task is to prove the H-theorem, i.e. the property dS(t) _- t)S(pq)>~0. dt St
(59)
Using (6), (9) and (14) we obtain t3S(pq) _ S(1) q- S(2) q- S(3) '
(60)
St where S (l) ----
dxl In ~l(rl;
. t3e(rl; t) S (3) = fJ d r l fl(rl;t) -~ •
it)
S (2) : -
dXl l n f l ,
St
(61)
Let us consider each term separately. Taking into account the identities (14), (54) and the kinetic equation (18) with the collision integral (46), we have for S (1) in the
I.P. Omelyan, M. V. Tokarchuk / Physica A 234 (1996) 89-107
103
explicit form
S(1)=-lidxlidx2F2(xi'xz;t){'i2vlcOi~<) 2 +8"9[~$(r12- "°v) "
+ ~
(b(rl2 - a{) - fi(rl2 - ak-))
])
(62)
.
k=l For the second term we find
S(2)= - I-2S dx,S
,/ i
= -2
dXl
dx2F2(Xl,X2;t)'F~In(fl(xl;t)fl(xz;t))
{o
dxz F2(xl,x2; t)grg
n*,m*
+Z
Z
(69)6(r12- o~) In Ba(gr)fl(xl(Xl, ;t)fl(xz;flt)fl(x2; t) t)
OPri(ag)6(ri2-aZri:') In^p fl(xl;t)fl(x2;t)
i,j=l p=b,c,d
aj
ali
B ri ( ~ r ) f l ( x 1 ; t ) f l ( x 2 ; t)
aj
(63) In order to calculate the third term it is necessary to obtain the equation of motion for the density of interaction energy. Comparing the first equation (s -- 1) of the BBGKY hierarchy (11) with the kinetic equation, multiplying the second (s = 2) equation of the hierarchy on 1/2~lz and integrating over dr1 dx2, we find, taking into account (5), (8), that
Oe(rl;t) cq---~ + ~i f
dr:
f dx2vl~-Vlr2tx:;x2;t ~,2.. )
= :'Id"Sd"O""'"'":":'S : dr: I dx2 F2(x:,x2; t)/~+CPl2 " .
2
i=1 m*
-~9(-6g)6(r12-
~r?)] + ~
As~[O(-~9 -
°(j)6(r12 - o"~.)
j=l
-- ~9(aO)a(r12 -- ~y~)] } "
(64)
104
LP. Omelyan, M.V. TokarchuklPhysica A 234 (1996) 89-107
Then the third term is of the form
S `3) ~---~
fdx2F2(x~,x2;t)
,221)10fl(rt; t)&_____~ +/~t218.gl
dXl
A~r[o(8.g - ~xr)~(rl2 - fir+) - O(-8.g)~(r12 - O-r~)]
×
L i=1 + ~
A 4 [ O ( - 8 . g - ~ ) a ( r , 2 - ~2j) -
0(8.~)6(r~2 -
~-)1
•
(65)
j=l Now we add the three terms (62), (63) and (65). Taking into account the inequality > 0 (the equality takes place only at h = z) and, since, the one-particle distribution function is positively defined, after some algebra one rewrites (60) in the form
ln(h/z) > 1 -z/h;h,z
n* m* d~(tt) ~Jhs ~- E Jri-Jr- E Jaj, i=1 j=l where
l/ /
Jhs = ~
dXl
(66)
dx2 (~(r12 - o'~){-8.g + 18.g [0(8"9)(1 - Ba(8.))}f2(xl,X2; t), (67)
Jri =
,/ dxl / dx2
"~
O'r; ) - ($(r12 - ar2[b))+
-
18.ol
o~(~a)
p=b,c,d
2p ^P ^ P / X~(rl2 - ~ri' )(1 - Bri(~r)exp(-flnEri)) f
Jaj = ~1 /
dxl
f
dx2
{-~g(5(r12
2p X ~(r12 - era) )(1 -
_
F2(Xl,X2;t ) ,
(68)
~)2 b ) _ 5(r,2 - (~'f))+
[Sg[ E oP(~g) p=b,c,d
^p . p } F2(xl,x2;t) Baj(cr)exp(-flnEaj))
(69)
,I and E~'~/= - ~ , ~ c = a~r, e g = 0, E~j = -A~), E:j = a ~ , E~ = 0. The term Jhs = 0 just as in the proof of H-theorem in the RET theory [25]. Considering the terms Jri,Jaj, one should use the fact that G2(r~,r2;t) is a discontinuous
105
I.P. Omelyan, M.V. TokarchukIPhysica A 234 (1996) 89-107
function. It is caused by the discontinuity of the initial potential ~12. It follows from (9), (12)-(14) that
/
drl2 F 2 ( x b x 2 ; t)O(r12 - O'+) r = exp(/312Ac r ) drl2 F2(Xl ,x2; t )6(rl2
- ~r+.) =
i
exp(-fllzAey)
drl2F2(Xl,X2; t ) ~ ( r l 2 - o-ri ) ,
fdr,2F2(Xl,X2; t)6(r12 - o-~). (70)
Then, one can show easily that each of the terms (68), (69) is equal to zero, i.e., Jri = 0, i = 1,n*; Jaj = 0, j = 1,m*. Then from (66) it follows that dS(t)/dt>,O and the H-theorem is proved. Let us determine the one-particle distribution function jit°)(xl;t) for which dS(t)/dt = 0. According to (66)-(69) such function must fulfil the following equations: ~(~g)6(r12 -- o-O~)[1 - - / T a ( r ~ ) ] f ( ° ) ( X l ;
t)f~°)(x2; t) = O,
OPri.(69)6 (rl2 - a2~ip) [exp (flizEPr}) -B~iy(~) 1 f(°)(xl;t)f~°)(xz;t) = O . ay
ayl
(71) By the straightforward calculation we can verify that the conditions (71) are fulfilled simultaneously at
l;°,(x,;,) _-,(,,.,)' k,('(')7 21t /#
exp
{-9(,)(o, - ,,))2}
,
where u = fvfl(x;t)dv is the hydrodynamic velocity. In the equilibrium spatially homogeneous case n,/7 = const., u = 0 and f(0) tends to the Maxwell distribution. In the presence of external field W(r) the kinetic equation transforms into the form (0
+L(1)
lOW(rl) ~ )
#
Orl
cGvi
fl(xl;t)
---- S dx2 l'12F2(Xl,X2;t),
(73)
i.e., only flow term is changed, while the collision integral remains the same (in the limit z ~ - 0 the external field does not change velocities of particles during collision). Putting fl - f(0) at n =- n(r), /3 = const., u = 0 and integrating (73) over drl, we obtain the integral equation for local density in an equilibrium state
± ~grl
[ln n(rl ) +/3 W(rl )] =
x
6(r12 - G~-) + Z ( 1
fdr2~n(,2
)Gz(rl, r2 )
-exp(-/312A~[)6(rl2
-- o-ri+)
i=1
+ Z j=l
where G 2 ( r l , r 2 )
(1 - exp(/312A~) ~(r12 -- O-+)
}
,
is the equilibrium distribution function.
(74)
106
I..P. Omelyan, M. II.. Tokarchuk/Physica A 234 (1996) 89-107
5. Conclusions We consider some special cases of choosing for the parameters of the model multistep potential (21). Choosing n* = m* = 0 corresponds to the hard-sphere potential. Here I[ =- 0 and I~ - 0, i.e., I -= Ihs and the kinetic equation is the equation of the RET theory [23]. In the case when n* = 0, m* = 1 and Ae~ ~ 0, the initial multistep potential transforms into the square-well one and I = Ihs + IF coincides with the collision integral of the RDRS theory [27]. It is interesting also to consider the case, when the initial potential transforms into a some smooth potential @t(r12) at r > ~0, i.e., LI~[ ~ O,
Ae r ~ O ,
i = 1..... n*; n*---*oo,
A~a ~ O,
Aej --~ 0,
j = 1,...,m*; m* --* c~,
(75)
AF.j--~ Cg@____L",2=~., t
(76)
at the condition -
-
-
-
--.)
O(IJt
0r12
r,2=ari
'
A¢~a
¢3r12
'
where Aa r = Aai, A ~ a = AM+,* and Aak = ak - ak-l. In such situation it can be shown that the collision integral transforms into the form I = Ihs + IMF, where IMF =
O~t ~F2(xl,x2; t)
dx2 t~rl
t~vt
(77)
coincides with the mean-field term of the KMFT theory [28]. The last case, however, does not fulfil the condition (22) for geometry of the interparticle potential ( d z becomes too low). As a result, the limit T --+ - 0 leads to the term IMF without collisions. At z ~ - 0 direct contribution of collisions into kinetics is provided only for singular potentials, when the distances between neighbouring walls are sufficiently big, and for high densities. On the other hand, besides the error due to the fact that Eq. (22) is not fulfilled for smooth potentials, there is also a model error caused by the fact that the multistep potential differs from real (as a rule smooth) potentials. As one can see, these two errors have opposite tendency. Thus, in order to use the obtained kinetic equation for real systems, it is necessary to find some compromise solution in the problem of geometry for the multistep potential. This problem can be solved as follows. Let the transport coefficients be as already obtained on the basis of the kinetic equation. Then for an optimal representation of a real potential, parameters of the multistep function can be accepted as those that give the best fit to the transport coefficients in a whole region of temperatures and densities. This will be the subject of a separate study. The condition (22) is not fulfilled when density decreases. In the limit of dilute gases the time of flee motion becomes too big: zf ~ T ~ d z and this fact does not give the possibility of the kinetic equation under consideration to transform into the Boltzmann one. Therefore, the kinetic theory proposed in this paper can be applied exclusively for dense gases and fluids.
I.P. Omelyan, M.V. Tokarchuk/Physica A 234 (1996) 89-107
107
Acknowledgements The authors warmly acknowledge stimulating discussions with Prof. D.N. Zubarev and Prof. V.G. Morozov concerning this work. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34]
N.N. Bogoliubov, in: Studies in Statistical Mechanics I (North-Holland, Amsterdam, 1962). J.R. Dorfman and E.G.D. Cohen, J. Math. Phys. 8 (1967) 282. J.G. Kirkwood, J. Chem. Phys. 14 (1946) 180. J.G. Kirkwood, J. Chem. Phys. 15 (1947) 72. M. Born and H.S. Green, Proc. Roy. Soc. A 188 (1946) 10. E.G.D. Cohen, Physica 28 (1962) 1045. E.G.D. Cohen, J. Math. Phys. 4 (1963) 183. L.S. Garcia-Colin and A. Flores, J. Math. Phys. 7 (1966) 254. G. Sandri, Ann. Phys. (N.Y.) 24 (1963) 332, 380. E.A. Frieman, J. Math. Phys. 4 (1963) 410. K. Kawasaki and J. Oppenheim, Phys. Rev. 139 (1965) 1763. J.R. Dorfman, in: Lectures in Theoretical Physics (Gordon and Breach, New York, 1967). E.G.D. Cohen, Physica A 194 (1993) 229. S. Chapman and T.G. Cowling, The Mathematical Theory of Non-Uniform Gases (Cambridge University Press, Cambridge, 1952). H.G.M. Hanley, R.D. McCarty and E.G. Cohen, Physica 60 (1972) 322. J.A. Barker and D. Henderson, J. Chem. Phys. 47 (1967) 4714. J.D. Weeks, D. Chandler and H.C. Anderson, J. Chem. Phys. 54 (1971) 5237. H.T. Davis, S.A. Rice and J.V. Sengers, J. Chem. Phys. 35 (1961) 2210. R.C. Castillo, E. Martina, M. Lopez, J. Karkheck and G. Stell, Phys. Rev. A 39 (1989) 3106. J. Karkheck and G. Stell, J. Chem. Phys. 75 (1981) 1475. R.C. Castillo and J.V. Orozco, Int. J. Thermophys. 11 (1990) 1025. M. Grmela, R. Rosen and L.S. Garcia-Colin, J. Chem. Phys. 75 (1981) 5474. J. Karkheck and G. Stell, Phys. Rev. A 25 (1982) 3302. H. Van Beijeren and M.H. Ernst, Physica 68 (1973) 437. P. Resibois, J. Star. Phys. 19 (1978) 593. P. Resibois, Phys. Rev. Lett. 40 (1978) 1409. J. Karkheck, H. Van Beijeren, I. de Schepper and G. Stell, Phys. Rev. A 32 (1985) 2517; see also H. Van Beijeren, J. Karkheck and J.V. Sengers, Phys. Rev. A 37 (1988) 2247. G. Stell, J. Karkheck and H. Van Beijeren, J. Chem. Phys. 79 (1983) 3166. D.N. Zubarev and V.G. Morozov, Teor, Mat. Fiz. 60 (1984) 270 [in Russian]. D.N. Zubarev, Nonequilibrium Statistical Thermodyamics (Nauka, Moscow, 1971) [English transl.: (Plenum Press, New York, 1974)]. D.N. Zubarev, V.G. Morozov, I.P. Omelyan and M.V. Tokarchuk, Teor. Matem. Fiz. 87 (1991) 113 [in Russian]. D.N. Zubarev, in: Itogi Nauki i Tekhniki, Sovr. Prob. Mat., Vol. 15 (VINITI, Moscow, 1980) [in Russian]. D.N. Zubarev and M.Y. Novikov, Teor. Mat. Fiz. 13 (1972) 406 [in Russian]. L. Boltzmann, in: Lectures on Gas Theory (University of California, Berkeley, 1964).