Papers in Physics, vol. 6, art. 060007 (2014)

Received: 14 June 2014, Accepted: 7 October 2014
Edited by: L. A. Pugnaloni
Reviewed by: R. Arévalo, School of Physical and Mathematical Sciences,

Nanyang Technological University, Singapore.
Licence: Creative Commons Attribution 3.0
DOI: http://dx.doi.org/10.4279/PIP.060007

www.papersinphysics.org

ISSN 1852-4249

Jamming transition in a two-dimensional open granular pile with rolling
resistance

C. F. M. Magalhães,1∗ A. P. F. Atman,2, 3† G. Combe,5 J. G. Moreira4‡

We present a molecular dynamics study of the jamming/unjamming transition in two-
dimensional granular piles with open boundaries. The grains are modeled by viscoelastic
forces, Coulomb friction and resistance to rolling. Two models for the rolling resistance
interaction were assessed: one considers a constant rolling friction coefficient, and the
other one a strain dependent coefficient. The piles are grown on a finite size substrate and
subsequently discharged through an orifice opened at the center of the substrate. Varying
the orifice width and taking the final height of the pile after the discharge as the order
parameter, one can devise a transition from a jammed regime (when the grain flux is always
clogged by an arch) to a catastrophic regime, in which the pile is completely destroyed
by an avalanche as large as the system size. A finite size analysis shows that there is a
finite orifice width associated with the threshold for the unjamming transition, no matter
the model used for the microscopic interactions. As expected, the value of this threshold
width increases when rolling resistance is considered, and it depends on the model used
for the rolling friction.

I. Introduction

Granular materials are ubiquitous either in nature
—desert dunes, beach sand, soil, etc.— or in indus-

∗E-mail: cfmm@unifei.edu.br
†E-mail: atman@dppg.cefetmg.br
‡E-mail: jmoreira@fisica.ufmg.br

1 Universidade Federal de Itajubá - Campus de Itabira,
Rua Irmã Ivone Drumond, 200, 35900-000 Itabira,
Brazil.

2 Departamento de F́ısica e Matemática, Centro Federal de
Educação Tecnológica de Minas Gerais, Av. Amazonas,
7675, 30510-000 Belo Horizonte, Brazil.

3 Instituto Nacional de Ciência e Tecnologia Sistemas
Complexos, 30510-000 Belo Horizonte, Brasil.

4 Universidade Federal de Minas Gerais, Caixa Postal 702,
30161-970 Belo Horizonte, Brasil.

5 UJF-Grenoble 1, Grenoble-INP, CNRS UMR 5521, 3SR
Lab. Grenoble F-38041, France.

trial processes as mineral extraction and process-
ing, or food, construction and pharmaceutical in-
dustries [1–3]. In fact, any particulate matter made
of macroscopic solid elements can be classified as
granular material. The vast phenomenology ex-
hibited by these systems combined with an incom-
plete understanding about the microscopic physical
mechanisms responsible for the macroscopic behav-
ior of these materials have motivated the increasing
interest of the physics community in the past years
[4, 5].

Although materials of this class are not sensitive
to thermal fluctuations, they can be found at gas,
liquid or solid phases [6]. The transition between
solid and liquid phases in granular matter, which
is commonly referred to as jamming/unjamming
transition, has been extensively studied from both
theoretical and experimental perspectives [4, 6–11].
Currently, a great effort is being made to under-

060007-1



Papers in Physics, vol. 6, art. 060007 (2014) / C. F. M. Magalhães et al.

stand the nature of this transition, which is still a
subject of debate [12]. The jamming/unjamming
transition is not a specific property of granular
matter, being observed in many kinds of materi-
als, such as foams [13], emulsions [14], colloids [15],
gels [16], and also in usual molecular liquids [17]
—glass transition. Liu and Nagel [9] proposed a
general phase diagram as an attempt to unify the
several approaches to study jamming/unjamming
in disordered materials. This work has motivated
several theoretical, experimental and numerical in-
vestigations, but a comprehensive understanding of
this transition is still lacking.

