Energy fluctuations in a randomly driven granular fluid
Abstract
We study the behavior of the energy fluctuations in the stationary state of a uniformly heated granular gas. The equation for the one-time two-particle correlation function is derived and the hydrodynamic eigenvalues are identified. Explicit predictions are subsequently determined for energy fluctuations. The results explain Monte Carlo numerical data reported in previous studies [Eur. Phys. J. B 51, 377 (2006)].
Full text
Molecular Physics Vol. 107, Nos. 4–6, 20 February–20 March 2009, 383–392 INVITED ARTICLE Energy fluctuations in a randomly driven granular fluid Marı´a Isabel Garcı´a de Soria a *, Pablo Maynar bc and Emmanuel Trizac a a Universite ´Paris-Sud, LPTMS, UMR 8626, Orsay Cedex F-91405 and CNRS, Orsay, F-91405, France; b Laboratoire de Physique The ´orique (CNRS UMR 8627), Ba ˆtiment 210, Universite ´Paris-Sud, 91405 Orsay cedex, France; c Fı´sica Teo ´rica, Universidad de Sevilla, Apartado de Correos 1065, Sevilla, E-41080, Spain (Received 7 November 2008; final version received 2 February 2009) We study the behavior of the energy fluctuations in the stationary state of a uniformly heated granular gas. The equation for the one-time two-particle correlation function is derived and the hydrodynamic eigenvalues are identified. Explicit predictions are subsequently determined for energy fluctuations. The results explain Monte Carlo numerical data reported in previous studies [Eur. Phys. J. B 51, 377 (2006)]. Keywords: energy fluctuations; granular gas 1. Introduction Recent years have witnessed ongoing interest in the microscopic and macroscopic properties of granular media. In such systems, a simple ingredient—energy dissipation resulting from collisions—has far-reaching consequences [1], with rich phenomenology: non- Gaussian velocity distributions [2,3], non-equipartition of energy [4–7], and spontaneous symmetry breaking [8–10], to name but a few. Theoretically, one of the tools used to understand this body of phenomena is kinetic theory, which is extended naturally to these systems by introducing an ‘inelastic collision rule’ in which the energy is not conserved. Although most of the work carried out thus far has focused on the oneparticle properties, and on the study of the corresponding Boltzmann equation, it has been shown that correlations are also important and, as a matter of fact, necessary to understand the behavior of the system when vortices or cluster are developed [11,12], or even in simpler situations where spatial homogeneity is enforced [13,14]. As a consequence of dissipation in collisions, the total energy of an isolated granular system decays monotonically in time. Under certain conditions, the system reaches a homogeneous cooling state in which the time dependence of the one-particle distribution function is entirely embodied in the kinetic (so-called granular) temperature, which evolves with time as t 2 [2,15]. It is experimentally difficult to probe such a regime (see, however, [16,17]). Nevertheless, it is possible to maintain a granular system in the fast-flow regime by injecting energy in such a way that a stationary state is reached. In these states, the energy injected by the thermostat is compensated by the energy dissipated in collisions. Several mechanisms can be introduced in order to obtain a stationary state. If, for example, energy is injected by a moving boundary such as a vibrating piston, the system reaches an inhomogeneous stationary state [18]. In this work, we will focus on a granular gas that is heated uniformly, coupling the velocity of each particle to a white noise, the so-called ‘stochastic thermostat’ [3,19–27]. For this kind of forcing, which is relevant for some two-dimensional experimental configurations with a rough vibrating piston [28], the system reaches a homogeneous stationary state after a transient regime. The advantage of such a driving mechanism is that it lends itself to theoretical progress. In this context, the single-particle distribution function has been characterized [3] and long-range correlations predicted from a hydrodynamic treatment [21]. More recently, the fluctuations of the total energy have been analysed [29] (see also [30] for a related numerical study in an inhomogeneous system). In Ref. [29], the second moment of the total energy fluctuations was evaluated by neglecting the correlations, which, by and large, could not explain the simulation results. The objective of the present work is to clarify and quantify the influence of the inelasticity induced correlations on the total energy fluctuations. The methods used bear some similarities to those reported in [31], where it was shown that, for *Corresponding author. Email: [email protected] ISSN 0026–8976 print/ISSN 1362–3028 online ß2009 Taylor & Francis DOI: 10.1080/00268970902794842 http://www.informaworld.com
the unforced system, the contribution coming from the correlations is of the same order as that coming from the one-particle distribution function itself. The paper is organized as follows. In Section 2, the equation for the two-particle distribution function is derived, taking due account of the thermostat, while in Section 3 the results are particularized to the homogeneous stationary state, which will play the role of our reference state in subsequent analysis. There we also summarize the main results already known and pertaining to the one-particle distribution function. In Section 4 we analyse the hydrodynamic equations for a homogeneous linear perturbation of the reference state, and obtain the corresponding modes and eigenvalues, which are finally used in Section 5 to meet our objective and obtain an explicit expression for the variance of the total energy. 2. Heated granular gas: two-body kinetic description We consider a gas of Nhard disks (dimension d¼2) or spheres (d¼3) of mass mand diameter that collide inelastically with a coefficient of normal restitution [1]. The system is heated uniformly by adding a random component to the velocity of each particle at equal times [21,22]. The driving is implemented in such a way that the time between random kicks is small compared with the mean free time. Then, between collisions, the velocities of the particles undergo a large number of kicks due to the thermostat. In addition, we will assume that the ‘jump moments’ of the velocities of the particles verify Bij, lim t!0hvi,vj,i t ¼2 0ij þ2 0 Nðij 1Þ ,i,j¼1, ...,N, ,¼1, ...,d,ð1Þ where we have introduced v i, v i, (tþt)v i, (t), v i, (t) being the component of the velocity of particle iat time t. We have also introduced the strength of the noise, 2 0, and hi, which denotes the average over different realizations of the noise. The non-diagonal terms (corresponding to i6¼jand ¼) are necessary in order to conserve the total momentum. In the dilute limit, assuming molecular chaos, i.e. that no correlations exist between colliding particles, and that the sizes of the jumps due to the thermostat are small compared with the velocity scale on which the distribution varies, the equation for the single particle distribution function in our system is the Boltzmann– Fokker–Planck equation [3,22] @ @tfðx1,tÞþLð0Þðx1Þfðx1,tÞ¼J½fjfþ2 0 2 @ @v1 2 fðx1,tÞ, ð2Þ where x i is short-hand for position-momenta coordinates {r i ,v i } and Lð0Þðx1Þ¼v1@ @r1 :ð3Þ The inelastic collision operator J[fjf] reads J½fjf¼d1Zdv2 T0ðv1,v2Þf1ðr,v1,tÞf1ðr,v2,tÞ,ð4Þ where T0ðv1,v2Þ¼Zd^ rð^ rgÞð^ rgÞ½2b1 ^ 1,ð5Þ with g¼v 1 v 2 the relative velocity, the Heaviside step function, ^ ra unit vector joining the centers of the particles at contact and b1 ^ ran operator replacing the velocities v 1 and v 2 appearing on the right by the precollisional values v 1b1 ^ rv1¼v11þ 2ðg^ rÞ^ r,ð6Þ v 2b1 ^ rv2¼v2þ1þ 2ðg^ rÞ^ r:ð7Þ The term 2 0 2 @ @v1 2 fðx1,tÞ is a diffusive Fokker–Plank term, and is a signature of the external noise. As we shall be studying fluctuations, it is convenient to introduce the two-particle distribution function, f 2 (x 1 ,x 2 ,t). The quantity f 2 (x 1 ,x 2 ,t)dx 1 dx 2 is defined as the number of pairs of particles in which one lies inside the differential volume dx 1 centred in x 1 and, likewise with dx 2 ,x 2 for the second particle. This definition is easily generalized to higher n-particle distribution functions, f n (x 1 ,...,x n ). The evolution equation for f 2 is [31,32] @ @tþLð0Þðx1ÞþLð0Þðx2Þ f2ðx1,x2,tÞ ¼ðr12Þd1 T0ðv1,v2Þf2ðx1,x2,tÞ þd1Zdx3½ðr13Þ T0ðv1,v3Þþðr23Þ T0ðv2,v3Þ f3ðx1,x2,x3,tÞþFTH,ð8Þ 384 M.I. Garcı´a de Soria et al.
where we have introduced F TH which accounts for the external driving. The evolution Equation (8) essentially contains three parts: the free streaming on the left-hand side, the two terms on the right-hand side corresponding to collisions, and the last term, F TH , due to the thermostat. The collisional contribution is split into one part corresponding to collisions of particles with velocities v 1 and v 2 , and the other that involves collisions of particles with velocities v 1 or v 2 with a third particle with arbitrary velocity, v 3 . The collisional contribution is identical to that appearing in the absence of forcing [31]. We concentrate now on the new term, F TH . Assuming that the sizes of the jumps due to the thermostat are small compared with the scale in which the distribution f 2 varies, we can expand F TH in the spirit of the Fokker–Planck description [33] FTH ’1 2X d ,¼1X 2 i,j¼1 Bij, @ @vi, @ @vj, f2ðx1,x2,tÞ ¼1 22 0 @2 @v2 1þ@2 @v2 22 N @ @v1@ @v2 f2ðx1,x2,tÞ,ð9Þ where we have taken into account Equation (1), and we have explicitly assumed that the jump moments B ij, do not depend on the magnitude of the velocities of the particles. Let us introduce the two-particle and three-particle correlation functions through the usual cluster expansion f2ðx1,x2,tÞ¼f1ðx1,tÞf1ðx2,tÞþg2ðx1,x2,tÞ,ð10Þ and f3ðx1,x2,x3,tÞ¼f1ðx1,tÞf1ðx2,tÞf1ðx3,tÞ þg2ðx1,x2,tÞf1ðx3,tÞþg2ðx1,x3,tÞf1ðx2,tÞ þg2ðx2,x3,tÞf1ðx1,tÞþg3ðx1,x2,x3,tÞ: ð11Þ The equation for the correlation function g 2 (x 1 ,x 2 ,t) can be obtained following the same lines as in Refs. [31,34]. Neglecting the three-body correlations, g 3 ,in Equation (8), we obtain @ @tþLð0Þðx1ÞþLð0Þðx2Þ g2ðx1,x2,tÞ ¼ðr12Þd1 T0ðv1,v2Þf1ðx1,tÞf1ðx2,tÞþ½Kðx1,tÞ þKðx2,tÞg2ðx1,x2,tÞ2 0 N @ @v1@ @v2 f1ðx1,tÞf1ðx2,tÞ, ð12Þ where we have introduced the linear operator K(x i ,t) defined as Kðxi,tÞ¼d1Zdx3ðri3Þ T0ðvi,v3Þð1þPi3Þf1ðx3,tÞ þ2 0 2 @ @vi 2 ,ð13Þ and where the permutation operator P ab interchanges the labels of particles aand bin the quantities on which it acts. As will become clear below, the 1/Nterm in Equation (12) is crucial for the calculation of the energy fluctuations. 3. The stationary state It has been shown numerically that, after a transient time, the system reaches a homogeneous stationary state [21] in which the energy input from the thermostat is compensated by the energy lost in collisions. In this section we will particularize the equations of the previous section to this state, summarizing the results that are already known concerning the one-particle distribution function and that are required for our theoretical analysis. The Boltzmann–Fokker–Planck Equation (2) for the distribution function, f H (v 1 ), in the stationary homogeneous state is 2 0 2 @ @v1 2 fHðv1ÞþJ½fHjfH¼0:ð14Þ It is convenient to introduce the scaled distribution function H fHðvÞ¼nH vd H HðcÞ,ð15Þ where n H is the homogeneous density, v H ¼(2T H /m) 1/2 is the thermal velocity defined from the granular temperature TH¼2 dnHZdv1 2mv2fHðvÞ,ð16Þ and c¼v/v H is the rescaled velocity. The distribution function has been studied in Ref. [3], where an approximate expression for H (c) was derived to second order in Sonine polynomials [35] HðcÞ¼ec2 d=21þa2ðÞS2 d=21ðc2Þ ,ð17Þ with S2 d=21ðc2Þ¼1 2c41 2ðdþ2Þc2þ1 8dðdþ2Þ,ð18Þ Molecular Physics 385
and a 2 () a coefficient related to the kurtosis of the function H (c) d dþ2hc4iH hc2i2 H¼1þa2ðÞ:ð19Þ An approximate expression for a 2 reads (see [23,36] for a discussion of various possible approximations) a2¼16ð1Þð122Þ 73 þ56d24d105þ30ð1Þ2:ð20Þ The expression for the temperature in the first Sonine approximation is TH¼md2 0ffiffiffi p ð12ÞdnHd1 2=3 ð1þOða2ÞÞ,ð21Þ where d ¼2 d/2 /(d/2) is the d-dimensional solid angle. We now turn to the equation for the correlation function, g 2,H (x 1 ,x 2 ). It is convenient to introduce the rescaled correlation function ~ gHvia g2, Hðx1,x2Þ¼ nH ‘dv2d H ~ gHðl12,c1,c2Þ,ð22Þ where ‘¼(n H d1 ) 1 is proportional to the mean free path and l¼r/‘. In these units, the equation for the reduced function ~ gHreads ðc1Þþðc2Þc12 @ @l12 ~ gHðl12,c1,c2Þ ¼ðl12Þ T0ðc1,c2ÞHðc1ÞHðc2Þ þe 2 0 nH‘d N @ @c1@ @c2 Hðc1ÞHðc2Þ,ð23Þ where we have introduced the linearized Boltzmann– Fokker–Planck operator (c) ðciÞhðciÞ¼Zdc3 T0ðci,c3Þð1þPi3ÞHðc3ÞhðciÞ þe 2 0 2 @ @ci 2 hðciÞ,ð24Þ with rescaled noise amplitude e 2 0¼2 0‘ v3 H :ð25Þ As can be seen from Equation (23), the correlation function, ~ gH, is determined by the properties of the linearized Boltzmann–Fokker–Planck operator, , and by the one-particle distribution function H . Consequently, it is important to study the spectral properties of , in particular the upper (hydrodynamic) part of the spectrum, in order to understand the fluctuations of global quantities. In the case of a granular gas in the homogeneous cooling state [37], and for a system under ballistic annihilation dynamics [38], it has been shown that it is possible to find the hydrodynamic eigenvalues and eigenfunctions of the linearized Boltzmann operator. Once these quantities are known, it becomes possible to evaluate the fluctuations of the relevant global quantities in the so-called ‘hydrodynamic approximation’ [31,34]. Below, we will see that we can evaluate the fluctuations of the total energy in an equivalent approximation, but without knowledge of the eigenfunction associated with the energy. The only information needed is the form of the linearized hydrodynamic equations and, in particular, the eigenvalues. 4. Hydrodynamic equations 4.1. Evolution of homogeneous perturbations The objective in this section is to consider the linearized hydrodynamic equations around a homogeneous perturbation in order to extract information about the behavior of a small perturbation of the total energy. The complete nonlinear hydrodynamic equations for a granular system heated by a stochastic thermostat are [21,25] @ @tn¼rðnuÞ,ð26Þ @ @tu¼uru1 mn rjPij,ð27Þ @ @tT¼urT2 dn ðrqþPijrjuiÞTþm2 0, ð28Þ where P ij is the pressure tensor, qis the heat flux and is the cooling rate, which is also a functional of the distribution function ¼ð12Þmðd1Þ=2d1 4d½ðdþ3Þ=2nkBT Zdv1Zdv2jv1v2j3fðr,v1,tÞfðr,v2,tÞ:ð29Þ Considering a homogeneous state, the previous equations reduce to @ @tn¼0, @ @tu¼0,@ @tT¼Tþm2 0:ð30Þ In the long time limit, the system is expected to approach a steady state with a constant temperature given by the equation HðfHÞTH¼m2 0:ð31Þ 386 M.I. Garcı´a de Soria et al.
Substituting the explicit form of the one-particle distribution function (15) in the equation above, we obtain the temperature given in Equation (21). Now let us consider a homogeneous state close to this homogeneous stationary state. We can write the hydrodynamic fields as n(t)¼n H þn,u(t)¼uand T(t)¼T H þT. We also define the dimensionless hydrodynamic fields ðÞ¼n nH ,wðÞ¼u vH ,ðÞ¼T TH ,ð32Þ where we have introduced the dimensionless time scale , proportional to the number of collisions per particle, defined as ¼Zt 0 dt0vH ‘¼vH ‘t:ð33Þ Assuming that the deviations are small, and taking into account Equations (30) and (31), we can write the linearized evolution equations for the dimensionless hydrodynamic fields in this new time scale: @ @ ¼0, @ @ w¼0, @ @ ¼0 3 20, ð34Þ where 0 ¼‘ H /v H is a dimensionless coefficient that is a functional of the one-particle distribution function in the stationary state. Its expression in the first Sonine approximation is [3] 0¼ð16 þ3a2Þðd1Þ=2ð12Þ 8ffiffiffi2 pdðd=2Þ:ð35Þ To obtain the equation for we have assumed that the perturbed distribution function scales as fðv,tÞ¼ n vðtÞdH v vðtÞ ,ð36Þ where vðtÞ¼½2TðtÞ=m1=2, and H is the same scaled distribution function as for the reference stationary state. This assumption has already been used and tested numerically in [21]. Then, the cooling rate for the state under scrutiny is proportional to T 1/2 (t) and we obtain the equation for the linearized energy written above. Equations (34) indicate that a perturbation in the total number of particles or total momentum does not decay, as a consequence of the fact that these variables are conserved, but a perturbation in the total energy will decay (exponentially in ) to the stationary value, as expected. Moreover, as the equation for the temperature can be rewritten in the following form: @ @ 2 3 þ ¼3 20 2 3 þ ,ð37Þ we can identify the hydrodynamic eigenvalues ¼0 and ¼3 20,being (dþ1)-fold degenerate. For the sake of clarity, it proves convenient to relabel these eigenvalues as 1¼0, 2¼0, 3¼3 20,ð38Þ where 2 is d-fold degenerate. The associated hydrodynamic modes, fygdþ2 ¼1, are y1¼,y2¼w,y3¼2 3 þ: ð39Þ The þsign in the last equation stems from the fact that an increased density leads to enhanced dissipation, and hence a lower temperature. 4.2. Enforcing consistency with the linearized Boltzmann–Fokker–Planck equation description We now turn our attention to the problem of finding the linearized hydrodynamic equations for a homogeneous perturbation, directly from the Boltzmann–Fokker–Planck equation. Enforcing consistency with the macroscopic considerations of Section 4.1, we will infer useful properties on the hydrodynamic part of the spectrum of (c). We first introduce the scaled deviation of the distribution function ðc,Þ¼vd H nH½fðv,tÞfHðvÞ:ð40Þ The evolution of the scaled distribution is governed by @ @ ðc,Þ¼ðcÞðc,Þ,ð41Þ where the operator (c) is the linearized Boltzmann– Fokker–Planck operator defined in (24). Let us also introduce the scalar product hfðcÞjgðcÞi Zdc1 HðcÞfðcÞgðcÞ,ð42Þ where fis the complex conjugate of f. Interestingly, the hydrodynamic modes introduced in (39) can then be written as y¼h ji,¼1, 2, 3, ð43Þ Molecular Physics 387
where 1ðcÞ¼HðcÞ, 2ðcÞ¼cHðcÞ, 3ðcÞ¼ c2 d1 6 HðcÞ: ð44Þ Taking the scalar product of the linearized Boltzmann– Fokker–Planck Equation (41) with the functions ,we obtain the linear Equations (34) (in the hydrodynamic time scale, that is if we wait long enough so that fast modes have vanished) only if the spectrum of admits the three eigenvalues written in (38), and the associated ‘hydrodynamic’ eigenfunctions, fgdþ2 ¼1, obey the orthogonality condition h 1j2i¼12,1,2¼1, 2, 3:ð45Þ In Appendix A, it is shown that the null eigenvalue is (dþ1)-fold degenerate, and the corresponding eigenfunctions, 1 and m 2 , are determined. Moreover, as a consequence of particle and total momentum conservation in a collision, 1and n2are the corresponding left eigenfunctions. We were not able to demonstrate that 3 is an eigenvalue of , but we have shown explicitly that h 3ji¼0, for ¼1, 2:ð46Þ In the following, we will assume that actually admits this third eigenvalue, with an unknown eigenfunction 3 . With the help of this assumption, we will see in the next section that it is possible to define a projector in the hydrodynamic subspace, which opens the way for evaluating the variance of the total energy fluctuations. 5. Energy fluctuations In this section, we study the fluctuations of the global energy for a system in the stationary state. As we are interested in global quantities, it is convenient to define a global correlation function H Hðc1,c2ÞZdr12 ~ gHðr12,c1,c2Þ:ð47Þ The energy fluctuations can be written as a functional of this correlation function and the one-particle distribution function, H , as [31,38] hðEÞ2iH¼m2 4Nv4 H"Zdcc4HðcÞ þZdc1Zdc2c2 1c2 2Hðc1,c2Þ#:ð48Þ In order to evaluate the integral over the correlation function H , we start from (23), integrating over the position variable. Assuming periodic boundary conditions, the spacial gradient terms disappear and we have the following equation for H : ½ðc1Þþðc2ÞHðc1,c2Þ¼ðc1,c2Þ,ð49Þ where ðc1,c2Þ¼ T0ðc1,c2ÞHðc1ÞHðc2Þ þe 2 0 @ @c1@ @c2 Hðc1ÞHðc2Þ:ð50Þ The solubility condition for Equation (49) is that does not have components in the subspace associated with the null eigenvalue. In our case, this subspace is generated by f 1, n2g. Due to the conservation of the number of particles and total momentum in a collision, and to the symmetry of the second term of , we have h 1ðc1Þ n2ðc2Þjðc1,c2Þi¼h n2ðc1Þ 1ðc2Þjðc1,c2Þi¼0, ð51Þ h 1ðc1Þ 1ðc2Þjðc1,c2Þi¼h 2,iðc1Þ 2,jðc2Þjðc1,c2Þi¼0, ð52Þ for i6¼j. The case i¼jis analysed in Appendix B, where it is shown that h 2, iðc1Þ 2, iðc2Þjðc1,c2Þi ¼ 0. In order to prove this property, the presence of the second term on the right-hand side of Equation (50) is essential. Hence, the solubility condition holds and the problem of finding H with Equation (49) is well defined. Let us also define a projector P 12 in the hydrodynamic subspace as P12hðc1,c2Þ¼X 3 1¼1X 3 2¼1h 1ðc1Þ 2ðc2Þjhðc1,c2Þi 1ðc1Þ2ðc2Þ,ð53Þ where fg3 ¼1are the right hydrodynamic eigenfunctions of the linearized Boltzmann–Fokker–Planck operator, f g3 ¼1the orthogonal set introduced in the previous section, Equation (44), and we have generalized the scalar product by hfðc1,c2Þjgðc1,c2Þi ¼ Zdc1Zdc21 Hðc1Þ1 Hðc2Þ fðc1,c2Þgðc1,c2Þ:ð54Þ Note that P 12 is a projector even if the function 3is not the true left eigenfunction of (c) (remember that 1and n2are the actual left eigenfunctions associated with the null eigenvalue). The fact that the set of functions f g3 ¼1and fg3 ¼1fulfil the orthogonality condition (45) is enough to guarantee that P2 12 ¼P12. Using this projector, we define the ‘hydrodynamic 388 M.I. Garcı´a de Soria et al.
part’ of H to be the function ðhÞ Hðc1,c2ÞP12Hðc1,c2Þ¼X 3 1¼1X 3 2¼1 a121ðc1Þ2ðc2Þ: ð55Þ The coefficients a 1 2 are the quantities we need to evaluate. As they are essentially the first moments of the correlation function H , they are directly related to the integral we have to calculate in (48). It is tempting to treat 3 as if it was an actual left eigenfunction of the linearized Boltzmann–Fokker–Planck operator, and we will in the following use the approximation P12ðciÞ¼P12ðciÞP12,ð56Þ which allows us to find a closed equation for ðhÞ H. This approximation has already been invoked in other systems such as the freely evolving granular gas [31], or the probabilistic ballistic annihilation model [34]. With the information available on the linearized Boltzmann– Fokker–Planck operator, it would seem to be the best that can be done technically. Let us also remark that the approximation is exact in the elastic limit. In the inelastic case, some exact results for the fluctuations of the total energy have been obtained for one dimensional Maxwell molecules [40]. Then, applying the projector P 12 to Equation (49) and taking into account the approximation (56), we obtain the following expressions for the coefficients a 1 2 : a12¼h 1ðc1Þ 2ðc2Þj T0ðc1,c2ÞHðc1ÞHðc2Þi 1þ2 , ð57Þ where it has been assumed that 1 þ 2 6¼0. The coefficients associated with the vanishing eigenvalue cannot be calculated using Equation (49), but are fixed by the boundary conditions. The coefficients a 1 2 are evaluated in Appendix C. The expression for ðhÞ His finally given by ðhÞ H¼a111ðc1Þ1ðc2Þþa13½1ðc1Þ3ðc2Þþ3ðc1Þ1ðc2Þ þa333ðc1Þ3ðc2Þ,ð58Þ where a 11 ¼1, a13 ¼1 3and a 33 can be obtained as a functional of the one-particle distribution function a33 ¼h 3ðc1Þ 3ðc2Þj T0ðc1,c2ÞHðc1ÞHðc2Þi 30 :ð59Þ An approximate expression is derived in Appendix C, and reads a33 ¼(15 þ7dþ14d23ð9þdð9þ2dÞÞ þ30ð1þdÞ26ð9þdÞ3Þ) 9dð19 þ2dð7þ3Þþ3ð9þ2ð1þÞÞÞ: ð60Þ Taking into account (48) and (58), the variance of the energy fluctuations can finally be written as 2 E¼NhðEÞ2iH hEi2 H¼ða2þ1Þdðdþ2Þ 4þd2a33 d25 36 : ð61Þ In Ref. [29], the value of 2 Ehas been measured by means of the Direct Monte Carlo simulation method (DSMC). In Figure 1 we compare our theoretical prediction (solid line given by (61)) with the DMSC simulations results (symbols). The agreement is satisfactory for the whole range of inelasticities, at variance with the theoretical attempt put forward in [29], which neglected velocity correlations. In particular, we note the non-trivial result in the elastic limit lim!1-2 EðÞ¼d=3 (i.e. 2/3 in Figure 1, which is well obeyed), while in the free cooling regime, this quantity vanishes [31]. We emphasize that the elastic limit is singular: the behaviour for elastic systems with ¼1 is not approached by taking the quasi-elastic limit !1 (we note that the divergence of the different moments of the velocity distribution as !1 is nevertheless indicative of the absence of a stationary state when ¼1). Such a singularity has already been reported in one dimension [39], but, to our knowledge, not for two-dimensional granular systems. It is also interesting to note that the singular nature of the quasielastic limit appears at the two-body level through the energy fluctuations, while, as far as rescaled distribution functions are concerned, the one-body level of description is regular, with a well-behaved velocity distribution approaching a Gaussian form [22]. 0 0.2 0.4 0.6 0.8 1 α 0.6 0.8 1 1.2 N=100 N=1000 σE 2 Figure 1. Scaled second moment of the energy fluctuations 2 Eas a function of the restitution coefficient . The solid line is the theoretical prediction and symbols are the two-dimensional Monte Carlo simulation results of Ref. [29]. Molecular Physics 389
6. Conclusions The problem of the fluctuations of the total energy of a granular (inelastic) gas maintained in a nonequilibrium stationary state by a random acceleration has been addressed. A numerical study of this quantity has been performed by means of Monte Carlo simulations and an argument assuming an uncorrelated non-Gaussian individual distribution function had been proposed in [29], without success. The main goal of this work was therefore to take due account of velocity correlations in order to study these fluctuations. To this end, the standard description at the single-particle level is not sufficient, and the twoparticle correlation function is needed. We have derived the evolution equation for such an object, and particularized the analysis to the homogeneous stationary state that is reached by the system in the long time limit. Our work shows that this equation is not a straightforward generalization of its counterpart arising in the context of the undriven granular gas (i.e. by only changing the linearized Boltzmann– Fokker–Planck operator into its driven form). A non-trivial non-diagonal term appears in the Fokker–Planck equation for the two-particle distribution function (contribution @ v 1 @ v 2 in Equation (9)), as a consequence of the coupling between velocities due to momentum conservation. We have seen that, for our purposes, exact knowledge of the hydrodynamic eigenfunctions is not needed. The important point is that we can construct a set of functions f g3 ¼1, which are linear combinations of 1, cand c 2 , that are orthogonal to the right eigenfunctions fg3 ¼1of the linearized Boltzmann–Fokker–Planck operator . This orthogonality property holds for the ‘real’ left eigenfunctions, that in our case correspond to the null eigenvalue (i.e. density and velocity fields associated with conserved quantities). The function 3is not a left eigenfunction of but it can be proved to be orthogonal to 1 and m 2 . In a subsequent step, the linear hydrodynamic equations around the reference state are derived and from that knowledge, the hydrodynamic eigenvalues are identified and the variance of energy fluctuations subsequently derived. Finally, our prediction has been successfully tested against the numerical results obtained by the Direct Monte Carlo simulation method for the complete range of values of the coefficient of normal restitution . This provides strong support for the theory developed here and assesses in retrospect the validity of our assumptions. Acknowledgements It is a pleasure to dedicate this work to Jean-Jacques Weis. We thank Paolo Visco for useful discussions, and for providing us with the Monte Carlo data of Figure 1. We would like to thank the Agence Nationale de la Recherche for financial support (grant ANR-05-JCJC-44482). M.I.G.S. and P.M. acknowledge financial support from Becas de la Fundacio ´n La Caixa y el Gobierno France ´s. E.T. acknowledges the support of Institut Universitaire de France. References [1] A. Barrat, E. Trizac, and M.H. Ernst, J. Phys.: Condens. Matter 17, S2429 (2005). [2] A. Goldshtein and M. Shapiro, J. Fluid Mech. 282,75 (1995). [3] T.P.C. van Noije and M.H. Ernst, Granular Matter 1, 57 (1998). [4] R.D. Wildman and D.J. Parker, Phys. Rev. Lett. 88, 064301 (2002). [5] K. Feitosa and N. Menon, Phys. Rev. Lett. 88, 198301 (2002). [6] A. Barrat and E. Trizac, Granular Matter 4, 57 (2002). [7] A. Santos and J.W. Dufty, Phys. Rev. Lett. 97, 058001 (2006). [8] H.J. Schlichtting and V. Nordmeier, MNU 49, 323 (1996). [9] J.J. Brey, F. Moreno, R. Garcı´a-Rojo, and M.J. Ruiz- Montero, Phys. Rev. E 65, 011305 (2002). [10] A. Barrat and E. Trizac, Mol. Phys. 101, 1713 (2003). [11] T.P.C. van Noije, M.H. Ernst, R. Brito, and J.A.G. Orza, Phys. Rev. Lett. 79, 411 (1997). [12] J.J. Brey, F. Moreno, and M.J. Ruiz-Montero, Phys. Fluids 10, 2965 (1998). [13] J.J. Brey, M.I. Garcı´a de Soria, P. Maynar, and M.J. Ruiz-Montero, Phys. Rev. Lett. 94, 098001 (2005). [14] J.J. Brey, A. Domı´nguez, M.I. Garcı´a de Soria, and P. Maynar, Phys. Rev. Lett. 96, 158002 (2006). [15] P.K. Haff, J. Fluid Mech. 134, 401 (1983). [16] B. Painter, M. Dutt, and R. Behringer, Physica D 175, 43 (2003). [17] C.C. Maaß, N. Isert, G. Maret, and C.M. Aegerter, Phys. Rev. Lett. 100, 248001 (2008). [18] J.J. Brey, M.J. Ruiz-Montero, and F. Moreno, Phys. Rev. E 62, 5339 (2000). [19] D.R.M. Williams and F.C. MacKintosh, Phys. Rev. E 54, R9 (1996). [20] A. Puglisi, V. Loreto, U.M.B. Marconi, and A. Vulpiani, Phys. Rev. E 59, 5582 (1999). [21] T.P.C. van Noije, M.H. Ernst, E. Trizac, and I. Pagonabarraga, Phys. Rev. E 59, 4326 (1999). [22] I. Pagonabarraga, E. Trizac, T.P.C. van Noije, and M.H. Ernst, Phys. Rev. E 65, 011303 (2002). [23] J.M. Montanero and A. Santos, Granular Matter 2,53 (2000). [24] S.J. Moon, M.D. Shattuck, and J.B. Swift, Phys. Rev. E 64, 031303 (2001). 390 M.I. Garcı´a de Soria et al.
[25] V. Garzo ´and J.M. Montanero, Physica A 313, 336 (2002). [26] A. Puglisi, P. Visco, A. Barrat, E. Trizac, and F. van Wijland, Phys. Rev. Lett. 95, 110202 (2005). [27] M.H. Ernst, E. Trizac, and A. Barrat, J. Statist. Phys. 124, 549 (2006). [28] A. Prevost, D.A. Egolf, and J.S. Urbach, Phys. Rev. Lett. 89, 084301 (2002). [29] P. Visco, A. Puglisi, A. Barrat, F. van Wijland, and E. Trizac, Eur. Phys. J. B 51, 377 (2006). [30] S. Aumaıˆtre, J. Farago, S. Fauve, and S. McNamara, Eur. Phys. J. B 42, 255 (2004). [31] J.J. Brey, M.I. Garcı´a de Soria, P. Maynar, and M.J. Ruiz-Montero, Phys. Rev. E 70, 011302 (2004). [32] M.H. Ernst and E.G.D. Cohen, J. Statist. Phys. 25, 153 (1981). [33] N.G. van Kampen, Stochastic Proccesses in Physics and Chemistry (North-Holland, Amsterdam, 1992). [34] P. Maynar, M.I. Garcı´a de Soria, G. Schehr, A. Barrat, and E. Trizac, Phys. Rev. E 77, 051127 (2008). [35] L. Landau and E. Lifshitz, Physical Kinetics (Pergamon Press, Oxford, 1981). [36] F. Coppex, M. Droz, J. Piasecki, and E. Trizac, Physica A329, 114 (2003). [37] J.J. Brey, J.W. Dufty, and M.J. Ruiz-Montero, in Granular Gas Dynamics, edited by T. Poeschel and N. Brilliantov (Springer, Berlin, 2003). [38] M.I. Garcı´a de Soria, P. Maynar, G. Schehr, A. Barrat, and E. Trizac, Phys. Rev. E 77, 051128 (2008). [39] A. Barrat, T. Biben, Z. Racz, E. Trizac, and F. van Wijland, J. Phys. A 35, 463 (2002); A. Barrat, E. Trizac, and M.H. Ernst, J. Phys. A 40, 4057 (2007). [40] G. Costantini, U.M.B. Marconi, and A. Puglisi, J. Stat. Mech.: Theory and Experiments P08031 (2007). Appendix A. Eigenvalue problem for , We consider here the eigenvalue problem for the homogeneous linear Boltzmann–Fokker–Planck operator , defined in (24) ðcÞðcÞ¼ðcÞ:ðA1Þ We are interested in the eigenfunctions and eigenvalues associated with linear hydrodynamics and, to perform the analysis, techniques similar to those in [31,37,38] will be required. Consider first the function 1ðcÞ¼HðcÞ:ðA2Þ When the linearized operator acts on H , we have ðc1ÞHðc1Þ¼Zdc2 T0ðc2,c3Þð1þP12ÞHðc2ÞHðc1Þ þe 2 0 2 @ @c1 2 Hðc1Þ:ðA3Þ Taking into account the equation for H , Equation (14), we obtain the following relation: ðc1Þ 1ðc1Þ¼e 2 0 2 @ @c1 2 Hðc1Þ:ðA4Þ Now let us considerer the function 2ðcÞ¼@ @cHðcÞ:ðA5Þ Taking the derivate in the equation obeyed by H (cw) with respect to w, and subsequently evaluating the result for w¼0, we obtain ðc1Þ 2ðc1Þ¼0:ðA6Þ Finally, we will consider the function 3ðcÞ¼c@ @cHðcÞ:ðA7Þ From the equation obeyed by 3 (c 1 ), we can take the derivate with respect to , and evaluate the result for ¼1. We arrive at an equation for 3 (c 1 ), ðc1Þ 3ðc1Þ¼ðdþ3Þe 2 0 2 @ @c1 2 Hðc1Þ:ðA8Þ From Equations (A4), (A6) and (A8), we can identify two eigenfunctions of . Making use of (A4) and (A8), it appears that ðcÞ1 3 @ @c½cHðcÞþHðcÞ ¼0:ðA9Þ Hence, from Equations (A6) and (A9) we can conclude that the null eigenvalue is (dþ1)-fold degenerate with the eigenfunctions 1ðcÞ¼1 3 @ @ccHðcÞ½þHðcÞ,n2¼@ @cHðcÞ:ðA10Þ Appendix B. Evaluation of the coefficient a 2, i2, i In this appendix, we show that h 2, iðc1Þ 2, iðc2Þjðc1,c2Þi ¼ 0. The integral corresponding to the second term of is simply Zdc1Zdc2c1, ic2, ie 2 0 @ @c1@ @c2 Hðc1ÞHðc2Þ¼e 2 0:ðB1Þ The other term can be written as Zdc1Zdc2c1,ic2,i T0ðc1,c2ÞHðc1ÞHðc2Þ ¼Zdc1Zdc2Hðc1ÞHðc2ÞZd^ rð^ rc12Þð^ rc12Þ½b^ r1c1,ic2,i ¼Zdc1Zdc2Hðc1ÞHðc2ÞZd^ rð^ rc12Þ ð ^ rc12Þ1þ 2ð^ rc12Þc12,i^ ið1þÞ2 4ð^ rc12Þ2^ 2 i ¼ðd1Þ=2 ½ðdþ3Þ=2 12 4dZdc1Zdc2Hðc1ÞHðc2Þc3 12 ¼e 2 0,ðB2Þ which is the desired result. Molecular Physics 391