Thermalization of an anisotropic granular particle

0.2 MB·English

by  P. Viot

Thermalization of an anisotropic granular particle
P. Viot1 and J. Talbot2
1Laboratoire de Physique Th´eorique des Liquides, Universit´e Pierre et Marie Curie,
4, place Jussieu, 75252 Paris Cedex, 05 France and
2Department of Chemistry and Biochemistry, Duquesne University, Pittsburgh, PA 15282-1530
4
0 Weinvestigatethedynamicsofaneedleinatwo-dimensionalbathcomposedofthermalizedpoint
0 particles. Collisionsbetweentheneedleandpointsareinelasticandcharacterizedbyanormalresti-
2 tutioncoefficientα<1. ByusingtheEnskog-Boltzmannequation,weobtainanalyticalexpressions
for the translational and rotational granular temperatures of the needle and show that these are,
n
in general, different from the bath temperature. The translational temperature always exceeds the
a
J rotational one, though the difference decreases with increasing moment of inertia. The predictions
of thetheory are in very good agreement with numerical simulations of themodel.
6

PACSnumbers: 05.20.-y,51.10+y,44.90+c
]
h
c
e I. INTRODUCTION masses of the two species, the moment of inertia of the
m
needle andthe normalrestitutioncoefficient. Boththese
- temperatures are smaller than that of the bath. The ro-
t Thedissipativenatureofthecollisionsingranularsys-
a tational granular temperature is usually lower than the
temsleadstofundamentallydifferentbehaviorfromtheir
t translationalone, except for very light bath particles for
s thermal analogues. A striking example is the lack of
. which the two are equal.
t energy equipartition between the degrees of freedom in
a We also report essentially exact numerical results ob-
m granular systems. For example experiments [1, 2] and
tainedusingastochasticsimulationmethod. Thetheory
computer simulations [3, 4, 5] haveshownthat in binary
- is in very good agreement with the simulation, which
d mixturesofisotropicinelasticparticles,thegranulartem-
validates the assumption of the Maxwellian shape of the
n peratures of the two species are different.
distribution functions.
o
These results have prompted theoreticians to investi-
c
gate some simple model systems. For example, Mar-
[
tin and Piasecki[6] examined the behavior of a spheri- II. MODEL
2 cal tracer particle immersed in a homogeneous fluid in
v
equilibrium at a temperature T. They showed that the
5 Weexamineatwo-dimensionalsystemconsistingofan
Enskog-Boltzmann equation of the tracer particle pos-
5 infinitely thin needle of length L, mass M and moment
sessesastationaryMaxwellianvelocitydistributionchar-
2 of inertia I immersed in a bath of point particles each
1 acterizedbyaneffectivetemperaturethatissmallerthan of mass m. The vector position of the center of mass of
1 T.
the needleandapointparticlearedenotedbyr1 andr2,
3
Althoughmanygranularsystemscontainparticlesthat respectively. Theorientationoftheneedleisspecifiedby
0
/ are manifestly anisotropic, most studies have been con- a unit vector u1 that points along its axis. Let r12 =
at fined to spherical particles. A notable exception is the r1 r2 and u⊥1 denote a vector perpendicular to u1. A
m work of Aspelmeier, Huthmann and Zippelius [7, 8] that col−lision between a needle and a point occurs when
examinesthe freecoolingofanassemblyofinelasticnee-
-
d dles with the aid of a pseudo-Liouville operator. They r12.u⊥1 =0 (1)
n predictedanexponentiallyfastcoolingfollowedbyastate
o with a stationary ratio of translational and rotational and |λ| < L/2. (see Fig. 1). The relative velocity of the
c energy. This two stage cooling was confirmed by event point of contact V is given by
:
v driven simulations.
i V=v12+λu˙1. (2)
X In this work we examine the breakdown of equipar-
tition in a steady state granular system containing an
r
a anisotropic particle. Motivated by the studies of Mar- The pre- and post-collisionalquantities (the latter are
tin and Piasecki[6] and Aspelmeier, Huthmann and labeled with a prime) obey the usual conservation laws:
Zippelius[7, 8], we consider an infinitely thin inelastic
Total momentum conservation
needle immersed in a bath of point particles. Start- •
ing from the collision rules of this model, we derive Mv′1+mv′2 =Mv1+mv2. (3)
the Enskog-Boltzmann equation. By assuming that the
points are thermalized and that the velocity and angu- Angular momentum conservation with respect to
lar velocity distributions of the needle are Maxwellian, • the point of contact
wederiveanalyticalexpressionsforthetranslationaland
rotational granular temperatures as a function of the Iω1′k=Iω1k+Mλu1 (v′1 v1), (4)
∧ −
2
netic energy,
M
∆E1T = 2 (v1′)2−(v1)2
(cid:0) 1 ∆p2 (cid:1)
=∆p.v1+
M 2
= (1+α)V.u⊥1v1.u⊥1. + M (1+α)2(V.u⊥1)2,
λu − 1 + 1 + λ2 2 1 + 1 + λ2 2
1 m M I m M I
(8)
(cid:0) (cid:1)
and rotational energy,
I
u⊥ r12 ∆E1R = 2 (ω1′)2−(ω1)2
1
= λ(cid:0)(1+α)V.u⊥1(cid:1)(ω1′ +ω1)
− 2 1 + 1 + λ2
m M I
= λ(1+α) V.u⊥1ω1 + λ2(1+α)2(V.u⊥1)2.
− 1 + 1 + λ2 I 1 + 1 + λ2 2
m M I m M I
(9)
(cid:0) (cid:1)
FIG.1: Geometryoftheneedleandapointintheplane: r12
denotesavectorjoining thepointlabeled 2andthecenterof III. PSEUDO-LIOUVILLE EQUATION
theneedle; u1 is a unit vectoralong the axisof theneedle, λ
isthealgebraicdistancebetweenthecenteroftheneedleand
For particles with hard-core interactions, the kinetic
thepoint of impact and u⊥1 is a unit vectorperpendicular to evolution of N-particle distribution f(rN,vN) is de-
the axis of the needle. For a collision to occur one requires
scribed by a pseudo-Liouville operator, where rN is a
that |λ|≤L/2 when thepoint lies on the line defined by the
short-hand notation for the positions (and internal de-
needle; i.e. when Eq.(1) is satisfied.
grees of freedom, i.e. angle θ1 with the x-axis for the
needle) of N particles and vN for their velocities (and
angular velocities for the needle). Originally derived
wherekis aunit vectorperpendicular tothe plane for spheres and elastic collisions[9], the formalism has
such that k=u1 u⊥1. beenextendedtosystemscomposedofinelasticparticles.
∧ Since collisions are instantaneous[10, 11], the pseudo-
LiouvilleoperatorL(rN,vN)=L0(rN)+ T isthe
i<j ij
As a result of the collision, the relative velocity of the sum of the free particle streaming operator L0(rN) and
P
contacting points changes instantaneously according to of the sum of two-body collision operators Tij. Since we
the following relations: areinterestedinthe homogeneousstate,the distribution
function f(v1,ω1) of the needle system is described by
∂f(v1,ω1) dθ1
V′.u⊥1 = αV.u⊥1 (5) =N dv2 dr2T12f(v1,ω1,v2)
− ∂t 2π
V′.u1 =V.u1 (6) Z Z Z (10)
where N is the total number of points, f(v1,ω1,v2) the
distribution function of the needle and a point, and T12
where α denotes the normal restitution coefficient. For
is the collision operator between a needle and a point.
the sake of simplicity we have taken the tangentialresti-
The assumption of molecular chaos yields a factor-
tution coefficientequalto one. This choiceis reflectedin
ization of the two body distribution f(v1,ω1,v2) =
the form of Eq. (6).
f(v1,ω1)Φ(v2), where Φ(v2) denotes the point distribu-
By combining Eqs. (2)- (5) one obtains, after some tionfunction. Sincethepointparticlesarethermalized,a
algebra, the change of the needle momentum ∆p = stationarystateforthesystem(tracerneedleandpoints)
M(v1′ v1) canbereached,andonefinallyobtainsfortheneedledis-
−
tribution
∆p.u⊥1 =−(11++α1)V+.uλ⊥12 . (7) Z dθ1Z dv2Z dr2T12f(v1,ω1)Φ(v2)=0 (11)
m M I where
mv2
The inelastic collision leads to a loss of translational ki- Φ(v2)∝exp − 2T2 (12)
(cid:18) (cid:19)
3
where T is the temperature of the bath.
To build the point-needle collision operator, T12, one
must include the change in quantities (i.e. velocity and ... dr2dv1dv2dω1dθ1Θ(L/2 λ)δ(r12.u⊥1 0+)
−| | | |−
angular momentum) produced during the infinitesimal Z Z
tfrimome iznetreorvoanllyofifththeecotwlliosiponar.tiTclheissaorpeeirnatcoornitsacdtiffaenrdenitf d|r1d2t.u⊥1| Θ − d|r1d2t.u⊥1| f(v1,ω1)Φ(v2)∆E1T =0.
(cid:12) (cid:12) (cid:18) (cid:12) (cid:12)(cid:19)
theparticleswereapproachingjustbeforethecollision[7]. (cid:12) (cid:12) (cid:12) (cid:12) (15)
(cid:12) (cid:12) (cid:12) (cid:12)
The explicit form of the operator is (cid:12) (cid:12) (cid:12) (cid:12)
Asimilarequationcanbewrittenforrotationalkinetic
T12 Θ(L/2 λ)δ(r12.u⊥1 0+) egniveerngyb.yTheneedledistributionfunctionf(v1,ω1)isthen
∝ −| | | |−
d|r12.u⊥1| Θ d|r12.u⊥1| (b12 1), (13) Mv12 Iω12
(cid:12) dt (cid:12) (cid:18)−(cid:12) dt (cid:12)(cid:19) − f(v1,ω1)∝exp −2γ T − 2γ T , (16)
(cid:12) (cid:12) (cid:12) (cid:12) (cid:18) T R (cid:19)
(cid:12) (cid:12) (cid:12) (cid:12)
(cid:12) (cid:12) (cid:12) (cid:12)
where b12 is an operator that changes pre-collisional where γT and γR are the ratios of the translational and
quantities to post-collisional quantities and Θ(x) is the rotationalneedle temperatures to the bath temperature,
Heaviside function. respectively.
The other terms of the collision operator corre- ByinsertingEq.(16)inEq.(15)andinthecorrespond-
spond to the necessary conditions of contact Θ(L/2 ing equation for the rotational energy one obtains, after
|λ|)δ(|r12.r⊥1|−0+), and approach Θ − d|r1d2t.u⊥1| . − sAopmpeenteddixioAus),btuhtestforalliogwhtinfogrwseatrdofceaqlcuualtaitoinosn(outlinedin
For an isotropic tracer particle in(cid:16)an(cid:12) isotropi(cid:12)c(cid:17)par-
(cid:12) (cid:12)
ticle bath, Martin and Piasecki[6] solve(cid:12)d the st(cid:12)ation- 1 √1+akx2 1+α 1 (1+akx2)3/2
dxb = dx ,
ary Enskog-Boltzmann equation, showing that the ve- 1+kx2 2 (1+kx2)2
0 0
locity distribution of the tracer particle remains gaus- Z Z (17)
sian in a thermalized bath. It is worth noting that
1 √1+akx2 1+α 1 x2(1+akx2)3/2
for finite dilutions, the velocity distribution function is dxax2 = dx ,
not Maxwellian[3], but the deviations from the gaussian Z0 1+kx2 2 Z0 (1+kx2)2
shape of the distribution function of the velocities (that (18)
can be calculated in a perturbative way) yields small
where
corrections to the estimate of the granular temperature.
More significant deviations are expected for finite dilu- L2
tion systems[12]. k = , (19)
4I 1 + 1
For a mixture of a needle and points, one can show m M
that a Maxwell distribution of the needle angular and (cid:0) (cid:1)
translational momenta cannot be a solution of the sta- (M +m)
a=γ , (20)
tionary Enskog-Boltzmannequation even for the case of R
M +mγ
T
an infinitely diluted needle. This is due to the fact that
the change of momentum depends on the location of the and
point of impact on the needle (the right-hand side of
M +m
Eq. (7) depends on λ). However, we impose gaussian b=γ . (21)
T
M +mγ
distributions for the translational and angular velocities T
of the needle. We show below that this is a very good
Explicit expressions for the integrals appearing in
approximation.
Eqs. (17) and (18) are given in Appendix B. Eq. (18)
In the stationary state, the loss of translational and is an implicit equation for a that, for a given value of
rotational kinetic energies is on average zero and can be α, can be solved with standard numerical methods. b is
expressed by means of the collision operator then easily obtained by calculating the ratio of integrals
of Eq. (17). Finally, from the values of a and b, γ and
T
∆E1T = T12f(v1,ω1)Φ(v2)E1T(v1) =0 γR can be obtained from Eqs.(20)-(21).
h i h i
∆E1R = T12f(v1,ω1)Φ(v2)E1R(ω1) =0 (14)
h i h i
IV. RESULTS
where the brackets denote an average over the indepen-
dent variables. This corresponds to taking the second We first note that for the elastic case, α=1, we have
moments of the distribution function of the needle. Af- the solution a = b = 1 which gives γ = 1 and since
T
ter substitution of the collision operator (Eq.(13)), one a/b = γ /γ , γ = 1 which corresponds to the equilib-
R T R
obtains explicitly for the translational kinetic energy rium case[13, 14].
4
1 1
0.8 0.8
0.6 0.6
γγR,T γγT,R
0.4 0.4
0.2 0.2
0 0
0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1
α α
FIG. 2: Ratio of the translational (full curve) γT and rota- FIG. 3: Ratio of the translational (full curve) γT and rota-
tional γR (dashed curve) granular temperature to the tem- tional γR (dashed curve) granular temperature to the tem-
perature of the bath versusthe normal restitution coefficient perature of thebath versusthe normal restitution coefficient
α for a homogeneous needle with M =m. α for an inhomogeneous needle(I =ML2/16)with M =m.
Figure 2 shows the ratio of the translational tempera- [6]. We conjecture that the result is independent of the
ture to the bath temperature and of the rotational tem- shape of the tracer particle and the dimension.
peraturetothebathtemperatureforahomogeneousnee- Figure 4 displays the variation with the mass ratio of
dle (I = ML2/12) whose mass is equal to the mass of a thetranslationalandrotationaltemperatureofthetracer
bath particle. As the normal restitution coefficient de- needle. Note that when the mass of the bath particles
creases from 1 to 0 both ratios decrease monotonically is very small compared to the needle particle, the two
from 1 to a strictly positive value, such that the trans- curves go to the limit given by Eq. (22).
lational temperature is always largerthe rotationaltem-
perature. An analogoussituation has been noted for the
free cooling of needles in three dimensions[7]. V. SIMULATION
Figure 3 shows the same ratios for an inhomogeneous
needle (I =ML2/16) corresponding to relatively lighter In order to assess the accuracy of the theory we de-
ends and a heavier center. Again the translational tem- veloped a Direct Simulation Monte Carlo (DSMC) code
perature is always larger than the rotational tempera- to obtain numerical solutions of the non-linear homoge-
ture, and the differences are greater than for a homo- neous Boltzmann equation. The simulation generates a
geneous needle. The difference between the two tem- sequence of collisions with a variable time interval be-
peratures vanishes in the limit of very large moment of tween each.
inertia. Contrary to the free cooling state of needles [7], At the beginning of each step, we first calculate the
thereisno“critical”valueofthemomentofinertiaabove total flux of colliding particles on both sides of the nee-
whichtherotationaltemperatureislargerthanthetrans- dle which is a function of the component of the velocity
lational temperature. of the center of mass normal to the length of the nee-
Wenowreturntoconsiderationofahomogeneousnee- dle v1.u⊥1 and the angular velocity ω. Next a waiting
dle and investigate the effect of varying the mass ratio,
timetothe nextcollisionissampledfromanexponential
m/M, at constant normal restitution coefficient. When distribution with rate constant equal to the total flux.
m/M 0,k 0,whichgivesa=b=(1+α)/2. Hence, The needle is rotated at its current angular velocity to
→ →
forverylightbathparticles,equipartitionbetweentrans- the point of collision. The location of the collision is
lational and rotational granular temperatures is asymp-
then selected with a probability proportional to the po-
totically obtained and we have sition dependent flux. The ratio of the flux on the left
hand side to the total flux is then compared to a uni-
1+α
γR =γT = . (22) form random number between zero and one. If the ratio
2
is greater than the random number, the collision occurs
Thissamelimitisobtainedforasphericaltracerparticle on the right; otherwise is occurs on the left. The normal
5
1 1
0.9
0.1
0.8
0.7
γγR,T P(v)x 0.01
0.6
0.5
0.001
0.4
0.3 0.0001
0 2 4 6 8 10 12 -4 -2 0 2 4
m/M vx
FIG. 4: Ratio of the translational (full curve) γT and rota- FIG. 5: Distribution of the normal velocity component of
tional γR (dashed curve) granular temperature to the tem- theneedlefortwovaluesofthenormalrestitutioncoefficient,
perature of the bath versus the ratio of the masses m/M for α=0.5andα=0.9(broadercurve). Thesolid(noisy)curves
a homogeneous needle (I =ML2/12) and α=0.9. are the simulation results and the dashed curves correspond
to Gaussian distributions of velocity evaluated for the corre-
sponding theoretical granular temperature.
component of the colliding point is then selected from a
gaussianflux. Finally, the collisionrules are applied and
thenormalvelocityandangularvelocityoftheneedleare
updated.
The simulation results for 107 collisions per run are
compared with the theory in Figs. 2- 3. Although the
theory is not exact, the agreement is excellent. The dis-
crepancy between simulationand theory is never greater
0.1
than1%whateverthe valueofthe restitutioncoefficient.
InFig.4,whenthemassoftheneedleincreases,thegran-
ular temperatures obtained by simulation are slightly
smaller than the analytical result.
0.01
Fig. 5 and Fig. 6 show the translational and angular
ω)1
velocity distributions, respectively, for two different val- P(
ues of the normal restitution coefficient (0.5 and 0.9).
On the scale of the plot, the simulated distributions are
0.001
in quantitative agreement with gaussian distributions at
thecorrespondingtheoreticalgranulartemperature,con-
firmingthe validityofemployingthelatterinthetheory.
0.0001
-10 -5 0 5 10
ω
VI. CONCLUSION 1
We have shown that the translational and rotational
granular temperatures of an anisotropic tracer particle FIG. 6: Same as 5, except the angular velocity is shown.
immersed in a bath of point particles depend on the ra-
tio of the masses of a point and the needle and the mo-
ment of inertia. They are, in general, not equal and dif-
fer from the bath temperature. In addition, we found
that equipartition is obtained asymptotically, regardless higher dimensions and for arbitrary shaped anisotropic
of the normal restitution coefficient, for very light bath particles. It should be possible to test this prediction
particles. We expect this to be a general feature, i.e. in experimentally.
6
APPENDIX A: NEEDLE AVERAGE ENERGY and
LOSS
s=(s1,s2,s3)
Asforbinarymixturesofspheres[3],weuseagaussian χ.u⊥1,ν.u⊥1,ξ (A8)
ansatz for the distribution functions and introduce two
differenttemperaturescorrespondingtothetranslational By inserting Eq. (8) in(cid:0)Eq.(A6), the(cid:1)averageenergy loss
and rotational degrees of freedom of the needle. The can be rewritten as
homogeneous distribution functions of the needle and of
the points are then given respectively by dλ dsexp s2 G.sΘ(pG.s)Θ λ L
− | | | |− 2
f(v1,ω1) exp Mv12γT−1 Iω12γR−1 , (A1) pX=±1Z Z (cid:0) (cid:1) (cid:18) (cid:19)
∼ − 2T − 2T (1+α)G.s 2T mγ
(cid:18) (cid:19) "−m1 + M1 + λI2sMγT +m(cid:18)γTs1+r MTs2(cid:19)
Φ(v2)∼exp(cid:18)−m2Tv22(cid:19), (A2) +2M(1+1 α+)2(1G+.s)λ22 2# (A9)
m M I
where T is the temperature of the bath, γ and γ the
T R
(cid:0) (cid:1)
ratioofthetranslational(androtational)temperatureof By defining a new coordinate system in which the z-axis
the needle to the bath temperature. is parallel to G, one find that the integrals of Eq. (A9)
We introduce the vectors χ and ν such that involve gaussian integrals of the form
1 π
χ= 2T(Mγ +m)(Mv1+mv2) (A3) dsexp −s2 (|G|sz)2Θ(±sz)Gisz = 8|G|2Gi (A10)
T
Z
(cid:0) (cid:1)
p mM and
ν = (v1+γTv2) (A4)
s2T(MγT +m)γT
π
dsexp s2 (Gs )3Θ( s )= G3 (A11)
z z
The scalar product V.u can be expressed as − | | ± 8| |
⊥1 Z
(cid:0) (cid:1)
which finally leads to Eq. (17). The equation for rota-
2T
V.u = (γ 1)χ.u tional energy is derived following exactly the same pro-
⊥1 smMγT T − ⊥1 cedure.
(cid:2)
m M
+√γT + ν.u⊥1 (A5)
M m
r r ! # APPENDIX B: INTEGRALS
Letusintroduceξ =ω1 2TIγR. Thetranslationalenergy The integrals appearing in Eqs. (17) and (18) can be
loss is given by the formqula evaluated explicitly. For completeness, we give below
their analytical expressions
dλ dθ1 dχ dν dω1
p= 1Z Z Z Z Z 1 √1+akx2
eXxp± χ2 ν2 ξ2 V.u Θ(pV.u )Θ λ L ∆ET. I1 =Z0 dx 1+kx2
− − − | ⊥1| ⊥1 | |− 2 1 a
(cid:18) (cid:19) = ln(√ak+√1+ak)
(cid:0) (cid:1) (A6) k
r
SinceEq.(A5)dependsonlyonχ.u andν.u ,onecan 1 a (1 a)k
⊥1 ⊥1 + − arctan − (B1)
freely integrateoverthe directionofu1 forthe vectorsχ r k r 1+ak !
and ν. The integration over θ1 can be easily performed.
If we introduce the three dimensional vectors G and s The integral I2 is defined as
with components:
1 (1+akx2)3/2
G=(G1,G2,G3) I2 = dx (1+kx2)2 (B2)
0
Z
2T 2Tγ m+M
T
= (γT 1), , and satisfies the equation
sMγT +m − sMγT +m √mM
λ 2TγR (A7) I2 =I1 (a 1)k ∂ 01dx√11++daxk2x2 (B3)
r I ! − − R ∂d (cid:12)(cid:12)d=k
(cid:12)
(cid:12)
(cid:12)
7
which gives Finally, one has
a3/2 (1 a)√1+ak 1 x2(1+akx2)3/2
I2 = √k ln(√ak+√1+ak)+ −2(1+k) I4 =Z0 dx (1+kx2)2 (B7)
1+a 2a2 (1 a)k
− arctan − (B4) which satisfy the algebraic equation
2 k(1 a) r 1+ak !
−
I2+kI4 =aJ1+(a 1)I1 (B8)
p −
1 x2√1+akx2 which gives
I3 = dx 1+kx2
Z0 √1+ak(ak 1+2a)
= 1 1dx 1+akx2 I1 (B5) I4 = 2k(k+−1)
k −
√a(cid:18)kZ+0 1 1p 2a (cid:19) + √a(3−4a)ln(√ak+√1+ak)
= − +ln(√ak+√1+ak) 2k3/2
2k √ak3/2
√1 a(1 4a) (1 a)k
√1−aarctan (1−a)k (B6) + −2k3/2− arctan r 1−+ak ! (B9)
− k3/2 r 1+ak !
[1] R. D. Wildman and D. J. Parker, Phys. Rev. Lett. 88, Rev. E60, 654 (1999).
064301 (2002). [9] M. H. Ernst, J. R. Dorfman, W. R. Hoegy, and J. M. J.
[2] K. Feitosa and N. Menon, Phys. Rev. Lett. 88, 198301 van Leeuwen, Physica 45, 127 (1969).
(2002). [10] T. van Noije and M. Ernst, in Granular Gases, edited
[3] A.Barrat and E. Trizac, Granul. Matter 4, 57 (2002). by S. Luding and T. Poschel (Springer, Berlin, 2001),
[4] A.BarratandE.Trizac,Phys.Rev.E66,051303(2002). Chap. KineticTheory of Granular Gases, p. 3.
[5] P.E. Krouskop and J. Talbot, cond-mat/0303263. [11] J.J. Brey,F. Moreno, and J.W. Dufty,Phys. Rev.E 54,
[6] P. A. Martin and J. Piasecki, Europhys. Lett 46, 613 445 (1996).
(1999). [12] T. Biben, P.A.Martin, and J. Piaseski, Physica A 310,
[7] T. Aspelmeier, T. M. Huthmann, and A. Zippelius, in 308 (2002).
Granular Gases, edited by S. Luding and T. Poschel [13] D. Frenkel and J.F. Maguire, Phys. Rev. Lett. 47, 1025
(Springer, Berlin, 2000), Chap. Free cooling of particles (1981).
with rotational degrees of freedom, p.31. [14] D. Frenkeland J. Maguire, Mol. Phys. 49, 503 (1981).
[8] M. Huthmann, T. Aspelmeier, and A. Zippelius, Phys.