O’Hern et al. [18, 19] have performed numeri-
cal simulations of granular materials approaching
to jamming in two and three dimensions. They
have explored the packing fraction axis of the gen-
eral phase diagram proposed by Liu and Nagel and
have demonstrated, by means of finite-size analy-
sis, the existence of a unique critical point in which
the system jams in the thermodynamic limit. The
authors have also shown some evidence that this
point is an ordinary critical point, indicating that
the jamming transition would be a second order
phase transition. These results were corroborated
later by experiments (c.f. Majmudar et al. [10])
and by simulations (c.f. Manna and Khakhar [20]).
In Ref. [20], the authors have revealed evidence of
self-organized criticality (SOC) by measuring the
internal avalanches resulting from the opening of
an orifice at the bottom of granular piles.

Experimental investigation of the jamming tran-
sition in granular materials under gravitational
field has been conducted in a variety of ways, ad-
dressing the role of many parameters, like the grain
shape, the friction coefficient, and the system ge-
ometry. However, a common feature of all these ap-
proaches is the analysis of the granular flow through
bottlenecks.

The jamming of three-dimensional piles seems to
be settled after the work of Zuriguel et al. [21].
They have demonstrated experimentally, for piles
composed of different kinds of grains, the diver-
gence on the mean internal avalanche size. It means
that, as the outlet size approaches a critical value,
the internal avalanche increases without limit and
a permanent flow is established. This critical outlet
is insensitive to the density, stiffness and roughness
of the grains, but shows a significant dependence on
the grain shape. For spherical grains, a critical out-

let width wc ∼ 4.94d was obtained, for cylindrical
grains wc ∼ 5.03d and for rice grains wc ∼ 6.15d.

Nevertheless, the jamming transition in two-
dimensional piles is still a question under debate.
In order to address it, To et al. [22] have carried out
experiments using two-dimensional hoppers in or-
der to find a critical outlet size for jamming events.
The jamming probability J has presented a rapid
decay from J = 1 to J = 0 close to the aperture
width w ∼ 3.8d, signaling a possible phase transi-
tion. The authors discussed, based on a restricted
random walk model, the connection between jam-
ming and the arch formation mechanism. Never-
theless, the point needs further investigation, espe-
cially at the limit of high hopper angle. There exist
several works [23–25] focused on the mechanisms of
arch formation, but none had explored its relation
to jamming probability.

Janda et al. [26] have made some progress
by simulating discharges of two-dimensional silos.
They have improved the definition of jamming so
that the internal avalanche size is considered. This
modification addresses the extremely long relax-
ation times associated with jamming. Within the
framework of a probabilistic theory concerning the
arch formation, the authors have tested two hy-
pothesis for the internal avalanche behavior: a
functional form that predicts a divergence in the
mean size, and a functional form where the mean
size exists for all values of the orifice width w. Since
the latter one is more compatible with the arch for-
mation model, the authors claimed that “no critical
opening size exists beyond which there is not jam-
ming” [26].

The results mentioned so far are related to fully
confined systems. Recently, a simulation study
on discharges of granular piles with open bound-
aries [27] provided new insights on the problem.
The piles are composed by homogeneous disks in-
teracting via elastic and frictional forces. Using
finite size analysis, the authors have shown that
a catastrophic regime, in which the pile is com-
pletely destroyed by the opening of the orifice, is
well defined. At the limit of infinitely large sys-
tems, the catastrophic regime coincides with the
unjammed phase, since it implies a divergent in-
ternal avalanche. Hence, the results indicate the
existence of the jamming transition. It is impor-
tant to note, however, that the pile geometry could
probably play a role in the causes for this distinct

060007-2



Papers in Physics, vol. 6, art. 060007 (2014) / C. F. M. Magalhães et al.

behavior, due to the absence of the Janssen effect,
but further investigations are necessary to confirm
this point.

