The total energy E(t) in a fluid of inelastic particles is dissipated through
inelastic collisions. When such systems are prepared in a homogeneous initial
state and evolve undriven, E(t) decays initially as t^{-2} \aprox exp[ -
2\epsilon \tau] (known as Haff's law), where \tau is the average number of
collisions suffered by a particle within time t, and \epsilon=1-\alpha^2
measures the degree of inelasticity, with \alpha the coefficient of normal
restitution. This decay law is extended for large times to E(t) \aprox
\tau^{-d/2} in d-dimensions, far into the nonlinear clustering regime. The
theoretical predictions are quantitatively confirmed by computer simulations,
and holds for small to moderate inelasticities with 0.6< \alpha< 1.Comment: 7 pages, 4 PostScript figures. To be published in Europhysics Letter

Brito, R.

Ernst, M. H.

English

arXiv

© EDP Sciences. It is a pleasure to thank D. Frenkel for an invaluable comment and helpful correspondence. We thank T. P. C. van Noije, J. A. G. Orza, I. Pagonabarraga and M. Hagen for stimulating discussions. The authors also acknowledge financial support from the Offces of
International Relations of Universidad Complutense and Universiteit Utrecht. One of us (RB) acknowledges support to DGICYT (Spain) number PB94-0265.The total energy E(t) in a fluid of inelastic particles is dissipated through inelastic collisions. When such systems are prepared in a homogeneous initial state and evolve undriven, E(t) decays initially as t^(−2) ~ exp[−2ετ ] (known as Haff's law), where “τ” is the average number of collisions suffered by a particle within time t, and ε = 1−α2 measures the degree of inelasticity, with α the coefficient of normal restitution. This decay law is extended for large times to E(t) ~ τ^(−d/2) in d dimensions, far into the nonlinear clustering regime. The theoretical predictions are quantitatively confirmed by computer simulations, and hold for small to moderate inelasticities with 0.6 < α < 1.Offices of International Relations of Universidad ComplutenseUniversiteit UtrechtDGICYT (Spain)Depto. de Estructura de la Materia, Física Térmica y ElectrónicaFac. de Ciencias FísicasTRUEpu

Brito, Ricardo

Docta Complutense

arXiv:cond-mat/9807224v1  [cond-mat.stat-mech]  15 Jul 1998EUROPHYSICS LETTERSEurophys. Lett., (), pp. ()Extension of Haff’s cooling law in granular flows.R. Brito1and M.H. Ernst21 Departamento de F́ısica Aplicada I, Universidad Complutense de Madrid28040-Madrid, Spain2 Instituut voor Theoretische Fysica, Universiteit UtrechtPostbus 80006, 3508 TA Utrecht, The Netherlands(received ; accepted )PACS. 05.20.y – Statistical mechanics.PACS. 05.20.Dd– Kinetic theory.PACS. 81.05.Rm– Porous materials; granular materials.Abstract. – The total energy E(t) in a fluid of inelastic particles is dissipated through inelasticcollisions. When such systems are prepared in a homogeneous initial state and evolve undriven,E(t) decays initially as t−2 ∼ exp[−2ǫτ ] (known as Haff’s law), where τ is the average number ofcollisions suffered by a particle within time t, and ǫ = 1−α2 measures the degree of inelasticity,with α the coefficient of normal restitution. This decay law is extended for large times to E(t) ∼τ−d/2 in d-dimensions, far into the nonlinear clustering regime. The theoretical predictions arequantitatively confirmed by computer simulations, and holds for small to moderate inelasticitieswith 0.6 < α < 1.A fluid of inelastic hard spheres (IHS) in d dimensions, prepared initially in a state ofthermal equilibrium with temperature T0, will remain for an extended period of time (measuredin average number of collisions suffered per particle) in a spatially homogeneous cooling state(HCS), provided the degree of inelasticity ǫ is small. In this state the average kinetic energyper particle E(t) = d2T (t) is decreasing like 1/t2 due to inelastic collisions, as first explained byHaff [1]. However, the HCS with a spatially uniform density and temperature, and a vanishingflow field is unstable against long wavelength spatial fluctuations, and leads finally to a statewith large scale clusters in the density field n(r, t) as well as in the flow field u(r, t) (vortices).This is illustrated in the MD simulations of fig. 1. In this state cooling slows down andlarge deviations from Haff’s law occur, which have not yet been explained in any quantitativemanner [2, 3].In general analytic results for nonlinear decay of rapid granular flows are rare [2]. Thegoal of this letter is to calculate the long time decay of E(t) from nonlinear hydrodynamics,using mode coupling methods. The basic idea is to decompose the quantities of interest, herethe total microscopic energy, or—in the study of long time tails of Green-Kubo formulas—thetotal microscopic momentum-, energy- or particle fluxes, into a superposition of products ofslow hydrodynamic fluctuations of density n(r, t), temperature T (r, t) and flow field u(r, t).Typeset using EURO-TEX2 EUROPHYSICS LETTERSFig. 1. –Snapshots (a) of flow field u(r, τ = 80) and (b) of density field n(r, τ = 160) in an IHSsimulation with N = 50000, α = 0.9, φ = 0.245, prepared in an initial thermal equilibriumstate. Crossover to the nonlinear regime occurs at τc ≃ 70.As one is dealing with fluctuations their amplitude is small in general, so their time decay canbe calculated from linear hydrodynamic modes Ref. [4], provided the fluctuations are stable.However, in the undriven IHS fluid the density fluctuations are growing exponentially likeexp(ǫτ) (clustering instability) [2]. The redeeming feature in the resolution of this paradox isthat the energy only depends on products of modes, in which exponentially growing densityfluctuations are balanced by exponentially decreasing temperature fluctuations.In the case of inelastic hard spheres two colliding particles loose a fraction, ǫ ≡ 1 − α2,of their relative kinetic energy, where ǫ is the degree of inelasticity and α the coefficient ofnormal restitution. Particles move more and more parallel after every collision, and theirrandom motion (temperature) decreases; i.e. the spatial fluctuations at shorter wavelength(which are rapidly damped out through viscosity and heat conduction) are not replenishedwith energy from the more microscopic degrees of freedom, as in the elastic case, wherefluctuations at all wavelengths are kept at thermal noise level, fed by the randomizing elasticcollisions. So, the dynamics selectively suppresses the shorter wave length components of theflow field, as is clearly illustrated in fig. 1a. After a time on the order of the mean free timet0 between collisions, the surviving fraction of the energy E(t) = (1/N)〈∑i12mv2i 〉 is totallystored in hydrodynamic modes (with wavenumbers k < 1, where k is the measured in units ofinverse mean free path l0) of the kinetic energy of the flow and the internal energy (granulartemperature),E(t) ≃ 1N∫dr{12〈ρ(r, t)u2(r, t)〉 + d2〈n(r, t)T (r, t)〉}. (1)The subsequent decay of E(t) is controlled by two different k-regimes [5]: (i) the elasticregime (ǫ < k < 1) at short times (homogeneous cooling regime) with t0 < t < t0/ǫ, and (ii)Extension of Haff’s cooling law in granular flows 3the dissipative regime (k < ǫ) at large times (t ≫ t0/ǫ).In the first regime the decay is controlled by the HCS, where n(r, t) = n, u(r, t) = 0 and ahomogeneous temperature, decaying according to Haff’s cooling law,T (t) =T0(1 + γ0t/t0)2≡ T0 exp(−2γ0τ), (2)where E0 =d2T0 is the initial (equilibrium) energy per particle, γ0 = ǫ/2d, and t0 = 1/ω(T0) isthe mean free time with collision frequency ω(T0) in the initial state. The second equality ineq. (2) defines the average number of collisions per particle within a time t. It is obtained byintegrating dτ = ω(T (t))dt, where ω(T ) ∼√T/l0 with a temperature independent mean freepath l0, given by the Enskog theory for a dense system of hard spheres (d = 2, 3). Haff’s law isillustrated in figs. 2a,b by the straight lines, labeled 1. The HCS is essentially an adiabaticallychanging equilibrium state, parametrized by {n,u = 0, T (t)}.To evaluate the contribution of the long wavelength fluctuations, we linearize the fields ineq. (1) around the HCS, where u(r, t), δn(r, t) = n(r, t) − n and δT (r, t) = T (r, t) − T (t),are small (at least their combined effects; see redeeming features below), and consider Fouriermodes. Retaining up to quadratic terms yields,E(t) =d2T (t) +12N∫dk(2πl0)d{ρ〈|u(k, τ)|2〉 + d〈δn(k, τ)δT (−k, τ〉}, (3)where T (t) is given by Haff’s law. The first term in the integrand is proportional to thestructure factors [(d − 1)S⊥(k, τ) + S‖(k, τ)], and the second one to SnT (k, τ), where thesubscript ‖ refers to the longitudinal component of u and the subscript ⊥ to the transverseones.According to the theory of ref. [3a,b], the fluctuations around the HCS are described byfluctuating hydrodynamic Langevin equations, and obey an ‘adiabatic’ fluctuation dissipationtheorem. The results from the structure factors, Sij(k, τ), obtained by numerically integrat-ing the fluctuating hydrodynamic equations of ref. [3b], could be integrated over the fullhydrodynamic range of wave numbers to obtain the decay of E(t) for all times t > t0.Here, however, we present an analytic method that enables us to calculate the asymptoticdecay for t ≫ t0/ǫ explicitly. For such times only modes in the dissipative regime with k ≪ ǫwill contribute to the integrands in (3). In this regime, the dispersion relations for the modes,relative to the HCS, are well known [5,3c].The structure factors for the (d − 1) transverse components u⊥α and the longitudinal oneul in the range k ≪ γ0 ∼ ǫ can be deduced from ref. [3b], and yield,V −1〈|u⊥α(k, τ)|2〉 = (T0/ρ)e−2∆⊥k2τ(1 + ∆⊥k2/γ0), (4)with a diffusivity ∆⊥ = ν/l20ω, where ν = η/ρ is the kinematic viscosity. The longitudinalstructure factor is given by a similar expression with ∆⊥ replaced by ∆‖. The explicit formof the longitudinal diffusivity (given in ref. [3b] as ∆‖ = γ0ξ2‖/l20) is not needed here. One canverify from the explicit expressions that ∆‖ > ∆⊥, and even ∆‖ ≫ ∆⊥ for γ0 → 0.It should be noted that the dimensionless diffusivities, ∆⊥ and ∆‖, do not depend on thelocal temperature, since transport coefficients and collision frequencies are proportional to√T .By inserting eq. (4) into (3) and performing the k-integrals, one easily finds that thecontributions Euu(t) of the term 〈|u|2〉 to eq. (3) decay as,Euu(t) ≃T02nld0{d − 1(∆⊥)d/2+1(∆‖)d/2}(18πτ)d/2. (5)4 EUROPHYSICS LETTERSNext, we consider the contributions EnT (t) in (3) from the n-T fluctuations. Calculationof the time dependence of δn(k, τ) and δT (k, τ) in the dissipative regime is in principlestraightforward, but technically much more involved. One solves the initial value problem forthe coupled linearized hydrodynamic equations by determining its eigenvalues (to O(k2)-termsincluded), and eigenvectors (to O(k)-terms). The result for k ≪ ǫ is,δn(k, τ) ≃ ikγ0C1e(γ0−∆‖k2)τul(k, 0) + · · ·δT (k, τ) ≃ ikγ0C2e−(γ0+∆‖k2)τul(k, 0) + · · · , (6)where the dots represent subdominant decaying terms, and Cn(n = 1, 2) are some constants.The technical details are available through ref. [6].The first equation in (6) exhibits the clustering instability. It shows that the densityfluctuation is growing at an exponential rate. We are not aware of any publication ona freely evolving rapid granular flows or IHS fluids, that explicitly demonstrate how theclustering instability is driven through the coupling to the long wavelength longitudinal velocityfluctuations—except for the closely related nonlinear slaving mechanism between density andflow field, discussed by Goldhirsch et al. [2]—, whereas all other fluctuations remain boundedby their values in the initial thermal equilibrium state, i.e. the fluctuations in the flow fielddecay diffusively at rates ∆⊥k2 and ∆‖k2, and those in the temperature decay rapidly at arate (γ0 + ∆‖k2).The above discussion also implies that an incompressible granular flow (∇ · u = 0) wouldnot exhibit a clustering instability, at least not in a linear stability analysis. Moreover, aninitial density fluctuation δn(k, 0) would not cause any clustering instability. Furthermore,we observe that 〈δn(k, τ)δT (−k, τ)〉 remains bounded, because the exponential increase of thefirst factor is balanced by an exponential decrease of the second. The resulting integral yieldsEnT (t) ∝ O(τ−d/2−1), which is a subleading asymptotic term when compared to (5).It can be shown similarly that the three mode contribution Enuu(t), neglected in thetransition from eq. (1) to (3) remains bounded, and decays as τ−d, again a subleadingcorrection. The dominant large time contribution Euu(t) is plotted in figs. 2a,b, where thelines labeled 2 and 3 represent the contributions of the transverse and longitudinal modes,respectively. The former contribution is a factor e3 ≃ 20 larger than the latter.The result (6) also shows that there is no balancing of unstable density fluctuations in thestructure factor Snn(k, τ) ∼ 〈|δn(k, τ)|2〉. Consequently, the predictions of the hydrodynamicLangevin equations for Snn breakdown at the crossover from the homogeneous cooling regime,as explained in ref. [3b], whereas those for S⊥, S‖ and SnT , used in this paper are expectedto remain valid until far into the nonlinear clustering regime.In summary, in the decay of the total energy E(t), one can distinguish two differentbehaviors: homogeneous cooling for τ < τc, where EHCS(t) ≃ d2T0 exp(−2γ0τ) is given byHaff’s law, and a diffusive decay for τ > τc, where E(t) ≃ Euu(τ) is solely determinedby the flow field, all long wavelength components of which decay diffusively. The equalityEuu = EHCS defines the crossover time τc, typical values of which are listed in the figurecaptions and table 1.Before comparing the preceding results with MD simulations, a final comment about therelation between the ‘internal’ kinetic time τ and the ‘external’ time t is needed: in thehomogeneous cooling state the relation between both times is correctly predicted by τ =(1/γ0) ln(1 + γ0t/t0) in eq. (2). In the nonlinear clustering regime the relation between bothtimes is not understood. The value of τ(t), as measured in the computer simulations for τ > τc,increases faster than those given by (2). For the cases corresponding to figs. 1a,b the time τExtension of Haff’s cooling law in granular flows 5becomes roughly linear in t for τ larger than 80.In fig. 2a the theoretical prediction for E(t) is compared with computer simulations in adense system (φ = 0.245) of N = 50000 IHS particles at low dissipation (α = 0.9). Thedensity or coverage is defined as φ = 14πNσ2/L2, where σ is the diameter of a disk. Theagreement remains good until far into the nonlinear clustering regime with τ ≫ τc ≃ 70 (seefig. 1a,b). The density field n(r, t) contains small (τ = 80) and large (τ = 160; see fig. 1b)spatial inhomogeneities. The flow field shows both at τ = 80 (see fig. 1a) and at τ = 160well-developed vortex patterns. To describe the crossover regime around τc one would have toevaluate E(t) in the total hydrodynamic time regime (t > t0), by performing the numericalcalculations described in the paragraph below eq. (3). The same good agreement is found atsmall inelasticities (α = 0.975, 0.9) for all densities (φ = 0.05, 0.11, 0.245, 0.4) in the largersystems with N = 50000 particles, as summarized in Table 1. When the number of particlesis reduced to N = 20000 and N = 5000, the deviations steadily increase, suggesting finite sizeeffects.0 40 80 120 160         τ (number of collisions per particle)-10-50ln[E(t)/E(0)]α=0.9φ=0.2452311+2+30 20 40 60 80        τ (number of collisions per particle)-10-50ln[E(t)/E(0)]1231+2+3α=0.6φ=0.05Fig. 2. –Energy decay E(t) vs. τ in simulations (◦), with N = 50000 IHS, (a) at low inelasticity(α = 0.9, φ = 0.245, τc ≃ 70), and (b) at high inelasticity (α = 0.6, φ = 0.05, τc ≃ 23),compared with theoretical prediction (labeled as 1+2+3, solid line), containing homogeneouscooling (1), vorticity diffusion (2) and longitudinal diffusion (3). At the latest simulation times(τl in Table 1) both systems are in the nonlinear clustering regime (see fig. 1b).Figure 2b shows a similar plot of E(t) in a low density system (φ = 0.05) with a rather largeinelasticity (α = 0.6), and consequently a short crossover time τc = 23. There is very goodagreement between theory and simulations over the whole time interval, until the time of thelatest measurement τl = 80. At still larger inelasticities and higher densities (α = 0.4, φ = 0.4)the theoretical prediction is an order of magnitude smaller than the simulated E(t).Table 1 also shows that the deviations at α = 0.6 between theory and simulationsincrease with increasing density (see φ = 0.05, 0.11, 0.4), suggesting that the increase of∆(m, α)/∆(m, 1) with m = (‖,⊥) as α → 0 is much larger at high than at low density. Apossible explanation might be that the Enskog theory for a dense IHS fluid, which is also basedon molecular chaos, breaks down at higher inelasticities. The effects of dynamic correlationsare expected to increase strongly with decreasing α, when frequent inelastic collisions force6 EUROPHYSICS LETTERSα φ τc τ‖ τl Esim/Etheor0.975 0.245 340 27 500 0.65∗0.9 0.245 70 240 160 0.80.9 0.4 68 56 160 10.6 0.05 23 104 80 10.6 0.11 19 3000 160 1.40.6 0.4 15 150 80 40.4 0.4 11 190 80 10Table I. – Comparison Esim(τl) with Etheor(τl) (∗N = 2 × 104; in all others cases N = 5 × 104).the particles to move more and more parallel. But so far ring kinetic theory [7] for calculatingtransport coefficients in IHS fluids has not yet been developed.The present theory is only applicable to thermodynamically large systems, because inderiving (3) the k-sums have been replaced by k-integrals. In the limit of small inelasticities(γ0 → 0), the smallest wave number is k ∼ 1/L, and the restriction k ≪ γ0 might be violatedat small L.For finite systems deviations from these predictions may occur because of interference effectsthrough the periodic boundaries. In elastic hard sphere fluids such effects can be expected fortimes larger than the acoustic traversal time of the system. Here the fastest spreading mode isthe longitudinal diffusion, and one may expect interference effects when the relevant diffusionlength, ξ‖ ≡√8∆‖τ , satisfies ξ‖ >∼ 12L. This implies τ >∼ τ‖ ≡ L2/32∆‖. For the simulationsin fig. 1 and 2a the time of the latest measurements, τl = 160, is of the order of τ‖ = 240, andthe ratio R(τl) ≡ Esim/Etheor ≃ 0.8. In fig. 2b, τl = 80 ≪ τ‖ = 104, and τ‖ = 23, and theratio R(τl) = 1.In summary, the mode coupling theory, applied to a freely evolving IHS fluid at smallinelasticity ǫ, explains how the hydrodynamic modes in the elastic range (ǫ ≪ k < 1)are responsible for Haff’s homogeneous cooling law (2) at times t <∼ t0/ǫ, and those in thedissipative range (k ≪ ǫ) for the diffusive decay, E(t) ∼ τ−d/2 valid for t ≪ t0/ǫ. In fact,the full range of hydrodynamic times (t > t0) can be evaluated numerically by applying thetheory of refs. [3a,b] over the full range of hydrodynamic wavenumbers (k < 1). The solidline in figs. 2a,b, representing the sum of homogeneous cooling and transverse and longitudinalvelocity diffusion, is simply an interpolation between the theoretical long and short time results.The agreement between theory and simulations is excellent over the large time interval fromhomogeneous cooling until far into the nonlinear clustering regime.We also observe that the motion of the IHS fluid at the largest τ -values (see fig. 1a,b and2a,b) is almost ‘frozen in’, as the kinetic energy has typically decayed by a factor of order10−3. Therefore, we do not expect a later crossover to truly nonlinear regime, which wouldimply that the decay law, derived here, is asymptotic indeed. At large inelasticities (α <∼ 0.4)the present theory, which is based on a separation of time scales for small ǫ, is not applicable,and the value of E(t), observed in the simulations, is an order of magnitude larger than theprediction of the present theory.Acknowledgements: It is a pleasure to thank D. Frenkel for an invaluable comment andhelpful correspondence. We thank T.P.C. van Noije, J.A.G. Orza, I. Pagonabarraga andM. Hagen for stimulating discussions. The authors also acknowledge financial support fromthe Offices of International Relations of Universidad Complutense and Universiteit Utrecht.Extension of Haff’s cooling law in granular flows 7One of us (R.B.) acknowledges support to DGICYT (Spain) number PB94-0265.REFERENCES[1] P.K. Haff, J. Fluid Mech. 134, 401 (1983).[2] I. Goldhirsch and G. Zanetti, Phys. Rev. Lett. 70, 1619 (1993); I. Goldhirsch, M-L Tan, and G.Zanetti, J.Sc. Computing 8, 1 (1993).[3] (a) T.P.C. van Noije, M.H. Ernst, R. Brito and J.A.G. Orza, Phys. Rev. Lett. 79, 411 (1997).(b) T.P.C. van Noije, M.H. Ernst and R. Brito, Phys. Rev. E 57, R4891 (1998).(c) J.A.G. Orza, R. Brito, T.P.C. van Noije and M.H. Ernst, Int. J. Mod. Phys. C 8, 953 (1997).[4] M.H. Ernst, E.H. Hauge and J.M.J. van Leeuwen, Phys. Rev. Lett. 25, 1245 (1970).[5] S. McNamara, Phys. Fluids A 5, 3056 (1993).[6] See web page: http://seneca.fis.ucm.es/brito/IHS.html.[7] T.P.C. van Noije, M.H. Ernst and R. Brito, Physica A 251, 266 (1998).

EPL (Europhysics Letters)

Extension of Haff's cooling law in granular flows

R Brito

M. H Ernst

Crossref

The total energy E(t) in a fluid of inelastic particles is 
dissipated through inelastic collisions. 
When such systems are prepared in a homogeneous initial state
and evolve undriven, E(t) decays initially as $t^{-2}
\sim \exp[ - 2\epsilon \tau]$ (known as Haff's law), where τ 
is the average number of collisions suffered by a particle within 
time t, and $\epsilon=1-\alpha^2$ measures the degree of inelasticity,
with α the coefficient of normal restitution. 
This decay law is extended for large times to $E(t) \sim 
\tau^{-d/2}$ in d dimensions, far into the nonlinear clustering 
regime. The theoretical predictions are
quantitatively confirmed by computer simulations, and hold for
small to moderate inelasticities with $0.6< \alpha< 1$

R. Brito

M. H. Ernst

EDP Sciences OAI-PMH repository (1.2.0)

http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.242.6764

Extension of Haff's cooling law in granular flows

Abstract

Similar works

Full text

Available Versions

Docta Complutense

Crossref

EDP Sciences OAI-PMH repository (1.2.0)