Table Of ContentGeometric origin of mechanical properties of granular materials
Jean-No¨el Roux
∗
Laboratoire Central des Ponts et Chauss´ees, 58, boulevard Lef`ebvre,75732 Paris cedex 15, France
Modelgranularassemblies,inwhichgrainsareassumedrigidandfrictionless,atequilibriumunder
some prescribed external load, are shown to possess, under generic conditions, several remarkable
mechanical properties, related to isostaticity and potential energy minimization. Isostaticity – the
9 uniquenessofthecontactforces, oncethelistofcontactsisknown–isestablished inaquitegeneral
0 context, and the important distinction between isostatic problems under given external loads and
0 isostatic(rigid)structures ispresented. Completerigidityisonlyguaranteed,onstabilitygrounds,in
2 thecaseofsphericalcohesionlessgrains. Otherwise,thenetworkofcontactsmightdeformelastically
n inresponsetosmallloadincrements,eventhoughgrainsare perfectlyrigid. Ingeneral, onegetsan
a upperboundon thecontact coordination number. Theapproximation of small displacements, that
J isintroducedanddiscussed,allowstodrawanalogieswithothermodelsystemsstudiedinstatistical
0 mechanics,suchasminimumpathsonalattice. Italsoentailstheuniquenessoftheequilibriumstate
2 (thelist ofcontactsitselfisgeometrically determined)forcohesionless grains, andthustheabsence
of plastic dissipation in rearrangements of the network of contacts. Plasticity and hysteresis are
] related to thelack of such uniqueness,which can be traced back, apart from intergranular friction,
t tonon-reversiblerearrangementsofsmallbutfiniteextent,inwhichthesystemjumpsbetweentwo
f
o distinctpotentialenergyminimain configuration space,or toboundedtensileforces, derivingfrom
s a non-convex potential, in the contacts. Properties of response functions to load increments are
.
t discussed. On the basis of past numerical studies, it is argued that, provided the approximation of
a
smalldisplacementsisvalid,displacementsduetotherearrangementsoftherigidgrainsinresponse
m
tosmallloadincrements,onceaveragedonthemacroscopicscale,aresolutionstoellipticboundary
- valueproblems (similar to the Stokesproblem for viscous incompressible flow).
d
n
PACSnumbers: 46.10.+z,05.40.+j,83.70.Fn
o
c
[
I. INTRODUCTION some grains or sets of grains carry but vanishing ef-
2 forts (“arching effect”). The histogram of contact forces
v spans a wide range. These phenomena have been ex-
A. Motivations
6
perimentally observed thanks to techniques like photoe-
4
lastic stress visualization [4, 5] and carbon paper print
2
A large research effort, both in the statistical physics analysis [6, 7]. They have also been studied in numeri-
1
0 and the mechanics and civil engineering communities, is cal simulations [8, 9], and some attempts of theoretical
0 currently being devoted to granular materials, aiming in descriptions have been proposed [10]. Such peculiar as-
0 particular at a better understanding of the relationships pects of granular systems render more difficult the ref-
/ between grain-level micromechanics (intergranular con- erence to existing models from other fields. Indeed, a
t
a tact laws) and macroscopic behaviours (global equilib- recent trend in the physics literature on static granular
m
rium conditions, constitutive relations) [1, 2, 3] systems [11, 12, 13, 14] insists on their difference with
- This aim –the traditional program of Statistical Me- ordinary, elastic solids, and suggests, instead of resort-
d
chanics – is far from fully achieved in dense granular ing to macroscopic displacement or strain variables, to
n
o systems near equilibrium, for one is facing at least two searchfor direct relations betwen the components of the
c fundamental difficulties. stress tensor.
:
v Firstly,thenon-smoothcharacterofcontactlaws,that The second basic difficulty stems from the incomplete
i involveunilateralityand, possibly,dry friction, is a com- knowledge of the mechanical properties of granular sys-
X
mon feature of granular assemblies that endows them tems,especiallythoserulingthedynamics. Whenagran-
r with a high level of disorder and a high sensitivity to ularsample is submitted to some prescribedexternalac-
a
perturbations. Tiny motions might significantly affect tions that are sufficiently slowly changing in time, its
the way forces are transmitted, since contacts between evolution is customarily described as an ordered set of
neighbouring grains might open or close (and the slid- equilibriumstatesthataresuccessivelyreached,withlit-
ing or non-sliding status of closed ones might change). tleornodependenceonphysicaltime. Thephysicalpro-
Hencethecharacteristicallyheterogeneousaspectofforce cessesbywhichkineticenergyisdissipatedare,however,
transportindense granulates: largeforcesarecarriedby mostoftensomewhatmysteriousorpoorlycharacterised.
a network of preferred paths (the “force chains”) while Theyare,intheframeworkofthequasi-static description
wehavejustmentioned,implicitlyregardedasirrelevant.
Onemightwishto assessthevalidity ofsuchanassump-
tion. Numerical simulations, that have to adopt some
∗Electronicaddress: [email protected] rule to move the grains, could in principle allow useful
2
investigationsoftheinfluenceofthedynamics. However, andconceptualimportance,sinceitallows,inparticular,
in view of the practical difficulty to obtain representa- ananalogywith problems of scalartransporton discrete
tive configurations close enough to equilibrium within a networks, as explained in section V.
reasonable computation time, they sometimes resort to
Once those essential ingredients made available, the
non-physicalparameters,andpickupthedynamicalrule
second part of the paper (sections VI to IX) establishes
among the restricted range of those that allow tractable
the main results and discusses their consequences, with
calculations.
referencetoprevioustheoreticalandnumericalwork,and
This paperaddressesboththosebasicconcerns,inthe
toknownaspectsofthemechanicalbehaviourofgranular
following way. Simplifying assumptions are introduced
materials.
(we consider,e.g., rigidfrictionlessgrains),thus restrict-
ing our attention to a certain class of model systems, Section VI is devoted to the generic isostaticity prop-
that are however argued to exhibit the same qualitative erty of equilibrium states in systems of rigid grains that
behavioursasmorerealisticones. Thosesystemsaresuit- mayonlyexertnormalcontactforcesononeanother. We
able candidates to test, most easily by numerical means, then prove and discuss (section VII) the uniqueness of
some recently proposed models and speculations, at the the equilibrium state in cohesionless systems within the
expenseofratherextensivenumericalcomputations. The
approximation of small displacements, and compare the
purpose of the present article is not, however,to present
determinatin of equilibrium states of such systems with
new results of numericalsimulations. We shall state and
other mechanical or scalar transport problems. Section
establish, rather, with a fair level of generality, some
VIII introduces the additional requirement of stability,
basic properties of such systems, and study their qual-
outside the approximation, which is dealt with, in the
itative consequences in terms of macroscopicmechanical
absenceoffriction,intermsofpotentialenergyminimiza-
behaviour. This analysis will shed some light on some
tion. In some restricted models, this allows to conclude
analogies and differences with other previously studied
to the isostaticity of the structure, a stronger property
problemsinstatisticalmechanics,suchasdirected‘poly-
than mere isostaicity of the problem under a given load.
mers’ in random environments and percolation models.
It is then possible to discuss the possible origins of plas-
Itwillalso,alongwiththeexploitationofpastnumerical
ticityinsystemsoffrictionlessgrainsandtheformofthe
resultsonasimplifiedmodel[9,15,16,17,18],allowusto
mechanicalresponsetosmallloadincrements. Thepaper
investigate the possible origins of some macroscopic fea-
endswithconcludingremarks(X)ontheroleofdisplace-
turesofgranularmechanics,thatareclassicallymodelled
ments and strains in granular materials and suggestions
withelastoplasticconstitutivelaws[3,19],andtodiscuss
for future research.
other recently proposed approaches [11, 12, 13, 14].
We will show that mechanics is to a large extent de-
termined by geometrical aspects (steric exclusion), thus
partiallyansweringconcernsaboutthe roleofdynamical
parameters. Finally we will discuss the status of dis-
placement and strain variables in quasi-static assemblies
II. BASIC DEFINITIONS AND PROPERTIES.
of rigid grains, and give perpectives for future investiga-
tions.
We are interested in the modelling of large packings
of solid bodies (grains), in equilibrium under some pre-
B. Synopsis.
scribed external forces. Grains are assumed to interact
via pointforcesmutuallyexertedontheirsurfaces,which
The paper is composed of two main parts. means that the distribution of stress on their areas of
First, sections II to V introduce useful definitions and contactor of influence caneffectively be viewed as local-
state basic properties that are necessary for the deriva- ized at a point, on the scale of the whole grain. Apart
tion of the main results. from this reservation, that excludes flat or conforming
Thus, section II presents useful definitions and me- surfaces [56], grains might have arbitrary shapes, and
chanical properties of static granular systems, i.e., col- our considerations apply to spatial dimension d equal to
lections of rigid bodies essentially interacting via point 2 or 3, although most examples will be taken with two-
forcesmutuallyexertedontheirsurfaces. Thosenotions, dimensionalsystemsofdiscs. Notethatwedonotrequire
that include the theorem of virtual power, generalized interactinggrainstotouchoneanotheratthisstage. We
forces and velocities for collective degrees of freedom, mostly restrict our attention here to frictionless bodies,
and the degree of indeterminacy of forces and of veloc- i.e., suchthatcontactforcesarenormaltothe grainsur-
ities, are not always familiar in the condensed matter faces. This might look like a severe limitation, but we
physics community. Section III introduces the potential shall argue that such simplified systems do possess the
energyminimizationproblemsforvarioussimplefriction- generic properties of granular media. We shall also as-
less contact laws. Section IV defines the approximation sume, unless otherwise specified, that the grains behave
ofsmalldisplacements,amodellingstepofbothtechnical as rigid undeformable objects.
3
A. System, external forces In order to enforce constraints of type 1, some external
efforts have to be exerted on the concerned bodies. On
requiring the power of such efforts to be balanced by
We consider a set of n grains, labelled with indices i,
that of internal forces (Fint) , one identifies the
wai‘tchen1te≤r’,iw≤hinch. mInigehatche.go.f,tchoeimnciwdee awribthitriatrsilcyenctheorosoef generalized force conjugatµe to1≤λ1µ≤aNsf
mass. In the case of spherical grains it is of course con-
venienttotakethegeometricalcenterofthesphere. The Q1 =− F(iin,αt)A(i,α). (3)
vg(dee-tlodhciemirtiewensistih(oΩntha)el)dv′e-dloimc,imteineasskiooenfuathplo(wstehitecheknditn′ee=rms,da((tVdi2−ci1)d)1)e≤gai≤rnenge,ustlaoor-f Wejustusedthepower(tio,Xαfi)∈nId1generalizedforces: thisisa
i 1 i n
freedom of the w≤h≤ole system, thus labelled by couples of manifestationoftheduality betweenforcesanddisplace-
mentsorvelocities,whichwillbe repeatedlyexploitedin
indices(i,α),with1 i nand1 α d+d = d(d+1).
We denote as I the s≤et o≤f such cou≤ples.≤If α>′ d, Vi,2α is ttherenaselqfuoercl.es,Tihse, bNyfc-doinmsternuscitoinoanl,vteoctboerrsepgaacredeFd aosf tehxe-
now a notation for Ω . Boundary conditions are of-
tenenforcedbyprescrii,αb−indgthemotion,ortheabsenceof dualspace,inthe ordinarysenseof linearalgebra,ofthe
N -dimensionalspace ofkinematicdegreesoffreedom.
motion, of walls. Those might be regarded as solid bod- f V
In general, it should be appreciated that the appro-
ies,or particular‘grains’themselves. Inthe followingwe
priatemathematicaldescriptionofconfigurationspaceis
shall sometimes write down large ‘velocity vectors’ that
not IRNf with its Euclidean structure, but, due to rota-
gatherallN kinematicdegreesoffreedomofthesystem,
f tional degrees of freedom, an N -dimensional manifold,
then denoted, with a single index, as (v ) . f
It mightalsobe convenientto keepsoµm1e≤gµr≤aNinfcoordi- doninwatheisc.h (xaµn)d1≤µ≤arNefriesspaescettivoefly(ltohceatl)ancugernvtiliannedarcoctoaonr--
natesfixed(thuschoosingoneparticularGalileanframe), V F
gent vector space at a given point, and depend on that
i.e., toimpose,forallcouplesi,αbelongingtosomesub-
point. Thus the definition of ‘constant velocities’, or of
setI ofI, V =0. Indices µarethenrenumbered,and
0 i,α ‘constant forces’ requires some care. However, these dif-
N is reduced accordingly,to label and to count the free
f ficulties areinessentialin our subsequenttreatment, and
kinematic parameters. Another classical way to impose
we shall assume ‘constant external forces’ are applied,
some boundary conditions is to require, for all i,α in
and derive from a potential energy:
some subset I of I, V to depend linearly on one or
1 i,α
several parameters, e.g.:
Nf
W = Fextx . (4)
(i,α)∈I1 Vi,α =Ai,αλ1, (1) −µX=1 µ µ
introducing some collective ‘generalized velocity’ λ .
1 It is easily checked that such a definition is devoid of
Once again, in such a case, N is reduced to count el-
f ambiguity in the following important cases.
ements of I (I I ), plus λ .
0 1 1
\ ∪
At least locally, it is possible to regard velocities and The set of grain center positions, as opposed to
generalized kinematic parameters (like λ1 in eqn. 1) as • grain orientations, define a ‘flat’ space, on which
timederivativesofspatialcoordinates,whichweshalldo constant vectors and covectors are unambiguous.
dX
in the following, thus writing, e.g., V = i,α. As we Whenever external efforts are not sensible to ori-
i,α
dt entational coordinates, as in the case of gravity (if
are only interested in those properties that do not de-
the grain ‘centers’ are their centers of mass), one
pend on dynamics, grain trajectories might as well be
may therefore ‘apply constant forces’.
described by any parameter, not necessarily by physi-
cal time. In the case of kinematic constraints of type
AnticipatingonpartIV,theapproximationofsmall
1, parameters A will be regarded as fixed, although •
i,α displacements assumes that the manifold might lo-
positions of the grains and the walls change. One then
cally be replaced by its flat tangent space.
defines a generalized coordinate Λ , such that dΛ1 =λ .
1 dt 1
Jreufsetrsliktoeftohrevwelhooclietiseest,tohfepcoosmitipoancatlncootoartdioinna(txesµ.)1≤µ≤Nf toTahsethceomlopalde.teSNomf-evteimcteosr, oitfiesxctoernnvaelnifeonrctetsoisderaelfewrriethd
External forces and torques may at will be exerted on parametrized sets of loads. When the direction of the
thegrainsthatarefreeofkinematicconstraints. Weshall loadisfixed,whileitsintensitymightvary,onehasaone-
use the same notations as for velocities, writing down parameter loading mode. Insucha situation,allexternal
large Nf-vectors of ‘external forces’ (some of their coor- force components are kept proportional to a single load-
dAitnaetqeusilsitbarniudmin,gt,haecytuaarlelyo,ffocrotuorrsqeuteos)b,easb(aFlaµenxct)e1d≤bµ≤yNinf-. iQn,gλp,acraanmbeteeirdQen,taifineddaogneenqeuraaltiiznegdtvheelopcoitwyecroonfjuexgtaetrentaol
ternal forces (Fµint)1≤µ≤Nf: forceswiththeproductQλ. λissomelinearcombination
ofthe kinematicdegreesoffreedom(v ) ,andthe
(1 µ N ) Fext+Fint =0. (2) time derivative of a generalized coordµin1a≤tµe≤ΛNf, equal to
≤ ≤ f µ µ
4
the same combination of coordinates (x ) . The These two slightly different boundary conditions (BC)
potential energy is then simply µ 1≤µ≤Nf are respectively abbreviated as BC1 and BC2 in the fol-
lowing.
W = QΛ. (5) System B (fig. 2) is a hexagonal sample of the same
−
material. It is submitted to external forces on the pe-
Let us now illustrate those notions with simple exam- riphery, which mimic hydrostatic pressure.
ples, that will be repeatedly used in the following. Sys-
tems A and B are packings of discs that are placed on
the sites of a regular triangular lattice. (Later on, we
shallallowfor aslightpolydispersityofthe grains. They
might move, gain or lose contacts with their neighbours,
and the lattice might be slightly distorted). System A
(fig.1) is a pile with slope inclined at 60 degrees with re- 1 2 3
spect to the horizontaldirection. Eachdisc is submitted
4 5 6 7
8 9 10 11 12
y 36
34 35 13 14 15 16
31 32 33
17 18 19
27 28 29 30
22 23 24 25 26
FIG. 2: System B : a hexagonal sample. Arrows depict ex-
ternal forces applied on peripheral discs.
16 17 18 19 20 21
SystemC(fig.3)isadisorderedcollectionofdiscswith
9 10 11 12 13 14 15
a larger polydispersity. It is embedded within a circular
x
wall the radius R of which might change. One controls
1 2 3 4 5 6 7 8
the generalized force conjugate to λ = dR, viz.
1 dt
Q = f , (9)
FIG.1: SystemA :apile undergravity. Thebottom bound- 1 iw
ary conditions are explained in thetext Xi
wherethesumrunsoverallparticlesiexertingforcesf
iw
normally onto the wall.
toitsownweight,exceptthoseofthe bottomrow,which
collectively set the boundary condition. One might keep
them fixed at regularly spaced positions, imposing, say B. The structure: a set of bonds.
(numbering them as onthe figure,and denoting asa the
lattice spacing)
The definitions we introduce here pertain to one spe-
cific configuration of the grains, with the positions and
x =(i 1)(1 Λ )a
(1 i 8) i − − 1 , (6) orientations fixed.
≤ ≤ (cid:26)yi =0 We call ‘bonds’ the pairs of neighbouring grains that
may exert a force on one another. We require this force
allowingforahorizontaldeformationparameterΛ . One
1 to be concentrated at the point of each grain which is
mayalsorequirethemtostayonthehorizontalaxisy =0
the closest to the other one, and directed normally to
and satisfy
the surface.[57] The more generalcase of arbitrary bond
forces will be briefly evoked later. Note that we neither
vy = λ (i 1)a, (7)
i − 1 − require the grains that are joined by a bond to be in
contact,norimposeanysignconstraintonthe force. We
withafreekinematicparameterλ . Accordingtoeqn.3,
1 thus define, somewhat arbitrarily at this stage, N such
the generalized force conjugate to λ is
1 bondsasdepicted onfigure4, alternativelylabelledwith
an index l, 1 l N, or with the pair of labels of the
8
Q1 =Xi=1Fii,nxt(i−1)a. (8) dtwenoogtreasinthsethuen≤yitjovi≤enc.toIfrbthoantdplocionntsnefrcotsmiiantodjj,,nnolromranlliyj
5
the form of internal forces ((Fint) ) and torques
i 1 i n
((Γint) ) in the system is speci≤fie≤d: they linearly
i 1 i n
depend≤on≤bond forces f , as
23 ij
25
26 24 22 Fint = f n
10 i − ij ij
27 12 11 9 Xj6=i . (10)
Γint = f R n
21 i − ij ij ∧ ij
28 13 4 3 Xj=i
6
8
5 20 Given the load (Fext) , equilibrium requires, in
29 14 1 2 view of eqns. 10 andµ 21t≤hµa≤tNtfhe bond forces (fl)1 l N
19 satisfy equations of the form ≤≤
6 37
7
30
N
15 16 18 36 (1≤µ≤Nf) Hµlfl =Fµext, (11)
31 17 Xl=1
32 35
33 defining a linear operator, H : IRN . Bond forces
34 → F
(f ) are then said to be statically admissible with
l 1 l N
the l≤oa≤d (Fext) . Bond forces that are statically
admissiblewµith1a≤lµo≤aNdfequaltozero(inequilibriumwith-
FIG.3: SystemC:adisorderedpackingsurroundedbyacir- out any external action) are the elements of a subspace
cularwallthatmightuniformlyexpandorshrink,asindicated S of IRN, the null space of operator H. Its dimension,
0
by thesmall arrows. that we denote as h, is the number of linearly indepen-
dentsuchself-balancedsetsofinternalforces,or,inother
nij words, the degree of indeterminacy of bond forces in the
system (also called the degree of hyperstaticity). If not
empty, the set of statically admissible bond forces is an
affine space of dimension h.
Therelativenormalvelocityofthegrainsiandjjoined
R by a bond is
R ji
ij
δV =n (V V +Ω R Ω R ), (12)
ij ij i j i ij j ji
· − ∧ − ∧
Grain j
with the convention that it is positive when the parti-
Grain i
cles are approachingeach other. Eqn.12 defines a linear
operator, G, acting on into IRN. The range of G is
V
the subspace of compatible relative normal velocities,
h
ij i.e., those N-vCectors for which one can effectively find
valuesforthe velocities,relations12being satisfied. The
null space of G is the vector space M of ‘mechanisms’,
alsocalled‘floppymodes’,i.e., motionsthatdonotalter
the lengths h of the bonds. Its dimension, denoted as
FIG. 4: Two grains i and j joined by a bond. hij is the l
k in the sequel, is the number of independent such mo-
minimum distance between their surfaces, measured where a
common normal unit vector is nij. Vector Rij (respectively tions, or, in other words, regarding the bonds as rigid,
Rji) points from the center of i (resp. of j) to the point of the degree of indeterminacy of velocities, also called de-
its surface that is closest toj (resp. to i) gree of hypostaticity. Imposing the condition δV = 0
ij
in all bonds of the structure restricts the possible values
of velocities (v ) to a vector space of dimension
tothesurfacesofbothgrainswherethedistancebetween k. Dependingµon1≤tµh≤eNftype of load and boundary con-
them , h , is the smallest. R is the vector joining the ditions, the whole set of grains might keep some overall
ij ij
center of graini (origin), to the point onits surface that rigid body kinematic degrees of freedom. System B, for
isclosesttograinj (extremity). Thiscontactzonemight instance, has 3 independent such motions, as any solid
transmitanormal force,alongn ,ofmagnitudef that body in 2D. If k d(d+1)/2 denotes the number of
ij ij 0
≤
will be counted positively when the grains repell each such particular motions allowed by the boundary condi-
other. Oncethissetofbondsisdefined,itisreferredtoas tions,thesystemissaidtoberigid whenitdoesnothave
the structure. The set of bonds defined by intergranular other mechanisms, i.e., when k =k .
0
contacts (h =0) will be called the contact structure. An important and useful result, the classical theorem
ij
As a consequence of the definition of a structure, of virtual power states the following. Let (δV ) be
l 1 l N
≤≤
6
any element of , corresponding to the velocity vector bonds according to the interaction law. One may, for
C
(v ) , and let (f ) be a set of bond forces example,declareabondtojointwograinswhenevertheir
staµti1c≤aµll≤yNafdmissible wilth1≤tlh≤eNload (Fext) . One surfaces are separated by a minimum distance smaller
µ 1 µ N
then has: ≤ ≤ than some threshold h > 0. The choice of a larger h ,
0 0
therebyincreasingN,willdecreasek and/orincreasethe
N Nf degree of hyperstaticity h.
f δV = Fextv . (13)
l l µ µ Let us remark that the properties we have just dealt
Xl=1 µX=1 with in the case of bonds that carry normal forces, are
very easily generalized to the case of arbitrary contact
Equality13,forarbitrary(‘virtual’)equilibriumsetofin-
forces, at the cost of minor modifications. Relative nor-
ternalforcesand velocities,stressesthe geometric mean-
mal velocities and normal contact forces are replaced by
ing of forces and the mechanical meaning of velocities.
d-vectors, IRdN replaces IRN, equalities 13 (with, now,
It is easily established in two steps: first use the force
a scalar product within the sum in the left-hand side)
balance equations in the right-handside; then transform
and14arestillsatisfied. Insteadof16, oneends up with
the sum over degrees of freedom into a sum over bonds.
dN +k = N +h. Adding friction increases h and/or
As a direct consequence of the theorem, one deduces f
decreases k.
that operator H is in fact (as one might check directly,
Returningtofrictionlesssystems,thecaseofspheresor
readingthematrixelementsineqns.12and11)thetrans-
pose of G : H = GT. This follows from the sequence of discsdeservesaspecialtreatment: nonormalforceisable
toexertanytorque,andallrotationaldegreesoffreedom
equalities
arethereforemechanisms. Itisconvenienttoignorethem
(f δV)=(f Gv)=(Hf v)= GTf v , altogether. Their number nd(d−1) (n is the number of
| | | | 2
(cid:0) (cid:1) particles) is then subtracted both from Nf and from k,
valid for arbitrary v (such that Gv = δv) and f (such and eqn.16 still holds. Such granular systems are then
thatHf =Fext), inwhich abracketnotationis usedfor analogous to ‘central-force networks’: networks of freely
scalar products. Consequently, S0, the null space of GT, articulated bars, or systems of threads tied together, in
isthe orthogonalcomplementaryto ,the rangeofG,in which only the translational degrees of freedom of the
IRN: C nodes matter. One should be aware, however, that the
presenceoffrictionreinstatesrotationsintotheproblem.
S0 =C⊥. (14) We now illustrate the notions and properties intro-
ducedinthis sectionwithexamplesofstructuresdefined
ThustocheckthatsomevaluesδV thatonemighttryto
l insystemsA,BandC,ignoring,asexplainedjustabove,
assign to the relative normal velocities are compatible,
disc rotations.
it is sufficient to ensure the orthogonality of N-vector
FirstconsidersystemB.Threedifferentstructuresare
(δV ) to all N-vectors of self-balanced bond forces
(orla1s≤pla≤nNning subset thereof): apparent on figure 2. The first one, that we denote as
SB1,isthesetofbondsthataredrawnasthicklines;the
second, SB2, contains all bonds of SB1, plus those that
(δV ) . (15)
l 1 l N 0
≤≤ ⊥S are drawn with thin continuous lines on the figure; and,
One thus uses forces (elements of S ) as cofactors in a finally, the third structure, SB3, comprises all possible
0
set of geometric compatibility conditions. bonds betweennearestneighboursinthe system,i.e., all
Recalling k (the number ofmechanisms) is the dimen- those of SB2 plus the dotted lines. Ignoring rotations,
sion of the null space M of G, one has one has Nf =2n=38.
Structure SB3 is a set of rigid triangles sharing com-
Nf =k+dim( ). mon edges with their neighbours. It is devoid of mecha-
C
nisms,exceptthe 3overallrigidbodydegreesoffreedom
As h=dim(S0), from 14, one also has: of the system. Thus k = 3. N = 42 bonds are present.
In view of eqn. 16, one has h=7. One can exhibit 7 lin-
N =h+dim( ).
C earlyindependentsystemsofself-balancednormalforces,
asfollows. Thesmallstructure,with12bonds,involving
Elimination of the dimension of from those two equal-
C 7discs,depictedonfig.5,allowstodefineonesuchsetof
ities yields the following relationship between the degree
forces. Noting that 7 such patterns are present on SB3
of hypostaticity, k, the degree of hyperstaticity, h, the
(centered on discs 5, 6, 9, 10, 11, 14 and 15), the right
number of bonds, N, and the number of degrees of free-
count is reached.
dom, N :
f
Structure SB2 is made of N = 35 bonds. It can be
N +k =N +h. (16) shown (on studying the properties of the corresponding
f
matrix G) to be devoid of self-balanced sets of forces,
As we will check on examples below, relation 16 holds h=0,andofmechanismsotherthanrigidbodymotions,
whatever the choice ofthe list ofbonds between objects, k =3. Thus N +k =Nf +h.
although it is of course desirable in practice to define Structure SB1, comprising N = 25 bonds only, still
7
contains, in addition, the two bonds drawn with dotted
lines(19-24and32-34). Dependingontheboundarycon-
dition, discs 1 to 8 either possess collectively one degree
of freedom (for BC2) and then N = 57, or none (for
4 f
f 3 BC1) and Nf =56.
SA2 has 57 bonds. It is devoid of mechanism (k = 0)
forwhateverBC.ForBC2,onealsohash=0andeqn.16
f -f -f f holds as an equality between the number of bonds and
1 the number of degrees of freedom. For BC1, one has
5 h = 1. Indeed, one may recognize, in the bottom left
-f -f 2
corner of the pile, with discs 1, 2, 3, 9 and 10, part of
the hyperstatic pattern of fig. 5. With BC1, one needs
f -f -f f not care about equilibrium of discs 1, 2 and 3 that are
perfectly fixed. A system of self-balanced bond forces is
thusfoundonattributingacommonvalue tothe normal
6 f 7 forces in bonds 1-9, 9-10, 10-3, and the opposite value
to the normal forces in bonds 2-9 and 2-10. In the case
of BC2, those forces do not balance, since the equilib-
rium equation for the collective degree of freedom of the
FIG. 5: A set of self-balanced normal forces. The 6 bonds bottomrow(a combinationofeqns.8 and10) is notsat-
of the regular hexagon perimeter carry some normal force f, isfied. As to SA1, it has the same properties as SA2,
whilethe6onesinvolvingthecentraldisclabelled1carrythe with 2 additional mechanisms (collective ones like that
opposite force.
of fig. 6).
Consider now structure SC that is shown, in system
C, on figure 3, with the lines connecting disc centers,
has h = 0. According to eqn. 16, it should possess 10
or joining discs to the wall, that define N = 70 bonds.
additional independent mechanisms. 2 of them are due
Takingintoaccountthedegreeoffreedomofthewall,one
todisc10,whichisnowcompletelyfree. 4othersinvolve
has N = 2n+1 =75. One may show h =0. Thus one
discs 5, 9, 12, and 19, which are still free to move in one f
has k=5. Twodiscs (10 and14)areentirely free,hence
direction. In the case ofa divalent disc like 5, this is due
4 mechanisms. The missing one is a globalrotation,as a
to the exact alignment, on the regular lattice, of bonds
solidbody, ofthe set ofallparticlesaroundthe centerof
4-5and5-6. Fourlesstrivialmechanismsaremorecollec-
thecircularcontainer,thewallremainingimmobile. Such
tive. One of them is shown on figure 6. Two structures,
a motion would not be possible if the same boundary
condition was used with another container shape.
C. The problem: the structure and the load.
1 2 3
Once a list ofbonds is chosen,thus defining the struc-
ture, we shall refer to the situation of the structure sub-
4 5 6 7 mitted to a given load as ‘the problem’.
Solvingtheproblemwouldmeanfindingthemotionor
equilibrium state of the system (determining, e.g., new
8 9 10 11 12 equilibrium positions and intergranular forces), once the
load, from an initial state of rest with no external force,
has been applied. We are not, of course, able to do that
13 14 15 16 at this stage, since no contact law relating the forces to
the relative motion of neighbouring particles has been
introduced. The only information available is that the
17 18 19 internalforcesarerequiredtobelongtosomevectorspace
that is known once the structure is defined, and to be
exerted on given points on the grain surfaces.
It is said that the load is supported by the structure
FIG. 6: A collective mechanism on structure SB1. Arrows if its application leads to an equilibrium state in which
represent disc velocities. internal forces, carried by the bonds of the structure,
balance the external ones.
SA1 and SA2, are defined, on fig. 1, in system A. SA1 is We can state a necessary condition for the load to be
made of all bonds drawnwith continuous lines, and SA2 supported: it must be possible to find statically admis-
8
sible intergranular forces. Necessarily, the N -vector of Analogously,structureSB1,alongwiththeloadshown
f
externalforcesmustlie inthe rangeofoperatorGT, i.e., on figure 2, defines an isostatic problem PB1 in spite of
it must be orthogonalto the null space M of G: the k =10 mechanisms. In particular,the loaddirection
(provided discs sit right on the regular lattice sites) is
(Fµ)(1≤µ≤Nf) ⊥M. (17) exactly orthogonal to the velocity vector represented on
figure 6. Structure SC, submitted to the following load:
This simply means that if the load is to be supported, it
must not set the mechanisms into motion. Such a load
is said to be supportable. All supportable loads are not (1≤i≤37) Feixt =0 (19)
always supported. (cid:26) Q1 =Q01,
By definition, the backbone of a structure is the set of
where a prescribed value Q0 is imposed to generalized
bonds l such that a list of statically admissible internal 1
0
force Q defined in eqn. 9, yields an isostatic problem.
forces (f ) exists with f =0. In the following we 1
shall alsol 1r≤elf≤erN, as ‘the backbl0on6 e’, to the set of grains Isostatic structures, on the other hand, are such that
allproblemsareisostatic,whateverthechoiceoftheload.
reached by such bonds.
More precisely, one requires all loads orthogonal to the
In general, a full mechanical characterization of the
overall rigid-body degrees of freedom to be supportable
equilibrium properties of the system requires some con-
with a unique determination of internal forces. Equiva-
stitutive law in the contacts. However, there are inter-
lently, both conditions h=0 and k =k are to be satis-
esting situations in which 0
fied. Both the degree of hyperstaticity and the degree of
condition17beingfulfilled,theloadissupportable; hypostaticity (excluding rigid-body motions) should be
•
equal to zero. This entails the well-known condition
if it is supported, then all intergranular forces are
•
uniquely determined by the equations of equilib- N =Nf k0, (20)
−
rium.
stating that the number of equilibrium equations (N
f
−
These two conditions define an isostatic problem. k0) is equal to the number of unknowns (N).
Further restrictions on internal forces are often en- Equality 20 is a necessary condition for the structure
forced in the form of inequalities. The definition of a to be isostatic, not a sufficient one. For example, in
supportable load is then modified accordingly, impos- the structure definedby the additionof the bondjoining
ing additional conditions, to be satisfied simultaneously discs19and24toSA1withthe firstboundarycondition
with17. Theirconsequenceswillbe discussedinsections (BC1), one has k0 =0, N =Nf =56, while h=k =1.
III and VII. Structure SA2, with BC2, is isostatic. SB2, with
N = 38 and k = 3, is isostatic. As to SC, it would
f 0
be isostaticuponremovalofgrains10and14,onlyifthe
D. Isostaticity: various definitions. global rotation of the set of grains with respect to the
wall were ignored. Of course, all those structures, as we
In section VI we shall see that equilibrium configura- are dealing with discs, are only isostatic if rotations are
tions of assemblies of rigid frictionless grains interacting ignored. Only problems with no external torque exerted
via contact forces only are generally such that the prob- on the grains are isostatic. This should be remembered
lem is isostatic. on comparing h and k with and without friction in such
Here, we first insist on the difference between an iso- systems.
static problem, as defined just above, and an isostatic As we shall see, isostatic problems, rather than iso-
structure, to be defined below. Once condition 17 is sat- static structures, naturally occur in some model granu-
isfied, the set of possible bond forces is an affine space lar systems. The distinction is relevant, for it accounts
of dimension h. One has an isostatic problem if both for disconnected or ‘dangling’ parts in disordered struc-
conditions 17 and h=0 are fulfilled. Some mechanisms tures like SC, and for the peculiarities of lattice models.
might still exist in the structure (k = 0), provided they Moreover, some systems can also spontaneously, as we
are orthogonalto the load direction.6 shall see, select a non-rigid (k > k0) equilibrium config-
Structure SA1 (figure 1), with discs exactly centered uration.
on the sites of a regular triangular lattice, is such that
the problem, denoted as PA1 in the following, defined
E. Generic versus geometric properties.
with BC2 and the following load: [58]
(9 i 36) Fext = pe The distinction between isostatic problems and iso-
(cid:26) ≤ ≤ Qi 1 = −52√943py, (18) static structures should not be confused with another
one: thatbetweengeometric andgeneric isostaticity. We
where p is the weight of one disc and e is the vertical have used a geometric definition of a structure, as asso-
y
upwardsunit vector,is isostatic,although 2 mechanisms ciated to one particular position of the system in con-
are present. figuration space, and accordingly the definition we gave
9
is that of geometric isostaticity. A topological one can importanceofgeometricalaspects. Thus wefirstpresent
be introduced which, irrespective of particle positions, the simplest case of rigid, frictionless and cohesionless
is only sensitive to the connectivity of the network of grains, in which contacts simply behave as struts. Then
bonds. In the case of spheres or discs, when rotations we introduce and briefly discuss other possible laws in
can be ignored, this amounts to regarding the structure whichunilateralityorrigidityconstraintsaremodifiedor
as a graph: a set of edges (bonds) joining at vertices relaxed. Most of those frictionless systems possess a po-
(grains). Operator G, spaces S , M and their dimen- tentialenergythatisstationaryatequilibriumstatesand
0
sionshandk smoothly dependonthe coordinatesofthe reaches then a minimum if they are stable. Throughout
grains, via vectors n and R . However, the rank of a this section, it is assumed that a one-parameter loading
ij ij
parameter-dependent matrix stays at its maximum ex- mode has been defined for varyingparticle positions and
cept for special values of the parameters. Equivalently, orientations, with constant external forces, and that the
the dimensionof the null space is genericallyequalto its potential energy of external forces,W, can be written in
minimum value. Applying this to both G and GT, one the forms of eqns. 4 and 5.
may define the generic degree of indeterminacy of veloc-
ities (with due account to the k rigid-body degrees of
0
freedom) k and the generic degree of indeterminacy of A. Rigid frictionless grains, no cohesion.
forces h as the respective generic (minimum) dimensions
of their null spaces. This allows to define a suitable iso- In this case, the contact law takes the form of the so-
staticity notion for topological structures: a generically called Signorini condition:
isostatic structure is one for which both numbers h and
k−k0 are equal to zero. fij =0 if hij >0 (21)
Itfollowsfromthe definitionsthatageometricallyiso- (cid:26)fij 0 if hij =0
≥
static structure is always, once regarded as a topolog-
ical structure, a generically isostatic one, but that the It shouldbe noted that this lawdoes notexpress a func-
reciprocal property is not true. Ref. [20] gives a coun- tionaldependenceoff onh . Letusstudythevariations
l l
terexample for a system of discs (like systems A and B, ofW nearequilibriumstates. First,considersuchastate,
equivalent to a network of articulated bars)on the regu- in which some non-negative contact forces f , in closed
l∗
lartriangularlattice. Inspecific configurations(likethat contacts(h =0)balancetheexternalloadQ. Letusap-
l
of a regular triangularlattice), one might exceptionnally ply the theoremof virtualpower with statically admissi-
have h=k k >0 on generically isostatic structures. ble force set (f ) , andarbitraryparticle velocities,
In two d−ime0nsions, there exists some powerful algo- correspondingl∗to1≤rel≤laNtive normal velocities δV = dhl
l − dt
rithms [21, 22] to evaluate the generic degrees of force and a value λ = dΛ for the kinematic parameter conju-
dt
and velocity indeterminacy in central-force networks (or gateto Q. Foranyl suchthatf >0,the Signorinicon-
l∗
systems of frictionless discs). Such computational meth- ditionrequiresthat h =0 andonemust haveδV 0 to
l l
ods only deal with connectivity properties, they do not complywiththeimpenetrabilityconstraints. Then≤,from
manipulate floating-point numbers and are therefore de-
void of numerical round-off errors. They were success- dW dΛ
= Q = Qλ= f δV
fully applied to systems of up to 106 nodes. However dt − dt − − l∗ l
theyareofcourseunabletocomputeposition-dependent Xl
quantities like force values.
it follows that any motion that does not lead to grain
interpenetrationcan only, to firstorder in t (any param-
eteronthetrajectoryinconfigurationspace)increase the
III. CONTACT LAW AND POTENTIAL potential energy. This non-negative first-order variation
MINIMIZATION. mightbeequaltozeroifδV =0foranyactivecontactl,
l
i.e.,ifamechanismexistsonthebackboneofthecontact
So far, the only restriction on intergranularforces was structure. Whether the equilibrium state correspondsto
that they should be normal to the grain surfaces.[59] In a minimum of W depends then on the sign of second
this section we consider some more specific cases of fric- or higher order variations. If the backbone of the con-
tionlessgrains,inwhichsome“contactlaw”,relatingnor- tact structure is rigid, then W is necessarily minimized
mal forces to relative positions, is known. This provides at equilibrium.
somelimitedadditionalinformation,thatisnotsufficient Conversely, let us assume that a configuration of the
in generalto predict the grain trajectories once they are grains has been reached, that locally minimizes W un-
submitted to external forces, for all dynamical aspects der the constraints hl 0. There must then exist some
≥
arestillunknownandthecharacterizationofequilibrium non-negative Lagrange multipliers fl, such that, for any
mightevenbeincomplete. Ouraimistodeduceasmuch coordinate xα,
aspossibleontheglobalpropertiesofthegranularassem-
∂ ∂h
bly from as little information as possible on the detailed (QΛ)= f l . (22)
l
mechanical laws of the contacts, in order to stress the − ∂xα −Xl ∂xα
10
Only for such indices l that h = 0 do the f take non- (bilateral) on the one hand, and cables (satisfying 23),
l l
vanishingvalues. Thepartialderivativeintheright-hand on the other hand. Their potential energy has the same
side of 22 is the opposite of matrix element G , while, properties as stated above.
l,α
from4,thatoftheleft-handsideistheexternalforcecon-
jugatetox . Thus,wehavejustwrittenthatparameters
α
f areinfactequilibriumcontactforcessatisfying21,and C. Systems with a smooth interaction potential.
l
reaction forces stem from geometrical constraints.
Wenowintroduceafewotherrelatedcontactlawsand
The model of perfectly rigid grains is physically rea-
mechanical models.
sonable when contact deformations (h < 0) are negligi-
l
ble in comparison with any other relevant length in the
problem. When this is no longer the case, or when one
B. Systems with tensile or bilateral forces.
wishes to model sound propagation, it is appropriate to
deal with contact laws that involve elastic deformations,
Networks of rigid strings or cables are analogous to e.g.,
frictionlessspheres(ignoringtheirrotations)ifthesignof
forces is reversedand if the distance constraint h 0 is f =0 if h >0
replaced by hl 0. The Signorini condition 21 belc≥omes (cid:26)fiijj =Kij hij m if hiijj 0, (26)
≤ | | ≤
f =0 if h <0 inwhichK isastiffnessconstantthatdependsonmate-
ij ij (23) ij
(cid:26)fij 0 if hij =0, rial properties and on the geometry of contact i,j. The
≤
exponent is m = 3/2 (Hertz law) for smooth surfaces
and the whole treament of the preceding subsection in 3D, and other values might model roughness and the
straightforwardlyapplies. presence of conical asperities [24, 25].
In the case of non-spherical grains, an analogous sys- Suchcontactforcesderivefromanelasticpotentialen-
tem supporting tensile forces is an idealized chain, in ergy:
which ‘grain’ -chain links- perimeters are free to cross.
Pairsofneighbouringlinks(interpenetrating‘grains’)ex- N
K
ert a force on one another, opposing their separation, Wel = w(hl),with w(hl)= l hl m+1. (27)
m+1| |
when their intersection reduces to a contact point. Xl=1
A bilateral contact law:
Likewise, rigid cables as introduced in subsection III.B
f =0 if h =0 couldbereplacedbyelasticones. Thatstableequilibrium
ij ij 6 (24)
(cid:26)fij unknown if hij =0, states correspond to minima, in the absence of frictions,
of the total potential energy
might model rigidcohesivegrains,that ‘stick’to one an-
other. Thestickingforcemightbelimitedbyanunequal- Wtot =Wel+W, (28)
ity:
sumoftheelasticpotential27andthepotentialenergyof
f =0 if h =0 external forces 4 or 5, is an extremely familiar property.
ij ij 6 (25)
(cid:26)fij f0 if hij =0, The Signorini condition might physically be regarded as
≥−
the limit of the interaction law expressedby equation26
When one simply uses the form 24, assuming the pairs whenthestiffnessconstantsbecomeverylarge,or,equiv-
thatarestuckincontactwillnotcomeapart,theconclu- alently,whenthe levelofintergranularforcesapproaches
sionsofsubsectionIII.Astillhold,ifunilateralconditions zero. Alternatively, it is mathematically possible to in-
on relative velocities and displacements are replaced by troduce a regularized contact law of the form 26 as an
bilateralones,andifallsignconstraintsoncontactforces approximation, when contacts are stiff enough, of the
are removed. Equilibrium configurations are character- ideal impenetrability constraint. Such a point of view
izedbystationarityofpotentialenergyW. Minimization is adopted in optimization theory: the procedure known
of W ensures stability. A sufficient, but not necessary as penalization of the constraints amounts, instead of
condition for minimization of W is the rigidity of the minimizing W subject to impenetrability constraints, to
backbone of the contact structure. searching for unconstrained minima of W +Wel.
Reciprocally, statically admissible normal contact Tensilecontactforcesoflimitedintensity,asincontact
forces naturally appear as Lagrange multipliers associ- law 25, might result from some attractive interaction of
ated with bilateral constraints h = 0 at a potential en- finite, but small, range, as depicted on fig. 7. It is in-
l
ergy minimum. teresting to note that the addition of an attractive tail
However, contact law 25 does not lend itself to a po- hasturnedthepotentialw(h)intoanon-convexfunction
tential energy formulation. of interstitial thickness h. At the inflexion point, A, the
Tensegrities[23](withrigidelements) areby definition attractive force reaches its maximum f . If one pulls,
A
mixednetworksofstruts(satisfyingcondition21)orbars withagrowingforce,ontwograinsincontactinorderto