In the present work, the investigation of jamming
in 2D open systems is extended in order to con-
sider a rolling resistance term in the grain interac-
tion model, following the prescription adopted by
Chevalier et al. [28]. The main objective of this
study is the verification of rolling resistance influ-
ence on jamming. Many factors contribute to the
appearance of rolling friction, including microslid-
ing, plastic deformation, surface adhesion, grain
shape, etc., but mainly, it is due to the contact de-
formation [29]. Here, it will be taken into account
only the effect due to the contact deformation by
implementing the micromechanical model proposed
by Jiang et al. [30]. The rolling friction produces
a resistance to roll which, among other effects, is
responsible for granting more stability to granular
piles [31], and for the occurrence of different types
of failure modes in granular matter [32,33]. Rolling
friction was also used to model a system of polygo-
nal grains by making a correspondence between the
rolling stiffness and the number of sides of the poly-
gon [34]. It was demonstrated that it is an essential
ingredient to reproduce experimental compression
tests in mixtures of two-dimensional circular and
rectangular grains [35]. These facts suggest that
jamming could be affected by rolling friction. Nev-
ertheless, most studies on granular materials based
on computer simulations do not deal with it.

The paper is organized as follows: after a review
in the Introduction, the next section is concerned
with the methodology. Then, we present the re-
sults and a brief discussion. Finally, the last section
gathers the conclusions and some perspectives.

II. Methods

The jamming transition is assessed by means of
discharges of granular piles, simulated using the
molecular dynamics method [36] with the Velocity-
Verlet algorithm [37]. In a few words, the molecu-
lar dynamics consists in integrating numerically the
equations of motion that governs the system dy-
namics. The system is constituted by N free grains
governed by Newton’s second law and by a finite
and horizontally aligned substrate made of fixed
grains. The free grains, which will form the pile, are

homogeneous bi-dimensional disks that are free to
translate and rotate around their center, and whose
radii are uniformly distributed around an average
value d, a small polydispersity of 5% was imposed
in order to avoid crystallization effects. All spatial
quantities will be expressed in terms of d. Since the
grains have all the same mass density, their masses
are proportional to their respective areas. Normal-
ized by the mass of the heaviest grain, the masses
are given by mi = d

2
i /d

2
max, where dmax stands for

the diameter of the largest grain. The finite sub-
strate of length L is composed by fixed grains of
the same kind but with a smaller and fixed diame-
ter ds = 0.1d. These grains are aligned horizontally
and are equally spaced in order to form a grid with-
out gaps.

The grains are subject to a uniform gravita-
tional field orthogonal to the substrate line, and to
short range binary interactions. Besides the visco-
elastic and coulomb friction interactions used in
past works [25, 27, 38–42], two grains in contact are
also subject to a rolling resistance moment due to
the finite contact length lc. The rolling resistance is
introduced through a micromechanical model of the
contact line between two grains [30]. This model
treats the contact as an object formed by a set of
springs and dashpots connecting the borders of the
two grains. As one grain rolls over the other, the
springs in one side of the contact line contract while
the springs in the opposite side stretch. This config-
uration generates an unbalanced force distribution
and a consequent moment with respect to the grain
center (see Fig. 1). This moment grows linearly
as the grain rolls, until the rolling displacement δr
reaches some threshold value at which the springs
located near one end of the contact line break up
and new ones emerge at the other end. At that
time, the moment saturates at some value that de-
pends on the properties of the grain. The rolling
displacement is defined by δr =

∑
(ωi − ωj)∆t,

where ωi refers to the angular velocity of grain i,
∆t is molecular dynamics time step, and the sum
runs over time during the whole existence of con-
tact. Based on these assumptions, the authors have
derived an analytic expression for the rolling resis-
tance moment as a function of rolling displacement.
They have also proposed a simplified version - the
one used in this study - in order to improve numer-
ical computations:

060007-3



