Energy transport in a one-dimensional granular gas
Ítalo’Ivo Lima Dias Pinto1 , Alexandre Rosas1, and Katja Lindenberg2
1
arXiv:0902.2401v1 [cond-mat.stat-mech] 13 Feb 2009
Departamento de Fı́sica, CCEN, Universidade Federal da Paraı́ba,
Caixa Postal 5008, 58059-900, João Pessoa, Brazil
2
Department of Chemistry and Biochemistry, and Institute for Nonlinear Science,
University of California San Diego, La Jolla, CA 92093-0340
We study heat conduction in one-dimensional granular gases. In particular, we consider two
mechanisms of viscous dissipation during inter-grain collisions. In one, the dissipative force is proportional to the grain’s velocity and dissipates not only energy but also momentum. In the other,
the dissipative force is proportional to the relative velocity of the grains and therefore conserves
momentum even while dissipating energy. This allows us to explore the role of momentum conservation in the heat conduction properties of this one-dimensional nonlinear system. We find normal
thermal conduction whether or not momentum is conserved.
PACS numbers: 45.70.-n,05.20.Dd,05.70.Ln,83.10.Rs
I.
INTRODUCTION
Energy transport in low-dimensional systems can be
pathological in the sense that Fourier’s law for heat
conduction might break down [1, 2, 3, 4]. In threedimensional solids, the energy (heat) transported is governed by Fourier’s law, which says that the heat flux J
is proportional to the gradient of the temperature,
J = −κ∇T,
(1)
where κ is the thermal conductivity. The thermal conductivity in Fourier’s law is independent of system size
and of time. This equation assumes that a local equilibrium is established at each time, so that one can define the local energy flux J(x, t) and temperature T (x, t),
when the temperature gradient lies along the x direction [4]. In one-dimensional systems one might expect
a similar equation, with the gradient replaced by the
derivative of T (x, t) with respect to x, and J simply
being the energy flux per unit time. However, it has
been observed that in many one-dimensional models the
thermal conductivity varies with the system size N as
κ ∼ N α , with α lying between 0.32 and 0.44 [4, 5, 6]. Selfconsistent mode coupling analysis for one-dimensional
nonlinear chains [7] predicts a universal value of α = 2/5,
while renormalization group analysis based on hydrodynamical models supports a universality class with α =
1/3. Recently, Mai et al. [5, 6] revisited the problem and
obtained the same α = 1/3 behavior for one-dimensional
nonlinear chains as obtained for the fluid-like systems.
They argue that the dynamic equations are those of a
fluid at all length scales even if the static order of the
chain extends to very large N . They explain the disagreement with most numerical results as a consequence
of the numerical difficulties in reaching asymptotic behavior, which requires extremely large N , but their conclusions have in turn been disputed [6]. Whether there
is a single universality class or there are two, there is a
connection between heat conduction and diffusion in onedimensional systems [8], that is, normal diffusion leads to
normal heat conductivity (κ ∼ N 0 ), whereas anomalous
transport is associated with κ ∼ N α , with α > 0 (< 0)
for super(sub)diffusion. Thus, superdiffusion is associated with a divergent thermal conductivity.
Another key concept in the understanding of the thermal conductivity in 1d systems is momentum conservation [4]. In fact, momentum conservation usually implies
the divergence of the thermal conductivity in 1d [1, 9],
and yet it is not necessary for the occurrence of anomalous diffusion [10]. An exception to the momentum conservation rule is a chain of particles interacting via the
nearest-neighbor potential V (qi − qi−1 ) = 1 − cos(qi −
qi−1 ), the so-called rotator model [11]. This momentumconserving model exhibits normal transport behavior to
very high numerical accuracy, and is said to occur because the rotator model “cannot support a nonvanishing
pressure, and thus infinite-wavelength phonons cannot
carry energy” [10]. In [4] this behavior is ascribed to the
periodicity of the potential and the associated independent jumps from valley to valley.
In this contribution we discuss heat conduction in onedimensional granular gases. Granular gases are dissipative systems, the associated energy loss usually being
modeled by introducing a coefficient of restitution as a
parameter in the description of granular collisions. Instead, we follow a more dynamical approach and introduce energy dissipation explicitly via viscous terms in the
equations of motion for the granules. There are a number
of different sources and descriptions of viscous effects in
the literature [12, 13, 14, 15, 16, 17, 18]. The different
ways in which they appear in the dynamical descriptions
suits our particular focus of interest, which is the role
of momentum conservation in the heat transport process. In particular, one way to introduce viscosity conserves momentum [19, 20, 21, 22], while the other does
not [18, 23, 24]. Interestingly, we establish that thermal
conduction exhibits normal behavior in both cases, that
is, that there is no divergence as the system size increases.
In Sec. II we describe the two dissipative models, and
in Sec. III we present simulation results to characterize
the transport of heat through the granular gas in both
cases. We conclude with a short summary in Sec.IV.
2
II.
DISSIPATIVE MODELS
We consider N identical granules constrained to move
on a line between two walls at different temperatures.
The granules move freely except during collisions either
with the walls or with one another. The system is a onedimensional “granular gas” because the distance between
the walls is much greater than the space occupied by the
granules. The walls act as heat baths, that is, whenever a
granule collides with a wall at temperature T , its energy
is absorbed by the wall and it acquires a new velocity
away from the wall according to the probability distribution [25]
P (v) =
|v|
exp −v 2 /2kB T .
T
(2)
Interparticle collisions are governed by the power-law potential
a |δ|n
V (δk,k+1 ) = n
δ ≤ 0,
k,k+1 ,
(3)
V (δk,k+1 ) = 0,
δ > 0.
Here δk,k+1 ≡ yk+1 − yk , a is a prefactor determined
by Young’s modulus and Poisson’s ratio, and the principal radius of curvature R of the surfaces at the point
of contact [26, 27]; and yk is the displacement of granule k from its position at the beginning of the collision.
The exponent n is 5/2 for spheres (Hertz potential), it
is 2 for cylinders, and in general depends on geometry.
We stress that the one-sided (only repulsion) granular
potential even with n = 2 is entirely different from a
two-sided (repulsion and attraction) harmonic potential.
The one-sidedness of the potential leads to analytic complexities even in the dissipationless case [24, 28, 29],
and even greater complexities in the presence of dissipation [18, 19, 20, 24].
In this paper we explore the low density limit, which
leads to enormous analytical and computational simplifications. The low density feature of these approximations
is implemented via the assumption that the collisions are
always binary [30], that is, that only two granules at a
time are members of any collision event, and that at any
moment of time there is at most one collision.
Our approach starts with the equations of motion of
the particles labeled by index k, and so we write yk as a
function of time τ . It is convenient to deal with scaled
position and time variables xk and t, related to the unscaled variables yk and τ as follows,
1/n
1/n
mv02
1 mv02
yk =
xk ,
τ=
t.
(4)
a
v0
a
Here m is the mass of the granules, and the velocity v0 is
an arbitrary choice in terms of which other velocities are
expressed. We also introduce a scaled friction coefficient
1/n
γ̃
mv02
γ=
,
(5)
mv0
a
where γ̃ is the friction coefficient for the unscaled variables. In the low density limit we only need to consider the equations of motion for two colliding granules
k = 1, 2. Furthermore, in this paper we only consider
cylindrical grains, which leads to considerable simplification while still capturing the important general features
of the system that we seek to highlight. We stress that
the one-sided granular potential (i.e., one with only repulsive interactions) even with n = 2 is entirely different
from a two-sided harmonic potential.
Consider a viscous force that is proportional to the relative velocity of the colliding granules. Such a force has
been considered not only theoretically [12, 19, 20] but it
has also been observed experimentally [21, 22]. While one
might be tempted to think of this force as arising because
one grain rubs against the other, the actual mechanism
is more complicated and involves the medium that surrounds the granules [21, 22]. The exact mechanism is
still a matter of conjecture. In any case, the appropriate
equations of motion in this case are
ẍ1 = [−γ(ẋ1 − ẋ2 ) − (x1 − x2 )] θ(x1 − x2 ),
ẍ2 = [γ(ẋ1 − ẋ2 ) + (x1 − x2 )] θ(x1 − x2 ),
(6)
where a dot denotes a derivative with respect to t. The
Heaviside function is defined as θ(x) = 1 for x > 0,
θ(x) = 0 for x < 0, and θ(0) = 1/2. It ensures that
the two particles interact only when in contact, that is,
only when the particles are loaded. The post-collisional
velocities (called u below) can be written in terms of the
velocities of the two granules at the beginning of the collision (called v) as
1
1
1 − e−γt0 v1 +
1 + e−γt0 v2 ,
2
2
1
1
−γt0
u2 =
v1 +
(7)
1+e
1 − e−γt0 v2 ,
2
2
p
where t0 = π/( 2 − γ 2 ) is the duration of the collision.
Since u1 + u2 = v1 + v2 , the momentum of the center of
mass is conserved.
A model in which the center of mass velocity is not
conserved is one governed by the equations of motion [18,
24]
u1 =
ẍ1 = [−γ ẋ1 − (x1 − x2 )] θ(x1 − x1 ),
ẍ2 = [−γ ẋ2 + (x1 − x2 )] θ(x1 − x1 ).
(8)
Here one might think of the viscosity arising from an interaction with the medium. However, such an interaction
would produce a damping term not only during a collision but also while the granules are moving independently
between collisions. In the low density limit the granules
would hardly collide before stopping entirely unless the
viscosity is extremely low, and one can then not talk of
heat transport along the one-dimensional system. Thus,
this model, while perhaps not realistic in any sense, is a
“toy model” in which momentum is not conserved while
still capable of supporting energy transport and therefore
3
u1 + u2 = [v1 + v2 ] e−γt0 ,
(10)
after the collision, so that the momentum is not conserved.
Caution must be exercised in the choice of the damping
coefficient γ. If it is too large, energy acquired from either
wall is simply dissipated before it crosses the system, in
which case the problem changes from one of energy transfer from one wall to the other through the granular gas
to that of two walls pumping energy into the granular
“sink.” To estimate the limit on the damping coefficient
we can imagine a sequence of collisions whereby a granule
starting from the hot wall with kinetic energy E0 collides
with the next granule, which in turn collides with the
next one, and so on, until the last granule, whose energy is EN , collides with the cold wall. The change in
the kinetic energy due to each collision in the momentum
conserving model is
1
1
1
∆E ≡ (u21 +u22 )− (v12 +v22 ) = − 1 − e−2γt0 (v1 +v2 )2 ,
2
2
4
(11)
and in the momentum non-conserving model it is
1
1
1
∆E ≡ (u21 +u22 )− (v12 +v22 ) = − 1 − e−γt0 (v1 +v2 )2 .
2
2
4
(12)
In both cases the kinetic energy before and after the collision are related by an expression of the form Eaf ter =
Ebef ore (1 − bγt0 ) where b is a velocity-dependent dimensionless quantity of order unity, and where we have assumed that γt0 ≪ 1. As an order of magnitude estimate
it is sufficient to write EN ∼ E0 (1 − γt0 )N . The average
energy of a hot granule is kB T1 /2 and in order for energy
to flow across the system the energy of the granule that
arrives at the cold wall must be greater than kB T2 /2.
Consequently, we estimate that a requirement for energy
transport is that
1/N
T2
γ .1−
.
(13)
T1
We now wish to explore the following specific questions. Is there a finite thermal conductivity for either
or both of these models? If so, which model has a
higher/lower thermal conductivity? How does the thermal conductivity depend on temperature? On viscosity?
On system size? These are the questions we address in
the next section.
III.
3
E (x 10 )
relevant to our question. Indeed, we can again calculate
the post-collision velocities,
i
i
1 h −γt0
1 h −γt0
− e−γt0 /2 v1 +
+ e−γt0 /2 v2 ,
u1 =
e
e
2
2
i
i
1 h −γt0
1 h −γt0
−γt0 /2
u2 =
v1 +
+e
− e−γt0 /2 v2 ,
e
e
2
2
(9)
p
where t0 = 2π/( 8 − γ 2 ) is again the duration of the
collision. It is easy to verify that in this case the momentum of the center of mass decays from v1 + v2 before a
collision to
NUMERICAL SIMULATIONS
45
40
35
30
25
20
15
10
5
0
-5
0
1
2
3
4
5
6
7
8
9
10
5
time (x 10 )
FIG. 1: (Color online) From top to bottom: energy injection
by the hot wall (red), energy absorption by the cold wall (magenta), energy injection by the cold wall (blue), and energy
absorption by the hot wall (green), for a momentum nonconserving system of N = 500 particles between walls maintained at Tl = 6 and Tr = 3 and with γ = 0.0002. The dashed
lines are the linear regressions of these curves (after the transient). They yield the rate of energy injection/absorption.
The initial temperature of the granular gas is T = 1.
The low density limit allows us to use an event driven
algorithm in our simulations. As indicated earlier, when
a granule collides with a wall at temperature T its energy is absorbed by the wall and it acquires a new energy as determined by the velocity distribution given in
Eq. (2). The initial temperature of the granular gas is arbitrarily set to T = 1, and we allow the system to arrive
at a steady state before beginning our “measurements.”
While ideally one would want the temperature gradient
in the steady state to be strictly linear as assumed in
the linear response theory that leads to the Fourier law,
the subtleties encountered in 1d systems are well-known
and observed pretty much no matter how the thermalization is implemented [6, 25] and the resulting gradients
are only strictly linear away from the boundaries (nonlinear behavior typically sets in near the boundaries). In
any case, we can address the question of whether or not
there are system size divergences in the transport of heat.
That we do in fact achieve a steady state can be seen, for
example, in Fig. 1, where we show the energy absorbed
and injected by each wall as a function of time for the
dissipating momentum model. After a short transient,
the rate of energy injection and absorption become constant, as they should in a steady state. We ran tests
to ascertain that the initial temperature of the granular gas is not important for this equilibration, that is,
we find that the rates of energy injection and absorption
are independent of the initial velocity distribution of the
granular gas. The momentum-conserving model equili-
4
brates equally well, also independently of initial condition. Henceforth we set the initial temperature of the
granular gas to unity.
κ (γ)
0.16
0.14
0.12
0.1
J
0.08
0.06
0.018
0.017
0.016
0.015
0.014
0.013
0.012
0.011
0.01
0.009
5
0.04
10
15
20
-5
γ (x 10 )
0.02
0
-0.02
0
2
4
6
8
10
12
∆T
FIG. 2: (Color online) Rates of energy transmission as a function of the temperature gradient. In this simulation the viscosity coefficient was set to γ = 0.0002 and the system is
composed of 500 granules. The (red) plus signs correspond to
the momentum conserving system, while the (blue) stars correspond to the momentum dissipating system. In this figure,
the cold wall temperature was set to 3.
FIG. 4: (Color online) Thermal conductivity as a function of
the viscosity coefficient. In this simulation the temperature
of the cold wall is set to 3 and that of the hot wall to 6.
The system is composed of 500 granules. The (red) plus signs
correspond to the momentum conserving system, the (blue)
stars to the momentum dissipating system.
0.05
0.04
κ (T)
J
0.03
0.02
0.019
0.018
0.017
0.016
0.015
0.014
0.013
0.012
0.011
0.01
0.009
0.02
0.01
100 200 300 400 500 600 700 800 900 1000
N
2
4
6
8
10
12
T
FIG. 3: (Color online) Thermal conductivity as a function of
the temperature. In this simulation the viscosity coefficient
was set to γ = 0.0002 and the system is composed of 500
granules. The (red) plus signs correspond to the momentum
conserving system, the (blue) stars to the momentum dissipating system.
Next, we discuss the dependence of the energy flux on
the temperature gradient. The energy flux J = dE/dt
was calculated as the slope of the transmitted energy
as a function of time, that is [(energy injection by the
hot wall - energy absorption by the hot wall) - (energy
injection by the cold wall - energy absorption by the cold
FIG. 5: (Color online) Rates of energy transmission as a function of the system size. In this simulation the viscosity coefficient was set to γ = 0.0001 and the temperatures of the
walls were set to 3 and 6. The (red) plus signs correspond
to the momentum conserving system, the (blue) stars to the
momentum dissipating system.
wall) - (energy dissipated during flow along the chain)]
per unit time. As seen in Fig. 2, both models lead to
behavior fairly well described by Fourier’s law. Repeating
this plot for different temperatures of the cold wall, we
obtain the dependence of the thermal conductivity on the
temperature. As observed in Fig. 3, κ is an increasing
function of T for both models, but it is larger for the
momentum dissipating system. On the other hand, for
a given temperature, the thermal conductivity decreases
(almost linearly) with the viscosity (see Fig. 4). Once
again, the dissipating system presents larger values of κ.
5
even when scaled in this way. In any case, we note the
very small scale of variation of the ordinate in Fig. 6.
0.0415
0.041
0.0405
J
0.04
IV.
CONCLUSIONS
0.0395
0.039
0.0385
0.038
0.0375
100
200
300
400
500
600
700
γN
FIG. 6: (Color online) Rates of energy transmission as a function of scaled system size. The temperatures of the walls were
set to 3 and 6. The (red) plus signs correspond to the momentum conserving system, the (blue) stars to the momentum
dissipating system.
Finally, we have analyzed the dependence of the rate
of energy transmission on system size, for a fixed temperature gradient and fixed viscosity. This is the crucial test
of normal vs anomalous behavior. In Fig. 5 we can clearly
see that as the size of the system increases, the rate of
energy transmission does not increase with system size.
In fact, it decreases. Therefore, in our model the thermal conductivity does not diverge with increasing system
size regardless of momentum conservation or dissipation.
One may be tempted to think of this behavior entirely
as a consequence of the energy dissipation because for
bigger systems the number of collisions necessary for the
energy to be transmitted from one end of the system
to the other increases. Hence more energy is dissipated
and less energy arrives at the cold wall. However, the
introduction of dissipation changes the dynamics more
profoundly. This can be seen in Fig. 6, where we plot the
rates of energy transmission against γN , thus taking into
account the “simple” effect of energy dissipation [see discussion preceding Eq.(13)]. The rate of energy transfer in
the momentum conserving system is essentially independent of system size in this scaled representation, but that
of the momentum dissipating system actually decreases
[1] O. Narayan,and S. Ramaswamy, Phys. Rev. Lett. 89,
200601 (2002).
[2] S. Lepri, R. Livi, and A. Politi, Phys. Rev. E 68, 067102
(2003).
[3] A. Dhar and K. Saito, Phys. Rev. E 78, 061136 (2008).
[4] S. Lepri, R. Livi, and A. Politi, Phys. Rep. 377, 1 (2003).
[5] T. Mai,and O. Narayam, Phys. Rev. E 73, 061202 (2006).
[6] T. Mai, A. Dhar, and O. Narayan, Phys. Rev. Lett. 98,
184301 (2007); L. Delfini, S Lepri, R. Livi, and A. Politi,
Phys. Rev. Lett. 100, 199401 (2008); A. Dhar and O.
Narayan, Phys. Rev. Lett. 100, 199402 (2008).
[7] S. Lepri, R. Livi, and A. Politi, Europhys. Lett. 43, 271
Heat transport in one-dimensional discrete systems
continues to be a problem of theoretical interest, uncertainty, and even controversy [1, 2, 3, 4, 5, 6]. Very
recently, experimental results in this arena have also
started to appear. In [31] the focus is the understanding
of transport of heat along colliding granular beads in a
liquid medium where the questions of interest involve the
nature of the contact regions (“liquid bridges”) between
granules. In [32] the issue is the breakdown of Fourier’s
law in nanotube thermal conductors in that the thermal
conductivity diverges with length of the nanotube. The
point is that even after many years of study the conditions that lead to the validity or violation of Fourier’s law
are not yet clear. While momentum conservation has often been featured as a condition closely associated with
the system size divergence of the 1d thermal conductivity, we have examined two dissipative granular gases, one
of which involves momentum dissipation along with energy dissipation, whereas in the other momentum is conserved. In both of these the thermal conductivity remains
finite and in fact decreases as the number of granules in
the system increases, as it must in a system where each
collision leads to energy dissipation. In the momentum
dissipating model, it decreases more rapidly than can be
accounted for by the simply heat loss mechanism of the
collisions, indicating a more profound change in the dynamics. This behavior points to the caution that must
be exercised when associating momentum conservation
with anomalous behavior in one dimension.
Acknowledgments
A.R. and I.L.D.P. acknowledge support from PronexCNPq-FAPESQ and CNPq. Acknowledgment is made to
the Donors of the American Chemical Society Petroleum
Research Fund for partial support of this research (K.L.).
(1998).
[8] B. Li, and J. Wang, Phys. Rev. Lett. 91, 044301 (2003).
[9] T. Prosen and D. K. Campbell, Phys. Rev. Lett. 84, 2857
(2000).
[10] T. Prosen and D. K. Campbell, Chaos 15, 015117 (2000).
[11] C. S. Giardiná, R. Livi, A. Politi, and M. Vassalli, Phys.
Rev. Lett 84, 2144 (2000).
[12] G. E. Duvall, R. Manvi, and S. C. Lowell, J. Appl. Phys.
40, 3771 (1969).
[13] C. Brunhuber, F. G. Mertens, and Y. Gaididei, Phys.
Rev. E 73, 016614 (2006).
[14] E. Arevalo, Y. Gaididei, and F. G. Mertens, Eur. Phys.
6
J. B 27, 63 (2002).
[15] N. V. Brilliantov, Phys. Rev. E 53, 5382 (1996); R. Ramrez, T. Pöschel, N. V. Brilliantov, and T. Schwager, Phys.
Rev. E 60, 4465 (1999); T. Pöschel and N. V. Brilliantov,
in Lecture Notes in Physics: Granular Gases, Vol. 564,
ed. by S. Luding (Springer, Berlin, Heidelberg, New York,
2000), pp. 203-212.
[16] M. Manciu, S. Sen, and A. J. Hurd, Physica D 157,
226 (2001); S. Sen et al., in AIP Conference Proceedings: Modern Challenges in Statistical Mechanics, Vol.
658, ed. by V. M. Kenkre and K. Lindenberg (AIP, New
York, 2003), pp. 357-379.
[17] S. McNamara and E. Falcon, Phys. Rev. E 71, 031302
(2005).
[18] A. Rosas, J. Buceta, and K. Lindenberg, Phys. Rev. E
68, 021303 (2003);
[19] A. Rosas, A. H. Romero, V. F. Nesterenko, and K. Lindenberg, Phys. Rev. Lett. 98, 164301 (2007).
[20] A. Rosas, A. H. Romero, V. F. Nesterenko, and K. Lindenberg, Phys. Rev. E 78, 051303 (2008).
[21] E. B. Herbold, V. F. Nesterenko, and C. Daraio, in
Shock Compression of Condensed Matter, 2005, AIP
Conf. Proc. No. 845 (AIP, New York, 2006), pp. 1523-
1526; See also cond-mat/0512367.
[22] E. B. Herbold and V. F. Nesterenko, Phys. Rev. E 75,
021304 (2007).
[23] A. Rosas, D. ben-Avraham, and K. Lindenberg, Phys.
Rev. E 71, 032301 (2005).
[24] A. Rosas and K. Lindenberg, Phys. Rev. E 68, 041304
(2003).
[25] R. Tehver, F. Toigo,J. Koplik, and J. R. Banavar, Phys.
Rev. E 57, R17 (1998).
[26] L. D. Landau and E. M. Lifshitz, Theory of Elasticity
(Addison- Wesley, Massachusetts, 1959), pp. 30.
[27] H. Hertz, J. reine u. angew. Math. 92, 156 (1881).
[28] V. F. Nesterenko, J. Appl. Mech. Tech. Phys. 5, 733
(1983).
[29] E. J. Hinch and S. Saint-Jean, Proc. R. Soc. London Ser.
A 455, 3201 (1999).
[30] A. Rosas and K. Lindenberg, Phys. Rev. E 69, 037601
(2004).
[31] J.-C. Géminard, D. Bouraya, and H. Gayvallet, Eur.
Phys. J. B 48, 500-517 (2005).
[32] C. W. Chang, D. Okawa, H. Garcia, A. Majumdar, and
A. Zettl, Phys. Rev. Lett. 101, 075903 (2008).