Papers in Physics, vol. 6, art. 060007 (2014) / C. F. M. Magalhães et al.

τrr =



−krδr, kr|δr| ≤ µrfel
−
δr
|δr|

µrfel, kr|δr| > µrfel ,
(1)

where kr is the rolling stiffness, µr is the coefficient
of rolling resistance, and fel is the compressive elas-
tic force, normal to contact line. As can be noted,
the expression for the rolling resistance momentum
possesses a striking resemblance with the Coulomb
static friction force, and as a matter of fact, the
rolling resistance interaction is implemented in the
molecular dynamics algorithm in the same way as
the static friction. The model predicts that the
rolling stiffness and the coefficient of rolling resis-
tance are related to the contact length lc by the
equations kr = knl

2
c and µr = µlc, where kn and

µ are respectively the normal elastic constant and
friction coefficient. The values used for these pa-
rameters were the same as in Ref. [27], µ = 0.5 and
kn = 1000 in normalized unities (see [40] for further
details). Two cases were investigated: systems in
which the rolling resistance parameters kr and µr
vary according to the above-mentioned equations,
and systems with fixed rolling resistance parame-
ters, assuming that all contacts have the same de-
formation value.

The simulation procedure consists of two steps:
(1) the formation of the granular pile with open
boundaries, by deposing grains from rest, under
gravity, over the substrate until an stationary state
is reached; (2) the discharge itself, which consists in
opening an orifice of a given width at the center of
the substrate. In the first step, the initial positions
of the free grains are randomly sorted along a hori-
zontal line located at height L from the substrate -
the releasing height is equal to the substrate length.
To avoid initial overlapping of grains, a 50% filling
ratio was imposed to each line of grains released,
and the time interval between successive rows is
the inverse of the frequency f. Each row was re-
leased after the predecessor had fallen a distance
equivalent to the maximum grain diameter. This
deposition protocol mimics a dense rain of grains.
During deposition, the grains may leave the sys-
tem through the lateral boundaries, so that the to-
tal number of grains in the pile fluctuates as the
process evolves. The release of grains ceases when
the number of grains in the pile reaches a station-
ary value. The deposition phase ends only when a

Figure 1: In panel (a), a perfectly rigid grain rolling
over another one is exhibited. The contact force is
concentrated on one point and is aligned through
the grain center, thus, not generating momentum
with respect to the center. Panel (b) shows the
same situation but with deformable grains. Now,
the contact force is non-uniformly distributed over
a segment of length lc (contact length). The re-

sultant force ~FR,c is dislocated from the center line
and an opposing moment with respect to the center
appears.

mechanical equilibrium state is attained, and the
configuration is recorded for later analysis. The
equilibrium state must satisfy the following crite-
ria: mechanical stability, absence of slipping con-
tacts, vertical and horizontal force balance, and
vanishingly small kinetic energy [43]. In the sec-
ond step, the configuration recorded is loaded and
an outlet of width w is opened at the center of the
substrate, allowing the grains to flow through it.
As the grains pass through the outlet, they are re-
moved from the system, in order to improve compu-
tational efficiency. The simulation runs until a new
equilibrium state is reached, which happens either
due to the formation of an arch above the outlet or
after the pile has been completely discharged. The
remaining pile configuration is then recorded.

The average height h of the resulting pile after
discharge, measured from the center of the sub-
strate, was taken as an order parameter to distin-
guish between the two regimes. If h does not change
significantly with respect to the original height, it
means that an arch does readily clog the flux af-
ter the orifice opening, but if h = 0, it means that
the pile has collapsed, and a jammed state was not
attained for the corresponding orifice diameter and
substrate length. As the orifice width w approaches

060007-4



Papers in Physics, vol. 6, art. 060007 (2014) / C. F. M. Magalhães et al.

Figure 2: Diagram of a possible equilibrium con-
figuration if rolling resistance is considered. The
grain on top has only one contact and, even so, it
can sustain a stable position. The maximum value
of the angle θ depends on the friction coefficient
and on the rolling resistance parameters.

a certain threshold value wt, the system suffers a
transition from jammed to unjammed state. This
threshold is defined as the orifice width for which
the height fluctuations is maximum. The connec-
tion to the jamming transition occurs when the
substrate length diverges since the collapse of an
infinite substrate pile implies a continuous flowing
state. Then, the jamming transition can be charac-
terized by a critical aperture width, as defined by
the following expression

wc = lim
L→∞

wt(L) . (2)

III. Results and Discussion

Two different models of rolling resistance interac-
tion were tested in the simulations: the original
model proposed by Jiang et al. [30] mentioned
earlier, with a rolling constant and a coefficient
of rolling resistance that depends on the contact
length (kr = knl

2
c and µr = µlc) and a simpler

derived version, in which these parameters are con-
stant over time assuming that lc = 0.05 d for all
contacts. This fixed lc model assumes a mean con-
tact length equivalent to 5% of the average grain
diameter, a scenario which could be associated to a
system composed by polygonal grains, for example,
a extreme rolling resistance regime.

For either models, the numerical simulations de-
scribed in the last section were carried out for

Figure 3: Images of the pile equilibrium states for
the three grain interaction models before the dis-
charge step. The piles were grown over a L = 100d
substrate and are representative samples from each
type of pile. As shown in the figure, the top image
represents a pile composed by grains with a fixed
kr rolling resistance interaction, the image in the
middle is a pile of grains with a kr ∼ l2c rolling re-
sistance interaction, and the bottom image is a pile
of grains without rolling resistance.

various system sizes and orifice widths. Figure
3 shows equilibrium configurations of typical piles
with L = 100d for the two tested rolling resistance
models and for the absence of rolling resistance
case. It can be seen that the inclusion of the rolling
resistance term modifies significantly the macro-
scopic features of the pile. It provides more stabil-
ity to the structure, which is reflected by a steeper
free surface, a fact also observed elsewhere [31]. In-
deed, it should be noted that the rolling resistance
makes possible some otherwise very unstable two-
grain configuration, as exemplified in Fig. 2.

The behavior of the order parameter h as func-
tion of w is presented in Fig. 4 for the three cases.

060007-5



Papers in Physics, vol. 6, art. 060007 (2014) / C. F. M. Magalhães et al.

0 2 4 6 8 10 12 14
w (d)

0

0.2

0.4

0.6

0.8

h
 (

n
o

rm
a
li

z
e
d

 u
n

it
)

k
r
 = const.

k
r
 ~ l

c

2

k
r
 = 0

Figure 4: Order parameter h as a function of the
orifice width w for all types of pile. The graphs were
generated from the simulation data of L = 100d
piles. While the symbols represent the parameter
h itself, the lines indicate the corresponding fluc-
tuations. The thick line is related to the fixed kr
curve, the medium thickness line to the kr ∼ lc2
curve, and the dashed line to kr = 0 curve.

Note that the transition region translates to the
right as the rolling resistance becomes more impor-
tant, which is an expected result since the increas-
ing of stability allows the formation of larger arches.
Figure 5 exhibits the dependence of wt (fluctuation
maximum) on the system size for the three mod-
els considered and, again, the curves are dislocated
along the wt axis. Nevertheless, they all share the
same functional aspect, which is an evidence that
the rolling resistance does not change the nature
of the jamming transition, only the critical value
of the threshold width. Despite the tendency to a
heavy tail distribution observed for large L and w,
in all scenarios, the data suggest the existence of
the jamming transition. It means that the features
described in Ref. [27] to characterize the transi-
tion were also observed in both scenarios tested
for the rolling resistance. The fitting values were
wc = (5.3 ± 0.1)d for contact length dependent
kr model and wc = (7.6 ± 0.2)d for the fixed kr
one, while for the model without rolling resistance,
wc = (5.0 ± 0.1)d [27]. These results indicate that
the rolling resistance only slightly changes the crit-
ical width when included in a system of disks, ex-
pressed by the contact length model. But, in other

0 0.01 0.02 0.03 0.04

1 / L (normalized units)

0

1

2

3

4

5

6

7

8

w
t 
(u

n
it

s 
o
f 

d
)

k
r
 = 0

k
r
 ~ l

c

2

Fixed k
r

Figure 5: Graphs of the threshold orifice width wt
as a function of the reciprocal system size 1/L for
models with and without rolling resistance. As the
caption indicates, the symbols represent the data
obtained from numerical simulations and the lines
are fitting curves, which provide the respective val-
ues of wc.

approaches, as in the fixed length model, it can al-
ter significantly the critical aperture width.

The probability density functions (PDF) of h for
the absent rolling friction and for the fixed rolling
friction models are exhibited in Fig. 6. For each
model, the PDFs are obtained for three characteris-
tic values of w: w = 1.0d , w ∼ wt, and w > wt. It
can be noted that the PDFs have an approximately
Gaussian peak around the initial height for all val-
ues of w, meaning that there is always a certain
amount of samples that are simply not disturbed
by the outlet. Apart from these samples, the great
majority is completely discharged for w > wt. This
result is in consonance with that obtained earlier
in a different context [25]. In that study, it was
shown that for large w, all blocked events occurred
after only few grains had passed through the outlet.
These facts indicate that the initial conditions may
play an important role in jamming experiments.
However, this issue needs further investigation.

IV. Conclusions

Evidence of the jamming transition was observed
in molecular dynamics simulations of open granu-
lar piles with rolling resistance, in consonance with

060007-6



Papers in Physics, vol. 6, art. 060007 (2014) / C. F. M. Magalhães et al.

0 0.2 0.4 0.6
h (normalized units)

0.01

0.1

1

10

100

h
 p

d
f

w = 1.0 d
w = 4.0 d
w = 6.0 d
Gaussian Fit

0 0.2 0.4 0.6 0.8 1
h (normalized units)

0.01

0.1

1

10

100

h
 p

d
f

w = 1.0 d
w = 7.0 d
w = 9.0 d
Gaussian Fit

No Rolling Friction

Fixed RF Parameters

Figure 6: Probability density functions of the order
parameter h for different values of w in the case of
absent rolling friction grains (top) and fixed rolling
friction grains (bottom).

the results found in granular piles without this in-
teraction. This result strengthens the expectation
that the transition exists in real granular piles and
is probably affected by the system geometry. We
observe that when there was rolling resistance, the
piles built were more stable, denoted for the large
mean height of the samples, and also more ro-
bust against perturbations, since the critical aper-
ture width increased. In future works, we plan to
present a detailed study of the arching statistics for
the two approaches considered here.

Acknowledgements - We are grateful for CNPq,
FAPEMIG Brazilian funding agencies. APFA and
GC thanks to CEFET-MG by the international in-
terchange which made this interaction possible.

[1] S J Antony, W Hoyle, Y Ding, Granu-
lar materials: Fundamentals and applications,
The Royal Society of Chemistry, Cambridge
(2004).

[2] J Duran, Sands, powders and grains, Springer,
Berlin (1997).

[3] T Halsey, A Mehta, Challenges in granular
physics, World Scientific Publishing, New Jer-
sey (2002).

[4] A Yu , K Dong , R Yang, S Luding, Pow-
ders and grains 2013: Proceedings of the 7th
international conference on micromechanics of
granular media, AIP Series. Vol. 1542, Sydney,
Australia (2013).

[5] GdR MiDi, Eur. Phys. J. E 14, 341 (2004).

[6] H M Jaeger, S R Nagel, Granular solids, liq-
uids, and gases, Rev. Mod. Phys. 68, 1259
(1996).

[7] P Cixous, E Kolb, N Gaudouen, J-C Charme,
Jamming and unjamming by penetration of
a cylindrical intruder inside a 2 dimensional
dense and disordered granular medium, In:
Powders and grains 2009, Proceedings of the
6th international conference on micromechan-
ics of granular media 1145, 539 (2009).

[8] A P F Atman, P Claudin, G Combe, R Mari,
Mechanical response of an inclined frictional
granular layer approaching unjamming, Euro-
phys. Lett. 101, 44006 (2013).

[9] A J Liu, S R Nagel, Jamming is not just cool
any more, Nature 396, 21 (1998).

[10] T S Majmudar, M Sperl, S Luding, R P
Behringer, Jamming transition in granular
systems, Phys. Rev. Lett. 98, 058001 (2007).

[11] C Mankoc, A Janda, R Arévalo, J M Pastor,
I Zuriguel, A Garcimart́ın and D Maza. The
flow rate of granular materials through an ori-
fice, Granul. Matter 9, 407 (2007).

[12] A P F Atman, P Claudin, G Combe, G H
B Martins, Mechanical properties of inclined
frictional granular layers, Granul. Matter 16,
1 (2014).

[13] G Katgert, M van Hecke, Jamming and geome-
try of two-dimensional foams, Europhys. Lett.
92, 34002 (2010).

[14] N D Denkov, S Tcholakova, K Golemanov,
A Lips, Jamming in sheared foams and emul-
sions, explained by critical instability of the
films between neighboring bubbles and drops,
Phys. Rev. Lett. 103, 118302 (2009).

060007-7



Papers in Physics, vol. 6, art. 060007 (2014) / C. F. M. Magalhães et al.

[15] A Kumar, J Wu, Structural and dynamic prop-
erties of colloids near jamming transition, Col-
loid. Surf. A 247, 145151 (2004).

[16] A Fluerasu, A Moussaid, A Madsen, A
Schofield, Slow dynamics and aging in col-
loidal gels studied by x-ray photon correlation
spectroscopy, Phys. Rev. E 76, 010401 (2007).

[17] B Duplantier, T C Halsey, V Rivasseau,
Glasses and grains: Poincaré Seminar 2009,
Springer, Basel (2011).

[18] C S O’Hern, S A Langer, A J Liu, S R
Nagel, Random packings of frictionless parti-
cles, Phys. Rev. Lett. 88, 075507 (2002).

[19] C O’Hern, L E Silbert, A J Liu, S R Nagel,
Jamming at zero temperature and zero applied
stress: The epitome of disorder, Phys. Rev. E
68, 011306 (2003).

[20] S S Manna, D V Khakhar, Internal avalanches
in a granular medium, Phys. Rev. E 58, R6935
(1998).

[21] I Zuriguel, A Garcimart́ın, D Maza, L A Pug-
naloni, J M Pastor, Jamming during the dis-
charge of granular matter from a silo, Phys.
Rev. E 71, 051303 (2005).

[22] K To, P-Y Lai, H K Pak, Jamming of granular
flow in a two-dimensional hopper, Phys. Rev.
Lett. 86, 71 (2001).

[23] A Garcimart́ın, I Zuriguel, L A Pugnaloni,
A Janda, Shape of jamming arches in two-
dimensional deposits of granular materials,
Phys. Rev. E 82, 031306 (2010).

[24] A Drescher, A J Waters, C A Rhoades, Arch-
ing in hoppers .2. Arching theories and critical
outlet size, Powder Technol. 84, 177 (1995).

[25] C F M Magalhães, A P F Atman, J G Moreira,
Segregation in arch formation, Eur. J. Phys. E
35, 38 (2012).

[26] A Janda, I Zuriguel, A Garcimart́ın, L A Pug-
naloni, D Maza, Jamming and critical outlet
size in the discharge of a two-dimensional silo,
Europhys. Lett. 84, 44002 (2008).

[27] C F M Magalhães, J G Moreira, A P F Atman,
Catastrophic regime in the discharge of a gran-
ular pile, Phys. Rev. E 82, 051303 (2010).

[28] B Chevalier, G Combe, P Villard, Experi-
mental and discrete element modeling stud-
ies of the trapdoor problem: Influence of the
macro-mechanical frictional parameters, Acta
Geotech. 7, 15 (2012).

[29] J Ai, J-F Chen, J M Rotter, J Y Ooi, Assess-
ment of rolling resistance models in discrete el-
ement simulations, Powder Technol. 206, 269
(2011).

[30] M J Jiang, H-S Yu, D Harris, A novel dis-
crete model for granular material incorporat-
ing rolling resistance, Comput. Geotech. 32,
340357 (2005).

[31] Y C Zhou, B D Wright, R Y Yang, B H Xu,
A B Yu, Rolling friction in the dynamic sim-
ulation of sandpile formation, Physica A 269,
536 (1999).

[32] K Iwashita, M Oda, Rolling resistance at con-
tacts in simulation of shear band development
by DEM, J. Eng. Mech. 124, 285292 (1998).

[33] X Li, X Chu, Y T Feng, A discrete particle
model and numerical modeling of the failure
modes of granular materials, Eng. Computa-
tion. 22, 894 (2005).

[34] N Estrada, E Azéma, F Radjai, A Taboada,
Identification of rolling resistance as a shape
parameter in sheared granular media, Phys.
Rev. E 84, 011306 (2011).

[35] E-M Charalampidou, G Combe, G Viggiani,
J Lanier, Mechanical behavior of mixtures of
circular and rectangular 2D particles, In: Pow-
ders and grains 2009: Proceedings of the 6th
international conference on micromechanics of
granular media, AIP Conf. Proc., Vol. 1145,
Pag. 821, (2009).

[36] D C Rapaport, The art of molecular dynamics
simulation, Cambridge University Press, Cam-
bridge (2004).

[37] W C Swope, H C Andersen, P H Berens, K R
Wilson, Computer simulation method for the

060007-8



Papers in Physics, vol. 6, art. 060007 (2014) / C. F. M. Magalhães et al.

calculation of equilibrium constants for the for-
mation of physical clusters of molecules: Ap-
plication to small water clusters, J. Chem.
Phys. 76, 637 (1982).

[38] C Goldenberg, A P F Atman, P Claudin,
G Combe, I Goldhirsch, Scale separation in
granular packings: Stress plateaus and fluctu-
ations, Phys. Rev. Lett. 96, 168001 (2006).

[39] S F Pinto, A P F Atman, M S Couto, S
G Alves, A T Bernardes, H F V Resende,
E C Souza, Granular fingers on jammed sys-
tems: New fluidlike patterns arising in grain-
grain invasion experiments, Phys. Rev. Lett.
99, 068001 (2007).

[40] A P F Atman, P Claudin, G Combe, Depar-
ture from elasticity in granular layers: Investi-
gation of a crossover overload force, Comput.
Phys. Commun. 180, 612 (2009).

[41] P A Cundall, O D L Strack, A discrete numer-
ical model for granular assemblies, Geotech-
nique 29, 47 (1979).

[42] R D Mindlin, Compliance of elastic bodies in
contact, J. Appl. Mech. 71, 259 (1949).

[43] A P F Atman, P Brunet, J Geng, G Rey-
dellet, G Combe, P Claudin, R P Behringer,
E Clement, Sensitivity of the stress response
function to packing preparation, J. Phys.:
Cond. Matter 17, S2391 (2005).

[44] H Hertz, On the contact of elastic solids, J.
Reine Angew. Math. 92, 156 (1881).

[45] M P Allen, D J Tildesley, Computer sim-
ulation of liquids, Clarendon Press, Oxford
(1987).

[46] O O’Sullivan, Computing quaternions, In: The
art of numerical manipulation, Eds. A Q Rista,
M Nadola, Pag. 132, North Holland, Amster-
dam (2003).

060007-9