University of Southern California Dissertations and Theses

An experimental study of the elastic theory of granular flows

(USC Thesis Other)

An experimental study of the elastic theory of granular flows

PDF

Download a page range

Download transcript

Copy asset link

Request this asset

Transcript (if available)

Content

An Experimental Study of the Elastic Theory
of Granular Flows
by
Tongtong Guo

In Partial Fulfillment of the Requirements
for the Degree of
Doctor of Philosophy

Faculty of The USC Graduate School
University of Southern California
December, 2017

Abstract

This paper reports annular shear cell measurements of granular flows with an eye towards
experimentally confirming the flow regimes laid out in the Elastic theory of granular flows. Tests were
carried out on four different kinds of plastic spherical particles for both constant volume flows and
constant applied stress flows. In particular, observations were made of the new regime in that model,
the Elastic-inertial regime, and the predicted transitions between the Elastic-Inertial and both the
Elastic-Quasistatic and pure Inertial regimes. Results from Discrete Element Method(DEM) computer
simulation using spring-dashpot-slider particles contact model were compared to experiment results
evaluating its performance. The results indicated that the particles’ friction coefficient decrease with
increasing normal force, and that including this effect is important for an accurate match between
simulation and experiment. It’s possible to fit a simulation to experimental data, although one cannot
use measured particle properties. Instead one must go through the laborious process of running many
simulations and trying to fit to the data. It shows that an accurate match requires complicated material
modelling.

Contents
1. Introduction .......................................................................................................................................... 1
2. Literature review ................................................................................................................................. 10
2.1 Quasi-static flow theory (slow flows) ......................................................................................... 10
2.2 Kinetic theory (rapid flows) ......................................................................................................... 17
2.3 Elastic theory ............................................................................................................................... 30
2.4 contact problem for two spheres ............................................................................................... 39
2.4.1 Constant normal displacement, varying tangential displacement ..................................... 41
2.4.2 Varying normal displacement, varying tangential displacement........................................ 44
2.5 Simulations .................................................................................................................................. 48
3. Experimental setup ............................................................................................................................. 52
3.1 Apparatus .................................................................................................................................... 52
3.2 Test Materials ............................................................................................................................. 54
3.3 Method and procedures ............................................................................................................. 60
4. Experimental results ........................................................................................................................... 61
4.1 Controlled volume results ........................................................................................................... 61
4.2 Controlled stress results ............................................................................................................. 66
5. Simulation results ............................................................................................................................... 70
6. Conclusions and Summary .................................................................................................................. 92
Acknowledgement ...................................................................................................................................... 94
Appendix ..................................................................................................................................................... 95
A.1 Review on the Inertial Number model ............................................................................................. 95
A.2 Comparison of three different frictional models .............................................................................. 98
References ................................................................................................................................................ 107

1

1. Introduction

The first mention of granular flows is probably by Lucretius (ca. 98-55 B.C.E.), “One can scoop up
poppy seeds with a ladle as easily as if they were water and, when dipping the ladle, the seeds flow in a
continuous stream.” (Duran, 2010)
Granular material is a bulk collection of discrete macroscopic solid particles. Those particles can
be of any size or shape although in colloquial usage, the term is limited to particles many orders of
magnitude larger than molecular scales and whose motion is not affected by temperature or Brownian
motion. Smaller particles, those with large surface area to volume ratios, are subject to interparticle
cohesion resulting from van der Wall’s forces, electrostatics or liquid bridges, but these effects are
negligible for larger particles. Examples of common granular materials are sand, rice, balls in a ball pit in
a playground, both coffee beans and ground coffee, and pills.

Figure 1. Some commonly seen granular materials.Picture is taken from
https://upload.wikimedia.org/wikipedia/commons/6/62/Granular_matter_examples.PNG

2

Granular materials play an important role in industries such as mining, agriculture, chemical,
cosmetic and pharmaceutical manufacturing. Coals and various minerals are mined, transported and
stored in granular form. It’s popular for chemical and pharmaceutical industries to prepare ingredients
in particulate forms both for better handling and because the increased surface area to volume ratio
improves chemical reaction rates. It’s also clear that natural processes such as landslides and the motion
of sand dunes are closely related to the flow of granular materials. As it is basic and cost-free, gravity is
the most common driving force for the flow of granular materials. As a result, hoppers and chutes
which utilize gravity are the most widely used tools for the storage and transportation of granular
materials. Examples of industrial hoppers and chutes are shown in Figure 2. Although it’s a most
common form of motion, the basic mechanisms governing the flow of granular materials are not well
understood and the tools used to handle granular materials are primitive.

(a)

(b)
Figure 2 Pictures of some widely used tools in storage and transportation of granular materials. (a)Hopper. Picture was taken
from http://www.claytonlambert.com/wp-content/uploads/2014/07/industrial.jpeg (b) Chute. Picture was taken from
http://www.thortsenmagnetics.com/Portals/0/5-blue-chute-06.jpg

Granular materials behave differently under different conditions. On one hand, they exhibit
liquid like characteristics. The bulk of the material can change its shape to accommodate themselves to
3

external boundaries and forces applied. Properties similar to buoyancy forces and hydraulic forces are
also found within granular materials. In flowing granular materials under gravity, heavier particles sink
to the bottom of the flow and lighter particles flow on top within the motion. If a rotating cylinder is
filled with granular particles, a surface with much larger curvature than particle diameter may be
observed. On the other hand, granular materials sometimes behave like solid and of course support the
weight of large buildings. In soil mechanics, engineers treat soil as a continuum object and classical
mechanics are applied to calculate the stress, strain and plastic deformation within the soil. The bulk
material can withstand shear forces like in an elastic object.
Jamming happens in many granular flows, especially in flows through hoppers, resulting in flow
stoppage that must be manually cleared. For uniformly sized, spherical and incompressible particles,
jamming usually happens at solid concentration 𝜈 ~ Random Close Packing (RCP) density. (The solid
concentration ν is the fraction of a unit volume filled with solid material, the rest of the volume making
up the space between the particles.) RCP is a loosely defined state that the material would assume if,
say, poured into a container. Closer packings are possible if the particles are regularly ordered. Figure 2
shows ordered close packing in two dimensions and three dimensions. In two dimensions, hexagonal
lattice packing achieves the highest possible solid concentration and Face Centered Cubic (FCC) packing
is the corresponding form in three dimensions. In two dimensions, 𝜈 ℎ𝑒𝑥
=
𝜋 2√3
≅0.91, while in three
dimensions, 𝜈 𝐹𝐶𝐶 =
√2𝜋 6
≅0.74. It’s found from simulations and experiments that for a large number of
particles, the RCP density is 𝜈 𝑅𝐶𝑃 2𝐷 ≅0.84 and 𝜈 𝑅𝐶𝑃 3𝐷 ≅0.64. There is no current mathematical
explanation for the RCP. Although randomly packed structure of particles seems fragile, RCP is
extremely stable. In 2 dimensions, any perturbation to the randomly close packed uniform disks may
rearrange the bulk into the hexagonal structure. Figure 3 shows the presence of both Random Close
Packing and Face Centered Cubic packing in the same pile of ball bearings.
4

Figure 3. Ordered Close Packing of the spherical particles in 2D and 3D. (a) hexagonal lattice packing (b) Face Centered Cubic
Packing

Figure 4 Random Close Packing and Face Centered Cubic Packing observed in ball bearings. Picture is taken from
http://www.pas.rochester.edu/~stte/BSSCM/lecture1.pdf

In general, the space between granular particles is always filled with some kind of interstitial
fluid. In most cases, granular materials are dry which means gas (usually air) is the only interstitial fluid.
When granular particles are large and heavy enough, it can be assumed that they are not affected by air
in between them; small light particles are exceptions. When interstices are filled with liquid, the
granular materials are wet and the bulk behavior changes dramatically. The difference between wet
sands that are saturated by sea water and dry sands that are further away from the sea is a perfect
5

example. When walking on the beach, one would immediately notice that it’s easier to walk on wet
sand than it is to walk on dry sand.
There are two internal stress transport mechanisms. For the stress tensor 𝝉 , each component 𝜏 𝑖𝑗

is the stress in 𝑖 -direction on a surface with outward pointing normal unit vector in 𝑗 -direction. One way
forces got transferred is through contacts with other particles. The contact stress tensor 𝝉 𝒄 , is defined
as an averaged dyadic product: 𝝉 𝒄 =<𝑭 𝒄 𝒍 >, where 𝑭 𝒄 is the contact force at the point of contact
between two particles and 𝒍 is the vector connecting the centers of the two particles in contact. The
brackets <> mean 𝝉 𝒄 is averaged over time and volume. In contrast to continuum solid mechanics, 𝝉 𝒄
is not necessarily symmetric, since due to friction, the contact force 𝑭 𝒄 does not have to be in the same
direction as 𝒍 , as shown in Figure 4. Furthermore, the averaging volume (control volume) cannot be
infinitesimal and has to be at least the size of two particles so that the common proof that the stress
tensor is symmetric, which involves shrinking the control volume to zero, does not apply.

Figure 5 Force generating by contacts may not be in the same direction of the line connecting the center of two particles.

Since the individual particles and the bulk of the material may travel at different speeds when in
motion as can be seen in Figure 5, momentum will be exchanged through a mechanism analogous to
6

kinetic theory or turbulence. If 𝒖 ′
is the fluctuating velocity (the difference between the instantaneous
particle velocity and the local average velocity), 𝜌 𝑝 is the particle density and 𝜈 is the solid
concentration, (note 𝜌 𝑝 𝜈 is the density of the bulk material) then streaming stress tensor 𝝉 𝒔 can be
defined as: 𝝉 𝒔 =𝜌 𝑝 𝜈 <𝒖 ′
𝒖 ′
>. It has a similar form to the Reynold Stress Tensor in a turbulent flow. In
most granular flows, the contact stress is the dominant stress and the streaming stress tensor can be
neglected. Exceptions are some high shear rate flow with small solid concentrations, where particles
travel some distance between collisions and the random velocity 𝒖 ′
is large, and one cannot neglect 𝝉 𝒔 .

Figure 6 Streaming stresses are generated by velocity difference 𝒖 ′
between individual particles and the bulk of material.

When two spherical particles are in contact, the contact force does not linearly with the
deformation due to the geometry of the contact. Figure 6 shows the model of the contact between two
particles. The contact force can be expressed as: 𝑓 𝑛 =𝐸𝐴𝘀 , where 𝐸 is the Young’s Modulus, 𝐴 is the
contact area and 𝘀 is the strain. For the strain, it is defined as 𝘀 =
𝛿 𝐿 , where 𝛿 is the deformation of the
particle at the contact point and 𝐿 is a carefully chosen length scale. When two particles are in contact
but with no force applied, 𝐴 is zero. With increasing the force pushing two together, both 𝐴 and 𝘀
change simultaneously, making a nonlinear force-deformation relation.
7

Figure 7 Model of contact between two particles. (a) The initial contact area is zero when no force is applied. (b) Both 𝜹 and A
change with the load, generating nonlinearity.

The contact problem of two objects with curved surface was first solved by Hertz in 1822,
providing an implicit solution. Assuming a small contact area 𝐴 , a force-deformation relation was
derived by Hertz, 𝑓 𝑛 =
4
3
𝑅 1
2
𝐸 1−𝑣 2
𝛿 3
2
, where in this case 𝑣 is Poisson’s ratio. When a large number of
particles are present, Bathurst and Rothenburg found that the elastic modulus of the bulk was:
𝐸 𝑏𝑢𝑙𝑘 ~𝑓 (𝑛 )
𝑘 𝑅 , where 𝑘 is the stiffness and 𝑛 is the coordination number (the average number of
contacts between a particle and the particles adjacent to it). The stiffness is defined as: 𝑘 =
𝑑 𝑓 𝑛 𝑑𝛿
=
𝑅 1
2
𝐸 1−𝑣 2
𝛿 1
2
. Note from the definition of 𝐸 𝑏𝑢𝑙𝑘 , the characteristics of the bulk material is associated with
the geometry of individual particles (R), material properties of individual properties ( 𝐸 and 𝑣 ) and
internal structure of the bulk material (through the coordination number 𝑛 ).
Leonard Da Vinci (1452-1519) was probably the first person who used a scientific approach to
study friction. Two hundred years before Isaac Newton (1643-1727) defined forces, Da Vinci found that
if loads were doubled, friction would be doubled. He stated that because the friction between
contacting rocks was proportional to the force pushing them together, the strength of a pile of sand was
8

associated to the pressure applied to it. This was the earliest quantitative statement related to the study
of granular materials.
The door into research on granular materials was opened by Charles de Coulomb (1736-1806),
in his paper “Essai sur un application de règles de maximis et minimis à quelques problèmes de staique
relatifs à l'architecture” which translated to “Essay on the rules of Maximus and Minimus applied to
some problems of equilibrium related to architecture”. Based on various tensile and shearing tests, he
investigated the failure planes with in soil and stone structures. He came up with two parameters to
determine the strength of soil structures, cohesion (𝑐 ) and shear resistance (or friction angle) (𝜙 ). His
findings was later summed to the famous Mohr-Coulomb criterion for a granular material to be statically
stable: 𝑇 ≤𝑐 +𝑁 𝑡𝑎𝑛𝜙 , where 𝑁 is the normal stress and 𝑇 is the shear stress. Within a granular
material, yield occurs when:
𝑇 =𝑐 +𝑁 𝑡𝑎𝑛𝜙 (1)
Erbst Chladni (1756-1872), did an interesting experiment using a thing vibrating plate covered
with sand as described in his book Entdeckungen über die Theorie des Klanges ("Discoveries in the
Theory of Sound") , which was published in 1787. He found that when the plate reached resonance,
sand moved to form small heaps. Those heaps outlined the nodal line of the resonating plate, where the
surface of the plate stayed still. The patterns of the formed heaps are called Chladni Figures. Similar
phenomenon was observed by Michael Faraday (1791-1867) and later named Faraday Waves.
Around 1885, Osborne Reynolds (1842-1912) described the dilatancy in granular material.
Imagine a balloon filled with sand saturated by water, slightly deform the balloon would cause the water
level to drop, contradicting common sense which would assume the water level to increase. He found in
order to deform a compact granular material, the bulk has to expand and reduce its density to a critical
value first. The same phenomenon can be observed as the drying of wet sand around a footstep on a
beach.
9

In the twentieth century, one man turned a new page in the research into granular materials.
His name was Ralph A. Bagnold (1896-1990). His research on the behavior of sand helped understand
the formation and motion of structures like sand dunes and ripples. He was the first person trying to
model granular flow from a perspective of individual particles. He found normal stresses within rapid
shear granular flows increased quadratically with the shear rate. The Bagnold number (Bagnold, 1954),
which is the ratio of particle collision stress to viscous stress in a granular flow with interstitial
Newtonian fluid, is now an important parameter characterizing granular flows. Since Bagnold, people
realized how little was known and thus the necessity of studying granular flows. As a result, more and
more scientists are devoted into the field.

10

2. Literature review

We have already discussed the contact model and basic transport mechanisms of granular flow
as well as some of the early history of granular research. In the following part, we will first review two
limiting granular flow models, the slow flow or quasi-static model, and the Rapid flow or kinetic theory
model. At last, we will examine the Elastic theory, which by including particle stiffness as a scaling
parameter, draws the entire flowmap for granular flows.
2.1 Quasi-static flow theory (slow flows)

Quasi-static flows, as the name implies, focus on flows that have relative small flow speed
where the developed stresses are independent of flow rate. Analyses of quasi-static flows are based on
metal plasticity theories that basically assume granular bulk behave solid-like until they yield at a certain
stress state. Despite the discrete nature of the granular materials, continuum models have been used
widely in solving both static and dynamic problems. Instead of investigating particles individually, in
these models, the granular bulk is treated as continuous media in that properties such as density and
shear rate are considered to be continuous functions of time and position. With such an approach,
numerous discrete particles are replaced by the infinite sets of infinitesimal control volumes.
Mathematically, as smooth continuous functions are differentiable, continuum models make it possible
to use calculus, making the problems significantly easier to manage. There have been many discussions
over the conditions under which continuum models are valid. Many scientists (Batchelor, 1967; Chung,
1988; and Fung, 1977) believe that as long as the length scale of the system is much larger than the
average spacing between particles (distance between centers of two adjacent particles), the continuum
models are accurate. However, most granular flow experiments conducted in laboratory do not meet
this criterion. The length scale of apparatus used in experiments often does not go beyond 100 particle
11

diameters, as in Brown and Richard (1970), Brennen and Pearce (1978) and Johnson et al (1990).
Chapman and Cowling (1964) suggested that the continuum model was valid if the time scale was large
enough comparing to the average time for a particle to move through control volume 𝑑𝑉 , even though
𝑑𝑉 was not much larger than particle diameter.
Once a granular material is considered a continuum, it is subjected to the balance law of
continuum mechanics. There have been numerous attempts to incorporate the behavior of granular
materials into continuum constitutive equations.
The first quantitative analysis of granular material’s yield was done by Coulomb as described by
equation (1). For a material yielding very slowly or about to yield, the components of the stress tensor 𝜎
satisfy the force balance law, which is described by the differential linear momentum balance equation
shown below:

𝜌 𝐷 𝐮 𝐷𝑡
=−∇∙𝛔 +ρ 𝐛
(2)
Where 𝜌 ,𝐮 ,𝑎𝑛𝑑 𝒃 are the density, velocity and body forces; 𝝈 is the stress tensor.
As the Coulomb yield condition in (1) predicts, yield happens on two planes which are called
𝑠𝑙𝑖𝑝 𝑝𝑙𝑎𝑛𝑒𝑠 . The normals to the slip planes are at ±(
𝜋 4
+𝜙 /2) relative to 𝜎 1
axis, where 𝜎 1
𝑎𝑛𝑑 𝜎 2

denote the major and minor principal stresses . The slip planes themselves are inclined at angles of
±(
𝜋 4
−𝜙 /2) to 𝜎 1
. According to Coulomb yield condition, yield only happen on slip planes.
On the 𝜎 1
−𝜎 2
plane, the yield condition (14) can be expressed as:

𝜎 1
𝜎 2
=
1+sin𝜙 1−sin𝜙 (3)
This is shown in Figure 8. The shaded area is the region where material behaves elastically, while
yielding occurs at the boundary.
12

Figure 8. Coulumb yield condition on 𝜎 1
−𝜎 2
plane. Material behaves elastically in the stress zone corresponding to the shaded
area.

It should be noticed that the yield condition only applies to static problems and does not
accurately explain the motion beyond the yield. For example, it tends to over predict the dilatancy
observed by Reynolds mentioned in the introduction. Supplemental information about how granular
materials deform is needed in order to explain the kinematics of motion at yield.
In order to understand the kinematics of the yield, we need to discuss the results of one kind of
shear tests designed to study soil mechanics. Figure 9. A schematic diagram of the Jenike shear cell.
shows a device called Jenike shear cell, which was first developed by Jenike (1961). It consists of a
cylindrical cell split into 2 halves. The bottom half is stationery, where the top part is free to move
horizontally. The two parts are initially amounted eccentrically. After the device is loaded with granular
material, it is covered with a lid, upon which a normal load 𝐹 𝑣 is applied. Notice that the lid is free to
move up and down according to the volume of the granular material. The top half is driven at a constant
speed, while the driving force 𝐹 ℎ
is measured. A Jenike shear cell is used to study the idealized shear
13

behavior between two infinite parallel plates under a constant normal load. Figure 9 shows an idealized
model of a Jenike shear cell. It’s assumed that at initial state, the granular material has a uniform density
𝜌 0
.

Figure 9. A schematic diagram of the Jenike shear cell.

Typical results of the test are shown in Figure 10. Depending on the initial bulk density 𝜌 0
and
normal stress 𝑁 , different behaviors are observed. On the shear stress vs. displacement (𝑇 −𝑋 ) plot,
when at a small displacement, the curve for a dense material or under a small load 𝑁 is almost linear. At
this stage, the behavior of the material is elastic and reversible. With increasing 𝑋 , it reaches a peak
value 𝑇 𝐹 , which depend on normal stress level 𝑁 and initial density 𝜌 0
. Beyond point 𝐹 , further
increasing 𝑥 , 𝑡 ℎ𝑒 shear stress 𝑇 decreases. At large 𝑥 , the shear stress 𝑇 reaches an asymptotic value
𝑇 ∞
(𝑁 ) , which is independent of 𝜌 0
𝑎𝑛𝑑 𝑜𝑛𝑙𝑦 𝑎 𝑓𝑢𝑛𝑐𝑡𝑖𝑜𝑛 𝑜𝑓 𝑁 . This behavior is often referred to the
process of failure (Atkinson & Bransby, 1982). The change in the height of the shear device, ∆ℎ,
associated with the shear motion is also shown in Figure 10a. Initially in the elastic region, the height
decreases slightly. After that, the material continues to dilate with an increasing ∆ℎ . Such a material is
called over-consolidated. Its behavior is denoted by the solid lines in Figure 10a.
14

When 𝜌 0
is small or the normal stress 𝑁 is large, the situation is quite different. The shear
stress- displacement curve displays an almost monotonic behavior. Similarly at small displacement,the
material responds elastically. Once it reaches a certain point 𝐶 , the process can no longer be retraced
reversiblely. Eventually as 𝑥 increases, 𝑇 approaches a asymtotic value 𝑇 ∞
(𝑁 ) , which is only a function
of normal stress 𝑁 . On the curve of ∆ℎ for a loosely packed material, it shows the bulk first gets
compressed significantly. After that, the material dilates. Such a material is under-consolidated. Its
behavior is denoted by the dash lines in Figure 10a.

Figure 10 The results from the Jenike cell tests. (a) Shear stress T and change in total height ∆ℎ vesus horizontal displacement X.
Solid line stands for the more condensed sample and dash line denotes the loosely packed sample. (b) schematic diagram of the
shear device. (results reproduced from Feda, (1982))

15

Eventually, both cases reach a state of isochoric deformation, where the shear continues with
no further density or volume changes. Roscoe et al. (1958) named this state of
deformation 𝑡 ℎ𝑒 𝑐𝑟𝑖𝑡𝑖𝑐𝑎𝑙 𝑠𝑡𝑎𝑡𝑒 . Commonly, a subscript 𝑐 is used to denote the 𝑐𝑟𝑖𝑡𝑖𝑐𝑎𝑙 𝑠𝑡𝑎𝑡𝑒 . It was
found that at the critical state for a certain type of material, the ratio of the shear stress 𝑇 𝑐 to normal
stress 𝑁 was nearly constant. Also, at the critical state, the solid fraction 𝜈 𝑐 was only related to the
normal stress 𝑁 .
For an over-consolidated material, following the failure which corresponds to the peak (point 𝐹 )
on Figure 10a, the shear stress decreases until it reaches its critical value. It’s obvious that the material
becomes weaker following failure. In that case, if the shear motion which causes the weakening
happened first within a narrow layer, both shear stress and density decrease in that rupture region.
Meanwhile, the structure of the adjacent areas remains mostly unchanged and therefore can withstand
a higher shear stress. As shown in Figure 11b, this causes the shear motion or the deformation to be
limited mainly to the narrow layer where the rupture initiates, while materials outside the rupture layer
stay rigid. Roscoe (1970)and Vardoulakis et al. (1991) found the average width of the narrow layer
within which shear occurs was about 10-20 particle diameter in width. For an under-consolidated
material, the deformation is uniformly distributed through the body (Figure 11c). Since the shear stress
increases as the material gets compressed, indicating strengthening within the material, this process is
often referred to the consolidation of granular material.

16

Figure 11. deformation in plane shear. (a) undeformed sample (b) deformation of dense sample (failure) (c) deformation of
loose sample (consolidation).

The Quasi-static models have been used widely in predicting the failure of soil structures in Civil
Engineering. Jackson (1983) examined mathematically some early plasticity theory based continuum
models used to describe granular flows. It was found these models over predicted the dilation of
granular materials. As mentioned earlier in the introduction, hoppers are probably one of the most
widely used tools in storage and transportation of granular materials. An important application of these
quasi-static models is to study the flow in hoppers. Figure 12 demonstrated two different flow situations
of hopper flows. It can be seen in Figure 12(a) that discharge happens only within a small channel
directly above the exit. The material outside this channel is stagnant. This kind of flows has been
referred to ‘funnel flow’ in hopper. The discharging channel is called a ‘rat hole’. While in Figure 12(b),
the material inside the hopper flows uniformly without stagnation area. This flow pattern is called ‘mass
flow’. It was found that the motion within the hopper high above the exit slot was highly complex and
determined by conditions at the surface. However, near the exit slot of the hopper, this complex motion
degenerated to a common motion. Early research of hopper flow was done by Jenike (1964, 1965), who
obtained the stress field near the exit using similarity solutions. By adding inertia terms to the
momentum equations, Brennen and Pearce (1978) determined the discharge rate at the bottom exit of
17

a hopper. It was found the discharge rate was independent of the depth of fill and conditions at the
upper surface, given the depth was large compared to the width of the exit slot. However, simulation
results of hopper flows from Potapov and Campbell (1996) suggested that quasi-static models could not
accurately describe hopper flow by showing the internal friction angle 𝜙 was not a constant inside the
hopper. This finding contradicted the basic assumptions of quasi-static models that the material was
always yielding in the hopper. However, the mechanics behind the variation of 𝜙 was not fully
understood.

(a)

(b)
Figure 12. Pictures of different flow profiles within hoppers.(a) Funnel flow. Picture was taken from
https://engineering.purdue.edu/~wassgren/research/gao_image.jpg (b) Mass flow. Picture was taken from
http://jenike.com/files/2012/10/Hopper-flow-pattern-mass-flow.jpg

2.2 Kinetic theory (rapid flows)

We have examined the slow, quasi-static flow of granular materials. Quasi-static flow is an
extreme case of all granular flows, for which the stress is not a function of the shear rate. One important
feature of slow flow is that friction is the dominant force. With increasing the shear rate, the
constitutive relations developed for slow flow is not valid once the inertia force cannot be ignored. In
18

the following part, we will examine the contrasting regime of rapid flow, which corresponds to high
shear rate flows. It will be shown that the stress varies as the square of the shear rate.
To simplify our analysis, we assume the granular material to be uniformly sized smooth spheres.
Instead of moving in narrow rupture layers, particles in rapid flow move freely and independent of its
neighbors. Particle concentration is generally smaller than seen in quasi-static flow and it is assumed
that particles never experience long duration contacts with their neighbors. Instead particles are
assumed to interact by instantaneous collisions. In rapid flow, the velocity of each particle consists of
two parts, the mean velocity of the shear flow and a random component of fluctuation motion. This
particle model is analogous to the gas molecules in the kinetic theory of gas. Thus, granular materials in
the rapid flow regime are also referred to as ‘granular gases’. Similar to the thermal motion of gas
molecules, the random velocity part of granular particles is a crucial parameter in describing the system.
The mean square value of the random velocity is commonly termed as ‘granular temperature’, which is
similar to temperature in a gas and determines the transport rate of mass, momentum and energy.
Granular temperature 𝑇 , in a three-dimensional case, according to Ogawa (1978), is defined as:

𝑇 =
1
3
<𝒖 ′2
>
(4)
Where 𝒖 ′ is the velocity component of random fluctuation, and <> denotes an average in both space
and time.
Within a rapid flow, there are two major ways that granular temperature is generated. Figure
13a shows the first mode, which is called ‘collisional mode’. Because of the existence of a velocity
gradient and the geometry of the particles, random velocities are generated by collisions between
particles even if they are moving with the mean velocities appropriate for their position. This type of
collision converts energy from the mean motion to the granular temperature. The magnitude of
granular temperature generated by this mode is proportional to the velocity gradient. The second way
19

to generate granular temperature is called the ‘streaming mode’, which is demonstrated by Figure 13b.
In this mode, granular temperature is generated when particles moving parallel to velocity gradient. The
local mean velocity difference between their initial and current positions appears as granular
temperature. For example, particles moving upward the velocity gradient produce a negative random
velocity and in the contrary, particles moving downward the velocity gradient produce a positive
random velocity. As a result, the ‘streaming mode’ generates random velocities only in the direction of
the mean flow field. As a result the granular temperature is generally anisotropic. At large solid
concentration, particles cannot move far between collisions and thus cannot pick up large velocity
differences. Thus, granular temperature is mainly generated by the ‘collisional’ mode. While at small
solid concentrations, the ‘streaming’ mode is the dominant mechanism.

20

Figure 13 granular temperature is generated in two major ways: (a) interparticle collisions (b) particle streaming.

As mentioned above, the two mechanisms of generating granular temperature are both results
of a velocity gradient. In another word, granular temperature is a product of a shear flow. The energy
flow pattern in a rapid granular flow is shown in Figure 14. A granular system will be at rest unless there
is work done on the bulk to agitate the system. Such an agitation could be motion of the bulk material
21

relative to the boundaries. Shear work converts kinetic energy of the mean flow to the random velocity
fluctuations which constitute the granular temperature. Because of the inelasticity of the particle
collisions this kinetic energy is dissipated through collisions into heat. Another possible path for energy
dissipation is by interstitial fluid drag. But here, since air is the only interstitial fluid and we are dealing
only with particles considered big and heavy comparing to air, the drag force can be neglected.

Figure 14. The way how energy flows in rapid granular flows.

Another interesting feature of granular temperature is that it may be conducted exactly
analogous to a molecular gas. The granular temperature of a particle may pass to neighboring areas
through collisions or simply carried with the particle as it moves through the material, which correspond
to collisional and streaming transport modes respectively. Campbell and Brennen (1985) observed the
diffusion of this kind of random fluctuation in particle velocity along the direction of decreasing granular
temperature.
Bagnold (1954) is one of the pioneers to study the constitutive relations of rapid flows. Through
the experimental research of suspended wax spheres in a glycerin-water-alcohol mixture, he found that
in a rapid flow, stress varied as the square of shear rate. Consider a situation where inelastic particles of
driving force
shear
motion of
granular
flow
granular
temperature
dissipated
by collision
or friction
into heat
22

diameter 𝑑 𝑝 , and density 𝜌 𝑝 are subjected to a steady and uniform shear flow of shear rate 𝛾 , under
solid concentration 𝜈 . From here forward, 𝛾 no longer refers to an angle, instead it stands for the shear
rate 𝛾 =
𝑑𝑢
𝑑𝑦
where u is the velocity in the 𝑥 -direction. Stresses within the flow can be expressed as:

𝜏 𝑖𝑗
=𝑓 (𝜈 ,𝛾 ,𝑑 𝑝 ,𝜌 𝑝 )
(5)

where 𝑑 𝑝 is the particle diameter and 𝜌 𝑝 is the density of the solid material from which the particles are
made.
Then, through dimensional analysis, it can be shown that the only solution is:
𝜏 𝑖𝑗
=𝑓 𝑖𝑗
(𝜈 )𝜌 𝑝 𝑑 𝑝 2
𝛾 2
(6)
Equation (6) is purely a result from dimensional analysis. It can be justified by how momentum is
exchanged within the flow. Particle momentum is transported in two different mechanism as mentioned
before: streaming mode and collisional mode. In the streaming mode, the momentum flux is
determined by 𝜌 𝑝 and the granular temperature which can be shown to be proportional to 𝑑 𝑝 2
𝛾 2
. While
in particles collisions, the momentum exchange is proportional to 𝜌 𝑝 𝛾 𝑑 𝑝 and the collision rate is
proportional to the relative velocity which is in turn proportional to 𝛾 𝑑 𝑝 .
Though Bagnold’s study provides a simple constitutive relation between stress and strain rate, it
is almost impossible to apply this result to situations where the strain rate tensor is more complex than
simple shear. Though there are various models derived from different theories. The basic mathematical
relations behind these models are almost identical, and result in equations similar to the Naiver-Stokes
equations.
The relevant hydrodynamic equations include balance equations of mass, momentum and
kinetic energy. The mass balance or equation of continuity is defined as:
23

𝜕𝜌
𝜕𝑡
+ ∇∙(𝜌 𝐮 )=0
(7)
The momentum balance equation is defined as:

𝜕 (𝜌 𝐮 )
𝜕𝑡
+∇∙(𝜌 𝐮𝐮 )=−∇∙𝝈 +𝜌 𝐛
(8)
Where 𝝈 is the total stress tensor and includes both collisional and streaming stress. 𝐛 is the
body force vector and 𝐮 is the velocity vector. Here, the density of the granular material 𝜌 , is defined as
𝜌 =𝜌 𝑝 𝜈 , where 𝜌 𝑝 is the density of the particle and 𝜈 is the solid fraction.
The internal energy can be separated into two parts: the transitional kinetic energy and granular
temperature which is the fluctuational kinetic energy. Here, we will only focus on the conservation of
granular temperature. The energy conservation equation is in the form of:

3
2
𝜌 D𝑇 Dt
=−∇∙𝐪 +Φ
̇ −Γ
(9)
Φ
̇ represents the granular temperature generated by shear work and generally has the form
Φ
̇ =𝝈 𝑻 :∇𝐮 . 𝐪 is the granular heat flux which describes the conduction of granular temperature. Γ is
the energy dissipation through inelastic collisions. We assume thermal temperature does not affect
particle properties, so the balance equation of thermal energy is not included.
To close the constitutive relations, expressions for 𝝈 , 𝐪 , and Γ are required. Haff (1983)
presented a heuristic solution of the problem, which was based on dimensional analysis. It was
intuitively correct in physics and simple in mathematics.
More rigorous analysis was derived from the kinetic theory of non-uniform dense gases
(Chapman & Chowling, 1964). Kinetic theory models make analogies between granular materials and
molecular gases. Constitutive relations were derived in similar ways to the Chapman-Enskog procedure
for gas molecules (Chapman and Cowling (1964)). However, there are fundamental differences between
24

the microscopic gas system and the macroscopic granular system. Since in a molecular gas system,
molecules are constantly in thermal fluctuations which balance the external forces. Macroscopic
particles’ thermal fluctuations are negligibly small and maintaining velocity fluctuations requires
external energy input. In order to apply the theory to macroscopic particles, it’s assumed that the
velocities of colliding particles are uncorrelated and independent of positions, which is called the
molecular chaos assumption. Using this assumption, it can be shown that the random velocities follow a
Maxwell-Boltzmann distribution. However, because of this assumption, applicability of kinetic theory is
limited to dilute granular system where particle contacts are instantaneous binary collisions only, which
means at any instant, one particle cannot be in contact with more than one other particle.
Despite the strong analogy between molecular gas and granular system, there are some basic
differences. The most obvious one is whether energy is conserved during collisions. For a molecular gas
system, collisions between molecules are elastic, which means there is no energy loss. The random
motion of gas molecules would maintain at equilibrium. However, in granular system, collisions are
inelastic. Energy is dissipated through interparticle collisions. Granular temperature decays unless the
system is constantly fed by an external energy input.
At collisions, energy is dissipated only through inelastic interparticle contacts, which is generally
described by the coefficient of restitution 𝑒 . In most analyses, a tangential force such as friction at the
point of contact is ignored. For elastic gas molecules, the temperature is determined by
thermodynamics and the distribution function can be calculated using Boltzmann equations. For
inelastic granules, the granular temperature is determined by the balance between dissipation rate and
energy supply rate
Early calculations done by Savage and Jeffery (1981) derived the collisional part of the stress
tensor for rapid granular flows of uniform sized, inelastic, smooth and sphere particles. Since the energy
25

conservation equation (28) was not included in the analysis, an explicit solution was not achieved.
However, they were able to show that the result depends on a dimensionless parameter 𝑆 , which is
defined as:

𝑆 =
𝑑 𝑝 𝛾 √𝑇 (10)
Over most range of solid concentration 𝜈 , it was shown by Campbell (1989) and Lun et al. (1984)
that 𝑆 was of the order of one. It indicates that the random velocity fluctuation is of the same order as
the relative velocity produced by the velocity gradient 𝑑𝛾 . Imposition of such a velocity gradient would
induce anisotropy and inevitably make the Maxwellian assumption break down. The anisotropy was
initially included in Savage and Jeffery (1981)’s study.
After that, Jenkins & Savage (1983) and Lun et al. (1984) were able to solve the constitutive
equations for slightly inelastic particles only, which meant 1−𝑒 2
≪1. It was made possible to use
Chapman-Enskog expansion on a small parameter 𝜖 =(1−𝑒 2
)
1
2
to approximate the collision integral.
Energy conservation equation was included in their analysis. Jenkins & Savage (1983) used a simplified
pair collision model which assumed anisotropy is linear in the direction of velocity gradient. Their
research gave fairly accurate results for cases corresponding to small 𝑆 values, which only occur at small
values of the solid fraction, ν.
Lun et al. (1984) first presented a method to calculate the streaming stress tensor. Despite the
anisotropy found in previous studies, their analysis abandoned the anisotropic collision distribution used
by Savage & Jeffery (1981) and Jenkins & Savage (1983). Instead, it was based on a perturbed
Maxwellian distribution function. Figure 15 shows stresses from Lun (1984) compared with results from
Jenkins & Savage (1983). His result shows a ‘U’ shape. It agrees with previous studies at high
concentration. But as the solid concentration goes to zero, Lun (1984)’s curve asymptote to ∞. In this
26

region, the streaming mode dominates the stress tensor and, as ν→0, the dissipation rate goes to zero
faster than the temperature production rate resulting in large granular temperatures (although this
effect would be mitigated by other dissipation mechanisms such as air drag). Lun’s analysis provides a
more reasonable solution at lower solid concentrations.

Figure 15 Dimensionless shear and normal stresses vesus solid fraction. Solid lines are results from Lun(1984), dash lines are
results from Jenkins & Savage (1983). (Diagrams are taken from Lun (1984)).

Other early calculations include Jenkins and Richman (1985). Instead of using an asymptotic
expansion on 𝜖 , they solved the constitutive relations using a Grad’s moment expansion model by
including higher moments of velocity distribution functions in the balance equations.
Up to this point, our focus has been on models for rapid flow of smooth particles only. In the
real world, granular materials are rough. Friction would bring tangential forces at particle contacts, and
add another way of dissipating energy. Including friction greatly complicates the calculation in the sense
that the slip/no-slip decision would create discontinuity in collision integrals.
27

There have been numerous attempts to incorporate friction into rapid flow model. Pidduck
(1922) was probably the first person presenting a kinetic theory for a dilute gas consisting of rough
molecules. In his model, it was assumed no slip occurred at the points of contact. McCoy et al. (1966)
further developed the theory for dense gases. In order to incorporate friction even approximately into
the models, one generally introduces a roughness coefficient 𝛽 , defined as the ratio of final to initial
tangential surface relative velocities. It can be understood as the coefficient of restitution in the
tangential direction. In this sense, 𝛽 ’s value is restricted to the rage −1≤𝛽 ≤1. When 𝛽 =−1, it
stands for the smooth particles case. When 𝛽 =0, there is no-slip between the particle surface on
recoil. When 𝛽 =1, it is referred to as “perfect rough” particle model, in which case the tangential
surface velocity is reversed with magnitude stay the same. For 𝛽 =±1, kinetic energy conserved.
However, 𝛽 is not an actual friction coefficient. It is an easy way to incorporate friction into various
models since slip condition related to real friction causes a discontinuity in collision integrals. Some
early studies assumed 𝛽 was constant (Lun & Savege, 1987), however, it should be pointed that in a real
situation, 𝛽 is a function of collision geometry (Goldsmith, 1960). Lun (1991) first presented a kinetic
theory for dense, nearly perfectly rough and slightly inelastic granular flows. Kinetic energy was divided
equally into transitional and rotational modes. The mean spin of particles was equal to half the vorticity
of the bulk. The resulting constitutive relations had the similar form as the one of smooth particles.
There are many recent work focusing on kinetic theory, however, not much progress has been
done. Goldhirsch et al. (2005) derived a set of hydrodynamic equations for nearly smooth particles. In
his analysis, transitional and rotational kinetic energy were individually conserved. The totally kinetic
energy was not equally divided into two halves. Kumaran (2004, 2006) presented a set of constitutive
relations including Burnett order terms for rough granular particles. Burnett order terms are referred to
second order terms of velocity and temperature gradients. By including Burnett order terms, the model
was able to incorporate normal stress differences.
28

A lot of studies (Kumaran, 2008; 2009) have focused on applying kinetic theory on chute flow
and dense granular flow at high volume fraction in the range 55%-59%. The problem examined
experimentally by Pouliquen (1999). Simulation study of chute flow (Silbert, et al., 2001) found within
steady dense chute flows, the velocity profile and stress-strain rate relations satisfied Bagnold’s law,
implying kinetic theories may be applied to model dense chute flows. By adopting corrected particles
distribution function (Torquato, 1995) for granular flow at high solid fractions, models derived from
kinetic theories (Kumaran, 2008; 2009) were able to reasonably predict stress and dissipation rate.
However, kinetic theories still could not explain some flow features such as the constant solid fraction
within the flows and the minimum thickness ℎ
𝑠𝑡𝑜𝑝 to maintain a flow. Some studies (Ertas & Halsey,
2002) linked ℎ
𝑠𝑡𝑜𝑝 to the presence of correlations in the flow and thus treated ℎ
𝑠𝑡𝑜𝑝 as a correlation
length scale.
There are several problems with rapid flow theory that should be noticed. First, friction still
cannot be incorporated perfectly into models. Up until now, only ‘perfectly rough’ and ‘nearly smooth’
particles can be handled by the theory. Second, the rapid flow model is based on kinetic theory of gases,
which assumes collisions between gas molecules are perfectly elastic and energy-conserving. When
applied to granular flows, inelastic collision and the resulting energy dissipation needs to be accounted
for and large energy dissipation greatly complicates kinetic theories. In order to solve the system,
particles within the flows are restricted to ‘nearly elastic’ ones, which have coefficient of restitution 𝘀 >
0.9. Third, all transport rates are assumed to be functions of granular temperature, which means it has
a dominant effect over velocity gradient or relative shear motion 𝑑 𝑝 𝛾 . However, as mentioned before,
computer simulation indicating that both mechanism has equally important roles by showing that over
most range of the flow, 𝑆 ≈1. Furthermore, relative shear motions generate anisotropies within the
flow (Lois, Lemaître, & Carlson, 2007). In the last, most granular flows happen at such a large
29

concentration, which particles interact with neighbors so frequently that put molecular chaos
assumption into question.
To finish, I would like to mention some recent studies related to dense granular flow or chute flow.
Other than the kinetic theory approach which was introduced earlier, GdR Midi (2004) presented a
useful model for dense granular flow, in what they called the Inertial regime where Bagnold rheology
applied. The rheology was based on the friction coefficient 𝜇 (𝐼 ) , which was found to be a function of a
dimensionless quantity 𝐼 - the 𝐼𝑛𝑒𝑟𝑡𝑖𝑎𝑙 𝑛𝑢𝑚𝑏𝑒𝑟 ( da Cruz et al, 2005; Pouliquen et al, 2006). The
I𝑛𝑒𝑟𝑡𝑖𝑎𝑙 𝑛𝑢𝑚 𝑏𝑒𝑟 is defined as
𝛾 ̇𝑑
√
𝑃 𝜌 𝑝 , where 𝑃 is the confining pressure or one-third of the trace of the
stress tensor. It should be noted that relations were discovered in situations where volume fraction was
not a fixed parameter. Instead, the pressure 𝑃 was controlled. Silbert et al (2003) and Takagi et al
(2011) pointed out that 𝜇 (𝐼 ) rheology failed to predict non-local behaviors. A non-local behavior is
defined as a behavior that is related to higher order derivatives of strain rate. Pouliquen & Forterre
(2009) developed a non-local 𝜇 (𝐼 ) rheology by integrating the constitutive relations over the entire
flow, while GDRMiDi (2004) and Ertas et al (2002) approach the problem by adopting the Prandtl mixing
length method. It was found their model was able to, at lest approximately, describe a wide range of
flow (from the citical state up to the inertial regime). However, 𝜇 (𝐼 ) rheology is highly limited when the
system becomes dilute. Simulations showed that 𝜇 (𝐼 ) diverged as 𝐼 increased beyond a certain value
(Forterre & Pouliquen, 2008). Non-local effect strongly suggested that complex correlation structures
existed.
Campbell (2011) studied the clusters in dense granular flows through computer simulation. A
cluster is defined as a network of particles, in which any two particles can be connected through
particles in contacts. It was shown clusters existed in dense flows for both elastic and inelastic reasons.
However, a really surprising result showed that the presence of clusters had no effects on stresses. This
30

finding created some troubles since one would normally expect a change in microstructure affects the
flow rheology. Later, this problem will be examined in detail.

2.3 Elastic theory

We have examined quasi-static and rapid flow models of granular materials. They are derived
for two extreme cases of granular flows. One question yet to be answered is what is happening in
between. Jenkins & Askari (1991) first approached the interface problem between quasi-static flow and
a granular gas using a rapid flow model. However, the theoretical results largely disagree with
experimental results. Zhang & Campbell (1992) pointed out that the transition from a solid behavior to
fluid behavior was a yield like process and could be described by Mohr-Coulomb failure criterion so that
the transition was found to be a function of a certain value of stress ratio
𝜏 𝑥𝑦
𝜏 𝑦𝑦
. However, Potapov &
Campbell (1996) and Campbell, Cleary & Hopkins (1995) found that the stress ratio was hardly constant
in simulations of hopper flow and large-scale landslides. Even though, both rapid flow theories and
quasi-static models, predict constant stress ratios that are not a function of shear rate.
Campbell (2002) introduced a new way to fill in the gap between quasi-static and rapid flow
regimes. By incorporating particle stiffness 𝑘 into the model, the entire granular flow was divided into
two global regimes, called the Elastic and Inertial regimes. In the elastic regimes, particles are generally
in constant contact with other particles and forces are transmitted through elastic deformation of the
particles contacts. The Elastic regime is further divided into two sub-regimes, the Elastic-Quasi-static
regime and the Elastic-inertial regime. In the Elastic-Quasi-static regime, which is the same as the old
Quasi-static regime, stresses are independent of shear rate. While in the Elastic-inertial regime, stresses
31

vary linearly with the shear rate. The inertial regimes contain two sub-regimes which are inertial-
collisional regime and inertial-non-collisional regime. They are separated by whether the flow is
dominated by binary collisions only (Collisional) or filled with clusters of particles (non-Collisional).
In the elastic regimes, granular flows are dominated by force chains, quasi linear structures of
particles experiencing large forces. With the help of photoelastic techniques, the structure of force
chains is visualized in Figure 16. Within a force chain the particles deform together and support the load.
Force chains are formed as shear motion bring particles together. Then, as the shear continues, the
chains rotate and eventually break. At small shear rates, when the flow is constrained to a constant
volume, the average force is determined by the degree of compression of the force chains as the
rotating force chains are forced to accommodate to the volume constraint. Force chains are created at a
rate proportional to the shear rate 𝛾 , and break apart after a duration proportional to
1
𝛾 . As a result, the
averaged stresses are independent of shear rate 𝛾 . This is the elastic-quasi-static regime and
encompasses the quasi-static flows described earlier.

Figure 16. A photoelastic picture of force chains within a two dimensional granular flow. Picture is taken from Howell et al
(1999).
32

At large shear rates and large solid concentrations, though force chains are still the dominant
structure inside the flow, inertial effects can no longer be neglected. The particle momentum is
proportional to the shear rate, 𝛾 and the force generated within the chains have the form 𝑎 +𝑏𝛾 .Here,
𝑎 is the force from quasi-static behavior and 𝑏𝛾 is the additional part generated by particles inertia.
Again, the chains were created at a rate proportional to 𝛾 and have a lifetime proportional to 1/𝛾 . So
the average stresses vary linearly with 𝛾 . This is called elastic-inertial regime. This was the new regime
which was observed in the landslide simulation by Campbell, Cleary & Hopkins (1995). In elastic-inertial
regime, stresses generated by particles inertia and elastic deformation are equally important.
In the elastic regime, it was shown that
𝜏 𝑑 𝑝 𝑘 is a proper way to scale stresses. Since 𝜏 ~𝐹 /𝑑 2
and
𝐹 =𝑘𝛿 , where 𝛿 is the particle deformation. The scaled stress can be understood as the ratio of elastic
deformation to particle diameter in the direction of the stress,
𝜏 𝑑 𝑝 𝑘 ~
𝛿 𝑑 . In the inertial regimes, force
chains disappear and flow follows the Bagnold’s Law. As a result, stresses are scaled inertially as
𝜏 𝜌 𝑑 𝑝 2
𝛾 2
.
A dimensionless stiffness parameter 𝑘 ∗
became significant in the theory. 𝑘 ∗
is defined as 𝑘 ∗
=
𝑘 𝜌 𝑑 𝑝 3
𝛾 2
. There are ways to interpret this parameter. Noting that
𝑘 𝜌 𝑑 𝑝 3
𝛾 2
=(
𝜏 𝜌 𝑑 𝑝 2
𝛾 2
)/(
𝜏 𝑑 𝑝 𝑘 ) , it is the ratio of
inertial stress scaling to the elastic stress scaling. Also, 𝑘 ∗
can be understood as (𝑑 𝑝 /𝛿 𝑖 )
2
,where 𝛿 𝑖 is
the deformation induced by particle inertia. In general, at small value of 𝑘 ∗
, elastic effects dominate and
at large 𝑘 ∗
, the system switches to inertial behavior.
Granular flows are examined under controlled stress and controlled volume conditions which
represent two important scenarios in soil mechanics, drained and undrained soil. In undrained soil, the
volume is kept constant by the water, which is incompressible. Even when the soil structure breaks
down, the water pressure continues to support the load and the particles are free to move like in the
33

controlled volume case. In drained soil, water is allowed to flow out before pressure could be built up.
Forces are supported by soil particles, so the soil acts similarly to a controlled stress behavior.

Figure 17. A flow map showing constant volume granular flows divided into four regimes as a function of solid concentration
and dimensionless stiffness. Diagram is taken from Campbell (2002).

Campbell (2002) studied the behaviors of granular flow under controlled volume condition.
Figure 17 is a flow map showing the entire flow divided into four sub regimes by parameter 𝑘 ∗
and solid
concentration 𝜈 . One would assume, at a fixed concentration, by increasing shear rate, the flow would
go through all four regimes, from quasi-static flow at low shear rate, eventually reaching the rapid flow
regime at high shear rate. However, that was found not to be true. Note that simply varying the shear
rate moves one horizontally in the flow map, and there is no transition between the quasi-static and
pure inertial regimes. A transition from elastic-quasi-static to inertial behavior requires force chains to
disappear which cannot be accomplished by simply varying the shear rate. At a large solid fraction, by
reducing shear rate (
𝑘 𝜌 𝑑 𝑝 3
𝛾 2
→∞), particles inertia becomes negligible, causing a shift from elastic-inertial
regime to quasi-static regime. It can be understood as in the 𝛾 →0 limit (slow flow). The dependence of
34

stresses on shear rate was used as an important indicator of the transition. Figure 18 shows the scaled
normal stress as a function of shear rate. The value of scaled normal stress decreases with 𝑘 /(ρd
3
γ
2
)
(increasing 𝑘 or decreasing 𝛾 ) until eventually reaches a constant value, which indicates quasi-static
behavior.

Figure 18 The scaled normal stress 𝜏 𝑦𝑦
, as a funtion of 𝑘 /(𝜌 𝑑 3
𝛾 2
) . Data were taken from a 1000 particle simulation at a fixed
solid concentration 𝜈 =0.6 for three different values ofparticle surface friction. Diagram is taken from Campbell (2002).

The boundary between collisional and non-collisional behaviors can be determined by the ratio
of particle contact time 𝑡 𝑐 to the binary collision time 𝑇 𝑏𝑐
. The binary collision time is defined as 𝑇 𝑏𝑐
=
𝜋 /√
2𝑘 𝑚 −
𝐷 2
𝑚 2
, where 𝑚 is the mass of a particle and 𝐷 is the dashpot coefficient used in a soft-particle
model. (Later, there will be a discussion of computer simulation techniques.) At a small solid
concentration, as
𝑘 𝜌 𝑑 𝑝 3
𝛾 2
→∞, a transition from non-collisional to collisional behavior occurs. This
35

transition should be thought of as in the hard sphere, 𝑘 →∞ limit(𝑇 𝑏𝑐
→0) . Thus, the flow behavior
shifts from non-collisional to collisional.
At small
𝑘 𝜌 𝑑 𝑝 3
𝛾 2
and relatively small solid fractions, increasing the shear rate, the flow enters the
elastic-inertial regime from the inertial regime. In the elastic-inertial regime, force chains and particle
inertia have comparable effect. Since
𝑘 𝜌 𝑑 𝑝 3
𝛾 2
~
𝑘 𝑚 𝛾 2
~
1
𝑇 𝑏𝑐
2
𝛾 2
, in the elastic-inertial regime, the flow
timescale
1
𝛾 is comparable to 𝑇 𝑏𝑐
. Particles are brought together at a rate comparable to a rate pushing
each other apart. So, at small
𝑘 𝜌 𝑑 𝑝 3
𝛾 2
, force chains may form within the flow at concentrations too small
for them to form under static conditions. Ironically, this indicates that rapid flow can only exist at low
shear rates.
Flow behaviors were also found closely related to the solid fraction and the particle surface
friction. Transitions were observed within only 1% change of solid fraction. In another paper, Campbell
(2003) studied the effect of surface friction on the stability of force chains. It was found that dissipation
helps stabilize particle contacts. Force chains are more stable when formed by particles with a larger
surface friction. Figure 18 shows the effect of particle surface friction on a scaled normal stress. Stress
within the flow is larger for material with a larger surface friction. For 𝜇 =0.1, the normal stress almost
disappears at the quasi-static limit. This indicates that the strength and structural integrity of the force
chains are strongly affected by friction.
Since most granular flows are not confined to be a certain volume, but rather under a constant
stress condition like chute flow, landslide and hopper flow. Then the material can expand or contract
locally in order to support an applied load. Campbell (2005) examined the rheology of granular flow
under stress controlled situations. In order to draw the complete flowmap, boundaries of different flow
regimes need to be identified. Unlike controlled volume case, in a controlled stress flow, the stress level
36

cannot be used as an indicator of flow regimes, since stresses are determined by the applied load not
the shear rate. Figure 19 shows the variation of solid concentration 𝜈 with scaled applied normal stress
𝜏 0
𝑑 /𝑘 for various k
∗
=k/(ρd
3
γ
2
) , where 𝜏 0
is the fixed normal stress. The solid concentration 𝜈 can be
used to find the transition from elastic-quasi-static to elastic-inertial regime. When shearing at quasi-
static regime, granular material is at critical state thus 𝜈 doesn’t change with shear rate 𝛾 . Once flow
enters the elastic-inertial regime, 𝜈 deviates from its critical value and becomes dependent on 𝛾 . At
large 𝜏 0
𝑑 /𝑘 , solid concentrations increase beyond the critical state line due to particle compressibility.
Normally, at a larger shear rate 𝛾 (corresponding to a smaller 𝑘 ∗
), the flow deviates from the critical
state at a larger confining stress (need a larger stress to confine the flow to the critical state).

Figure 19 Solid concentration 𝜈 as a function of scaled applied normal stress 𝜏 0
𝑑 /𝑘 for various 𝑘 ∗
=𝑘 /(𝜌 𝑑 3
𝛾 2
) , where 𝜏 0
is the
fixed normal stress. 𝜇 =0.5. Diagram is taken from Campbell (2005).

37

The shear to normal stress ratio
𝜏 𝑥𝑦
𝜏 𝑦𝑦
can also be used to distinguish quasi-static and elastic-
inertial behaviors. The intersection of the sloping lines and horizontal lines in Figure 20 indicates the
transition between elastic-quasi-static and elastic-inertial regimes. In the quasi-static regime, 𝜏 𝑥𝑦
/𝜏 𝑦𝑦

remains constant, while in elastic-inertial regime,
𝜏 𝑥𝑦
𝜏 𝑦𝑦
increases with shear rate. Again, the ratio of actual
contact time to binary collision time 𝑡 𝑐 /𝑇 𝑏𝑐
is an indicator of transition between inertial regime and
elastic-inertial regime.

Figure 20 shear ratio 𝜏 𝑥𝑦
/𝜏 𝑦𝑦
as a function of 𝑘 /(𝜌 𝑑 3
𝛾 2
) for vairous 𝑡 ∗
=𝜏 0
𝑑 /𝑘 , where 𝜏 0
is the fixed applied normal stress.
Data were take for 𝜇 =0.5. Diagram is taken from (Campbell, 2005).

Figure 21 shows a flow map of granular material under constant stress constraint. By increasing
the shear rate from zero, a granular flow under constant applied stress would experience the expected
progression from elastic-quasi-static to elastic inertial regimes and eventually end up in the pure inertial
regimes. One thing to notice is that rapid flow (inertial-collisional) regime is now at 𝑘 /𝜌 𝑑 3
𝛾 2
→0 limit,
38

which corresponds to 𝛾 →∞ (rapid flow). This limit can also be explained as the soft particle limit,
where particle stiffness 𝑘 →0. This may seem contrary to the hard particle assumption in the rapid flow
theory. However, particle stiffness is also in the ordinate 𝜏𝑑 /𝑘 . Reducing 𝑘 with all other quantities fixed
would cause the flow behavior to move diagonally from lower-right to upper-left, avoiding the rapid
flow regime.

Figure 21 A flow map of granular flow under constant normal stress. Diagram is taken from Campbell (2005).

It should be noticed that there are obvious differences between controlled volume and
controlled stress flowmaps. First, under controlled stress, increasing the shear rate 𝛾 at a fixed stress
level 𝜏𝑑 /𝑘 , the flow would experience transitions from quasi-static regime to elastic-inertial regime, and
eventually end up in the pure inertial regime. However, in the controlled volume cases, by simply
varying shear rate, transitions were only observed between quasi-static/elastic-inertial regimes or
elastic-inertial/inertial regimes and a transition between quasi-static/inertial regimes can only occurs by
changing solid the concentration 𝜈 . The fixed volume flow approaches rapid flow regime as 𝑘 /𝜌 𝑑 3
𝛾 2
→
∞, which correspond to the hard sphere limit 𝑘 →∞. While the controlled stress flow approach rapid
39

flow regime as 𝑘 /𝜌 𝑑 3
𝛾 2
→0, which correspond to the high shear rate limit 𝛾 →∞. Second, it had
already been observed in earlier studies (Hvorslev, 1937; Schofield & Wroth, 1968) that controlled stress
flow could support a wide range of stress under the same shear rate with little concentration change
(the critical state). In controlled volume flows, a single stress value was found at certain shear rate and
solid fraction. It was also found that elastic-inertial behavior is more accessible in a controlled stress
case. The Elastic-inertial regime must always appear in between the quasi-static regime and inertial
regime as there must always be a region where elastic and inertial effects are of the same magnitude. In
a controlled stress flow, elastic forces are determined by the applied stress and it’s relatively easy to
generate a comparable inertial force. While in controlled volume flows, elastic forces are generated by
the kinematics of rotating force chains to meet the geometrical constraints. The elastic forces could be
large and it is difficult to generate comparable inertial force. These findings indicate that there are
fundamental differences between fixed concentration and fixed stress flows.
Both quasi-static theory and rapid flow theory are developed from other branches of science:
metal plasticity and kinetic theory of molecular gases. Particle stiffness governs how particles see each
other and thus it’s one of the most important parameters in describing particle contacts. It’s surprising
to see the particle stiffness had been neglected for so long. By incorporating 𝑘 into the model, quasi-
static flow and rapid flow are unified in a single theory and a better view of the entire flow is provided.

2.4 contact problem for two spheres

The Hertz contact model has already been mentioned in the introduction. However, it will be
shown in the computer simulation section, that the frictional models are key to matching simulated and
experimental data. Thus, an examination of the friction modelling is in order. The most detailed such
40

model for two frictional elastic spherical particles’ contact was developed by Mindlin and Deresiewicz
(Mindlin & Deresiewicz, 1953), when a Hertzian model was employed for the direction normal to the
contact point and friction is taken into account in the tangential direction. This contact model deals
differently with increasing and decreasing normal displacements, and the tangential displacement at any
point along the loading path depends on the whole loading history due to the presence of micro slip in
the tangential direction (Johnson K. L., 1985). Microslip occurs when the frictional limit is exceeded on
part, but not all of the contact leading to a case where tangential stresses are relaxed over that portion
but there is no gross tangential motion between the two spheres. Consider two identical spherical
particles in contact, as shown in Figure 22. Let 𝛿 𝑛 and 𝛿 𝑡 denote the normal and tangential
displacements of the centers of the spheres. Assuming there is no torsion or any relative rolling, the
normal and tangential displacements at the center of the contact relative to the center of either spheres
will be equal to 𝛿 𝑛 and 𝛿 𝑡 , respectively. For a general oblique contact, it is impossible to derive a direct
force-displacement relationship. Following the rules set by Mindlin and Deresiewicz, the change in
tangential displacement can be solved with an incremental procedure and the problem is discussed in
two different loading scenarios.
41

Figure 22 the contact between two spherical particles

2.4.1 Constant normal displacement, varying tangential displacement

42

Figure 23 loading curve for the case of constant normal displacement

An incremental procedure is used in general, relating the force-displacement relationship in the
tangential direction of the contact. The actual tangential force 𝑓 𝑡 is calculated based on previous
tangential force 𝑓 𝑡 0
, the change of the tangential displacement (𝛿 𝑡 −𝛿 𝑡 0
) and an incremental stiffness
𝐾 𝑡 .
𝑓 𝑡 =𝑓 𝑡 0
+𝐾 𝑡 (𝛿 𝑡 −𝛿 𝑡 0
) (11)

A typical loading- unloading- reloading curve is shown in Figure 23. A limit on the tangential
force exists as described by Coulomb’s Law 𝑓 𝑡 ≤𝜇 𝑓 𝑛 . Gross sliding happens when 𝑓 𝑡 =𝜇 𝑓 𝑛 throughout
the whole contact area. However, in reality, equality may only occur on part of the contact area and
only on that portion may slipping occur. Under the constant normal displacement case, the contact
43

force in the tangential direction can always be related directly to displacement in the tangential
direction. When a normal force is held constant and a tangential displacement` is applied for the first
time, the tangential stiffness can be calculated as ：

𝐾 𝑡 =𝐾 𝑡 0
(1−
𝑓 𝑡 𝜇 𝑓 𝑛 )
1
3

(12)
Where 𝐾 𝑡 0
=8𝐺 √𝑅 𝑐 𝛿 𝑛 is the initial stiffness. It is a function of shear modulus 𝐺 , radius of the contact
area 𝑅 𝑐 and displacement in the normal direction 𝛿 𝑛 .
When the tangential displacement changes direction and begins to decrease, the unloading
force-displacement relationship follows a different path and a hysteretic behavior is observed. The
tangential stiffness for unloading is:

𝐾 𝑡 =𝐾 𝑡 0
(1−
𝑓 𝑡 𝑇𝑃
−𝑓 𝑡 2𝜇 𝑓 𝑛 )
1
3

(13)
Where 𝑓 𝑡 𝑇𝑃
denotes the tangential contact force at the turning point.
Following successive incremental steps, an unloading curve can be built. Let 𝛿 𝑡 𝑇𝑃
denote the
tangential displacement at the turning point. It can be shown that the point (−𝛿 𝑡 𝑇𝑃
,𝑓 𝑡 𝑇𝑃
) falls on the
unloading path, and beyond that point all previous first loading history has no more influence. Once a
second turning point is encountered, the stiffness constant for reloading can be expressed as ；

𝐾 𝑡 =𝐾 𝑡 0
(1−
𝑓 𝑡 −𝑓 𝑡 𝑇𝑇𝑃
2𝜇 𝑓 𝑛 )
1
3

(14)
Where 𝑓 𝑡 𝑇𝑇𝑃
is the value of tangential force at the second turning point.

44

2.4.2 Varying normal displacement, varying tangential displacement

Figure 24 Loading path for the case: increasing 𝛿 𝑛 , increasing 𝛿 𝑡

When normal displacement and tangential displacement are changing simultaneously at the
contact, the analysis follows the basic idea that constant normal displacement 𝑓 𝑡 −𝛿 𝑡 curves are built
for both previous and current states of the system, then a series of incremental steps are created to
move between those curves. In order to do that, a basic assumption has to be made that any state of
the system, independent of load path, can be reached through a constant normal displacement 𝑓 𝑡 −𝛿 𝑡
curve.
Figure 24 shows the case of increasing normal displacement, increasing tangential displacement.
A two-step procedure is used to solve the displacement problem in the tangential direction. The system
45

goes from the initial state 0 to final state 2 through an intermediate state 1. The normal displacement
𝛿 𝑛 0
,𝛿 𝑛 2
can be computed independently of the tangential displacement problem. The tangential
displacement problem then is solved by first constructing the constant normal displacement 𝑓 𝑡 −𝛿 𝑡
curves though initial and final steps, which correspond to the dotted and solid lines on the graph
respectively. Following the two-step procedure, state 1 is defined as:
𝑓 𝑡 1
=𝑓 𝑡 0
+𝜇 (𝑓 𝑛 2
−𝑓 𝑛 0
) (15)

𝛿 𝑡 1
=𝛿 𝑡 0
+
𝜇 (𝑓 𝑛 2
−𝑓 𝑛 0
)
𝐾 𝑡 0

(16)
Where 𝐾 𝑡 0
is calculated as 𝐾 𝑡 0
=8𝐺 √𝑅 𝑐 𝛿 𝑛 2
.
Then the second step from state 1 to state 2 is characterized by the constant normal
displacement 𝑓 𝑡 −𝛿 𝑡 curve introduced in 2.4.1.

𝐾 𝑡 1
=𝐾 𝑡 0
(1−
𝑓 𝑡 1
𝜇 𝑓 𝑛 2
)
1
3

(17)
𝑓 𝑡 2
=𝑓 𝑡 1
+𝐾 𝑡 1
(𝛿 𝑡 2
−𝛿 𝑡 1
) (18)
Above is the procedure that should be followed for first loading case. The unloading and
reloading 𝑓 𝑡 −𝛿 𝑡 curves can be derived in similar ways. Figure 25 shows the unloading path in the
tangential direction while both normal and tangential displacements are decreasing. On the graph, state
0 denotes the initial state and 3 is the final state. State 1, state 2 are two intermediate states, while
1,2,3 are on the same constant normal displacement 𝑓 𝑡 −𝛿 𝑡 curve and the tangential force at 2 satisfies
𝑓 𝑡 0
=𝑓 𝑡 2
.
46

Figure 25 Loading path for the case: decreasing 𝛿 𝑛 , decreasing 𝛿 𝑡

As system goes from state 0 to state 1, normal displacement decreases (𝛥 𝑓 𝑛 <0). In the
tangential direction:
𝑓 𝑡 1
=𝑓 𝑡 0
−𝜇𝛥 𝑓 𝑛 (19)

𝛿 𝑡 1
=𝛿 𝑡 0
−
𝜇𝛥 𝑓 𝑛 𝐾 𝑡 01

(20)
Where 𝐾 𝑡 01
=8𝐺 √𝑅 𝑐 𝛿 𝑛 2
. In order to proceed, the turning point 𝑓 𝑡 𝑇𝑃 (1)
corresponding to the new
constant normal displacement 𝑓 𝑡 −𝛿 𝑡 curve (solid line on Figure 25) need to be calculated. The new
turning point is defined as:
47

𝑓 𝑡 𝑇𝑃
|
1
=𝑓 𝑡 𝑇𝑃
|
0
+𝜇𝛥 𝑓 𝑛 (21)
Next, the intermediate state 2 and final state 3 can be solved as:

𝐾 𝑡 12
=𝐾 𝑡 01
(1−
𝑓 𝑡 𝑇𝑃
|
0
−𝑓 𝑡 1
2𝜇 𝑓 𝑛 )
1
3

(22)

𝛿 𝑡 2
=𝛿 𝑡 1
+
𝜇𝛥 𝑓 𝑛 𝐾 𝑡 12

(23)

𝐾 𝑡 23
=𝐾 𝑡 12
(1−
𝑓 𝑡 𝑇𝑃
|
0
−𝑓 𝑡 2
2𝜇 𝑓 𝑛 )
1
3

(24)
𝑓 𝑡 3
=𝑓 𝑡 2
+𝐾 𝑡 23
(𝛿 𝑡 3
−𝛿 𝑡 2
) (25)
The remaining reloading case can be treated in a similar way. One thing to notice is that the
above procedure to solve the contact problem of two spheres in the tangential direction is displacement
driven, which means the tangential force is calculated based on change of tangential displacement.
Experimental studies of the contact frictional behavior were also performed (Bowden & Tabor,
2001; Mullier, Tüzün, & Walton, 1991). The results were compared with prediction with the Elastic
Contact Theory (Mindlin & Deresiewicz, 1953). However, Mullier, Tüzün, & Walton (1991) found the
frictional behavior was more complicated than the theory predicted. Their research experimentally
showed a result that when gross sliding happened, the friction coefficient is a function of the normal
load. Pre-compression to a higher normal load appears to reduce the friction coefficient. Some
researchers (Singer, Bolster, Wegand, & Fayeulle, 1990) argued that since the friction force was
proportional to the contact area 𝑓 ~𝐴 , and according to Hertzian model, the contact area 𝐴 ~𝐹 𝑛 2
3
, was a
function of the normal force at the contact. As a result, the friction coefficient 𝜇 =
𝑓 𝐹 𝑛 ~𝐹 𝑛 −
1
3
. However,
Archard (1957) found 𝐴 ~𝐹 𝑛 𝛼 and thus 𝜇 ~𝐹 𝑛 𝛼 −1
, where 𝛼 takes a value between 0.67 and 1 depending
upon the form of the surfaces .This was because the friction coefficient which reflects the dependence
48

of the friction force on the normal load, is closely related to the asperities on the contact surface and
the state of deformation at contact interface. The surface asperities on a rough contact surface are
often in the micron size range. The real contact area is determined by the number of asperities brought
into contact by the normal load. As a result, when two particles are in contact, the friction coefficient
becomes a function of normal load applied and the surface topography. (Archard, 1957; Tüzün, Adams,
& Briscoe, 1988).
It’s of great importance to understand the contact mechanics of particles, as the results may
have significant implications in the simulations of granular flows.

2.5 Simulations

Finally, it would be necessary to review some different simulation techniques used in studying
granular flows. Computer simulation is extremely useful and important in granular flow research due to
the difficulty of making measurements in real granular flows. There are two general types of computer
models, rigid-particle and soft-particle models.
The rigid-particle model (Campbell, 1989; Campbell, 1997) assumes particles are perfectly rigid
and cannot deform at all. The wave speed through such a material is infinite which implies that particle
interactions are instantaneous. As a result, only binary collisions are possible. Time in the system is
updated from collision to collision (a rigid-particle model is also called an event-driven model). In such a
way, the simulation time step is proportional to the collision frequency. Simulations adopting the rigid-
particle model are highly efficient at low solid concentrations where collisions are infrequent. However,
49

at high particle density, this technique fails due to its inability to model long-duration contacts in
structures such as clusters or force chains.
The soft-particle model (Campbell, 2002; Campbell, 2005) which is first introduced by Cundall
(1974) is more realistic in the sense that it allows particle contacts to have finite durations. The particle
contact is described by a spring-dashpot model in the normal direction of the contact, and a slider in the
tangential direction of the contact, which is illustrated in Figure 26. This technique generally assumes
the shape of particle stays undeformed. Instead, an overlap is allowed when two particles come into
contact. The time step of a soft-particle simulation is determined by the particle stiffness (binary
collision time). Each particle’s motion is determined by applying Newton’s second law at every time
step. In order to increase the time step, soft particles with small stiffness are often used. However, a
problem with the soft-particle model is that when system is at large concentrations, the particles may
deform enough to push through the stress singularity that occurs at high concentrations. Typically,
overlap is limited to 1% of particle diameter.
50

Figure 26 A schematic of the spring dash contact model used for soft particles.𝜇 is the friction coefficient, 𝑘 is the spring stiffness
and 𝐷 is the dashpot coefficient. In such a way, the coefficient of restitution is expressed as 𝘀 =𝑒𝑥𝑝 [−𝜋𝐷 /√2𝑚𝑘 −𝐷 2
], where
m is the particle mass.

Although the spring-dashpot model in Figure 26 is the most commonly used model in simulating
behaviors of granular materials, computer models following simplified (Tsuji, Tanaka, & Ishida, 1992;
Walton & Braun, 1986; Walton, 1993; Zhang & Vu-Quoc, 1999) or complete (Zhang & Vu-Quoc, 2007)
contact law developed by Hertz and Mindlin & Deresiewicz (1953) were also used in the past.
Comparisons were conducted between different kinds of models. Some (Elata & Berryman, 1996) argues
that certain types of simplified Mindlin & Deresiewicz models (Tsuji, Tanaka, & Ishida, 1992) and linear
spring model are thermodynamically inconsistent in that they allow energy to be generated at no cost.
However, a more recent and comprehensive study (Di Renzo & Di Maio, 2004) shows that the spring-
dashpot-slider model, simplified Hertz-Mindlin-Deresiewicz and the complete Hertz-Mindlin-Deresiewicz
models are equally accurate in describing the contact behavior of two spheres (Zhang & Vu-Quoc, 2007).
51

Looking back at all the granular flow theories (quasi-static, rapid flow and elastic theories) we
have reviewed, particle stiffness 𝑘 is shown an important parameter in granular system. It’s reasonable
to choose soft-particle model over rigid-particle model in studying granular flows. The soft-particle
model is more capable and accurate in simulating flows especially at relative high concentrations where
force chains are formed.

52

3. Experimental setup

3.1 Apparatus

The major goal of present research is to further study the elastic behavior of granular flows.
Direct experimentation is carried out in an annular shear cell that is similar to the ones in many previous
studies (Savage & Sayed, 1984; Wang & Campbell, 1992; Daniels & Behringer, 2005, 2006). Figure 27
shows the cross section of the annular shear cell used in tests. The whole apparatus consists of two
concentric circular aluminum disks mounted on the same vertical shaft. The lower disk assembly is
mounted on the shaft, so that it rotates with the shaft. The annular trough on this bottom disk has an
inner diameter of 351.4 mm and a width of 142 mm. The depth of the trough is 55.7 mm. The upper disk
is mounted along the shaft through a combined radial and linear ball bearing assembly, which allows the
upper disk to move freely in both rotational and axial direction although the rotation of the upper disk is
restrained mechanically by a stop bar. The center shaft is driven by a 10 HP AC motor through a 10:1
worm gear speed reducer. The speed of the motor is adjusted by computer through a Dura Pulse GS3-
2010 Variable Frequency Drive (VFD).
The apparatus was designed to perform both controlled volume and controlled stress
experiments. In controlled volume cases, the position of the top plate is locked by a collar mounted on
the center shaft. In controlled stress situations, the top plate is allowed to move freely in the axial
direction. The mass of the top plate is 62kg, which is balanced by a set of counterweights. By adjusting
the counterweights, various stresses can be applied. A pair of pneumatic cylinders is used to exert extra
loads when necessary. A Kistler type 9317B three component (𝑥 ,𝑦 ,𝑧 ) force sensor is connected to a
sensor plate amounted at the center of the top surface of the shear cell. The area of the sensor plate is
48.57 𝑐 𝑚 2
, over which stresses are averaged. A displacement transducer is used to measure the vertical
53

position of the top plate thus determine the height of the annular shear cell. A Koyo TRD-N2500-RZWD
encoder connected to the shaft measures the rotational speed of the cell. The data are transmitted to a
computer through a Kistler type 5010 dual mode amplifier and a Measurement Computing’s PCI
DAS6014 data acquisition card.

Figure 27 A schematic of the annular shear cell used in this study.

A small number of the test particles were glued onto the contact rings in order to create the
rough wall boundary conditions (Campbell, 1993). Figure 28 shows the pattern of glued particles for
both top and bottom contact surfaces in the test channel. The left part is half of the top surface, where
Stopping bar
Sensor
Test region
Steel contact plate
Upper disk
Lower disk
Center shaft
Linear ball bearing
assembly
54

each row contains the same number of particles. At the center of the ring, particles are spaced by two
diameters from center to center. The distance between two adjacent particles in different rows is two
particle diameters as well. The pattern follows such a rule and grows from the center row until reaching
edges of the ring. On the bottom surface, glued particles are arranged in similar ways. The only
difference is that center-to-center spacing is four particle diameters in the center row and between
adjacent particles in different rows. The different spacing is the result of computer simulation studies
that showed that centrifugal forces produced a stronger coupling between the particles and the rotating
lower disc. The tighter spacing on the upper surface helped balance this by creating an effectively
rougher surface thus increasing the coupling with the upper surface.

Figure 28 Patterns of glued particles on bottom and top contact rings. All units are in mm.

3.2 Test Materials

55

The choice of possible test material is very limited. Since the elastic nature of granular flow is
the subject of present study, soft spherical particles are preferred. As the particles deform elastically,
they must be softer than the material from which the shear cell is constructed, otherwise the elastic
deformation of the shear cell walls will bias the measurements. This precludes using stiff materials like
glass beads, in favor of softer materials like plastics. They also had to be relatively cheap as these
experiments are hard on the particles and the test materials had to be replaced often. Furthermore, to
get to small values of 𝑘 ∗
=𝑘 /𝜌 𝑑 3
𝛾 2
, particles of large diameter d, are preferred. Four different kinds
of particles were tested. They were of 2 different sizes, 6mm and 8mm in diameter and three different
densities. All the materials are pellets used in Air-soft games. Pictures of the test particles are shown in
Figure 29 and their properties are listed in Table 1. The test particles’ stiffness was measured in a
INSTRON 5567 uniaxial compression test device with an insulated cage within which the temperature
can be carefully controlled. Measurements were done at various temperatures. Figure 30 shows the
compression test result at room temperature for all four kinds of test particles. For all these materials
the contact force can be described by 𝑓 =𝑘 𝑛 𝛿 𝑛 , where 𝑓 is the normal contact force; 𝑘 𝑛 ,𝑛 are two
constants and 𝛿 is the deformation at the contact. Ideally for elastic spheres, at a small deformation 𝛿 ,
the interparticle contact force in the normal direction can be described by a Hertzian Law, 𝑓 =𝑘 𝐻 𝛿 3/2
,
where 𝑘 𝐻 is the Hertzian pre-coeffcient. Notice that 𝑘 𝐻 ~ 𝑘 𝛿 0.5
, where 𝑘 =𝜕𝑓 /𝜕𝛿 is the stiffness. It can
be seen in Figure 4 that the 0.12g 6mm, 0.2g 6mm and 0.27g 8mm particles followed Hertzian behavior.
However, the force-deformation relation for 0.3g 6mm particles followed a different power law which
was found to be 𝑓 =𝑘 1.15
𝛿 1.15
, and as a result 𝑘 1.15
~𝑘 𝛿 0.15
We will see later that this difference in
materials elastic property will be reflected in the final results.

56

(a)

(b)
Figure 29 Particles used in the test. (a) 6mm particles with singe particle mass 0.12g, 0.2g and 0.3g (b) 8mm particles with
single particle mass 0.27g

Table 1 Different particles used in experiment.
Mass of
Particle (𝑔 )
Diameter
(𝑚𝑚 )
Measured
density
(𝑘𝑔 /𝑚 3
)
Friction
coefficient
Behavior at contact
Measured
coefficient
𝑘 𝑚
Manufacturer
0.12 6 1033 0.17
Hertzian
𝑘 𝑛 =𝑘 𝐻 ,𝑛 =1.5
2.63E+08
𝑁 /𝑚 1.5

CHIMEI,
Taiwan
0.2 6 1841 0.26
Hertzian
𝑘 𝑛 =𝑘 𝐻 ,𝑛 =1.5
2.83E+08
𝑁 /𝑚 1.5

CHIMEI,
Taiwan
0.3 6 2593 0.3
Power Law with
𝑘 𝑛 =𝑘 1.15
,𝑛 =
1.15
1.02E+07
𝑁 /𝑚 1.15

CHIMEI,
Taiwan
0.27 8 1034 0.24
Hertzian
𝑘 𝑛 =𝑘 𝐻 ,𝑛 =1.5
1.85E+08
𝑁 /𝑚 1.5

MARUSHIN,
Japan

57

Figure 28 Compression test results for all four kinds of particles at room temperature. 0.12g 6mm, 0.2g 6mm and 0.27g 8mm
particles exhibited Hertzian behavior, 𝑓 =𝑘 𝐻 𝛿 3/2
, while the test result for 0.3g 6mm particles followed a power law with an
exponent of 1.15, , 𝑓 =𝑘 1.15
𝛿 1.15
. 𝑘 𝐻 has the unit of 𝑁 /𝑚 1.5
and 𝑘 1.15
has the unit of 𝑁 /𝑚 1.15
.

Since particle properties are greatly affected by temperature. Compression tests were repeated
by increasing the temperature in increments of 10°C. Tests result of 0.12 g particles is shown in Figure
31. At 30°C, a particles’ force-deformation data could still be described by Hertzian contact model.
However, when temperature exceeds 40°C, the relationship between contact force and deformation
became linear. It’s obvious that at higher temperatures, particles become softer. Similar results were
also observed for other three kinds of particles. Thus, during the experiment, temperature inside the
test cell was carefully monitored and recorded using a type K thermocouple attached to the sensor
plate. The stiffness was interpolated using measured temperature data so that the data could be
corrected for the stiffness changes. If the temperature reached 30 °C, the motor would be stopped and
the experiment would not resume until the test cell cooled – a process that took several hours. This
assured that the stiffness model of the particles did not change. But as at large concentrations and shear
58

rates, the 30 °C limit could be reached after taking only a single data point, this meant the progress of
the experiments was very slow.

(a)

(b)

(c)

(d)

(e)

(f)
59

(g)

(h)
Figure 29 compression tests of 0.12g particles at different temperatures. (a) 30°C (b)40°C (c) 50°C (d) 60°C (e) 70°C (f) 80°C (g)
90°C (h) 100°C. Test data is in blue. Red line is the curve fitting 𝑓 𝑛 =𝑘 𝐻 𝛿 1.5
to the measurement using shown 𝑘 𝐻 ’s value.
Pink line is the curve fitting using a linear relationship 𝑓 𝑛 =𝑐𝛿 with shown value of c.

From Table 1, 0.12g (6mm) and 0.27g (8mm) particles have the same density and it is
reasonable to assume they are made from the same material. However, it is not clear why 8mm
particles have a smaller stiffness (according to Hertz theory, it should be larger). It was also noted during
the experiments that 0.2g (6mm) and 0.3 (6mm) particles have larger friction coefficients than the other
two kinds of particles. Static friction coefficients were roughly measured for all four kinds of particles by
gluing particles of each type to the corners of two identical square plates. One plate was placed on one
another, making only four point-contacts at the tips of 4 pairs of particles. The plates were then tilted
until sliding begins and the static friction coefficient is computed as the tangent of the angle at which
slide commenced.
Another problem is the generation of dust. As the experiment progressed, dust generated due
to surface abrasion would accumulate. While this did not noticeably change the particle properties, the
dust would fill the interstitial area and disturb the test results. This required periodically discarding the
test material and replacing it with fresh particles.

60

3.3 Method and procedures

The basic goal of the tests was to obtain curves of normal and shear stress as functions of shear
rate and solid concentration for both fixed concentration and fixed stress flows.
For fixed concentration tests, a typical experiment started with filling the trough in the lower
disk by a certain mass of particles. Then the upper disk is lowered and locked to a preset height, which
separates the top and bottom surfaces of the trough by six particle’s diameters. It was shown as long
as 𝐻 /𝑅 ≥10, where 𝐻 was the vertical spacing between top and bottom surfaces and 𝑅 was the
particle radius, the vertical spacing didn’t present significant influence on the result (Campbell, 1993).
The test begins by starting the motor. At various shear rates, stresses and displacement of the top disk
were recorded and averaged for a long enough period of time to get a stable result.
For fixed stress experiments, the procedures are similar. It starts with filling the test cell with a
certain amount of particles (corresponding to 60% at 𝐻 =6𝐷 ). Then the top disk was lowered but not
locked in the vertical direction, so the material could expand or contract in order to balance the applied
load. At small stress levels, the load was applied by removing weight from the counterweight. At large
stresses, the load was generated by adding extra weights on top of the upper disk or by using a pair of
pneumatic cylinders to add further downward force on the disk. Stresses and displacement of the upper
disk were recorded at a set of pre-selected values of 𝑘 ∗
. A threshold value for displacement was set that
would stop the motor before the top cleared the wall of the trough and allowed particles to escape.

61

4. Experimental results

4.1 Controlled volume results

In Figure 30, unscaled raw normal stress 𝜏 data from the experiment was illustrated against
shear rate 𝛾 on log-log plots. Lines in red were added on the graph in order to help identify the flow
behavior. At small concentrations, as shown in Figure 30(a), flows of all four particles stayed in the
inertial regime, where 𝜏 ~ 𝛾 2
. By increasing the solid concentration, flows for 0.3g particles entered the
elastic-inertial regime first at 𝜈 =0.54, where 𝜏 ~𝛾 . In Figure 30(b), with an increasing shear rate, the
slope of curve for 0.3g particles changed from 2 to 1, indicating a transition between the inertial to the
elastic inertial regimes. At 𝜈 =0.58 in Figure 30c, while flows of 0.12g particles and 0.27g particles
were in elastic-inertial regime, 0.2g particles exhibited some kind of transitional behavior between
quasi-static and elastic-inertial regimes, which corresponded to a straight line with a slope between 0
and 1 on the plot. A complete transition between quasi-static and elastic-inertial regimes happened for
the flow of 0.3g particles at 𝜈 =0.58. At 𝜈 =0.6, at low shear rates, all four kinds of particles exhibit
quasi-static behaviors that stresses were independent of shear rate. By increasing the shear rate, inertial
effects gradually appeared in the flows of 0.12g, 0.2 and 0.27g particles and eventually these flows enter
the elastic-inertial regime. One thing to be noticed was that transitions for 0.3g (6mm) particles
happened earlier, at smaller concentrations when compared to other three kinds of particles. This would
be expected if the 0.3g particles were more frictional, which agreed with results from static friction
coefficient measurement.
62

(a)

(b)

63

(c)

(d)

(e)
Figure 30 Normal stress 𝜏 as a function of 𝑠 ℎ𝑒𝑎𝑟 𝑟𝑎𝑡𝑒 𝛾 at various concentrations. (a) ν=0.48 (b) ν=0.54 (c) ν=0.56 (d) ν=0.58
(e) ν=0.60.

Since the contact properties are different for different particles, we cannot simply scale the
stress and shear rate as 𝜏𝑑 /𝑘 and 𝑘 /𝜌 𝑑 3
𝛾 2
. From Figure 28, it can be seen that the contact force is
non-linear and may be expresses generically as 𝑓 =𝑘 𝑛 𝛿 𝑛 or 𝑘 =
𝑑𝑓
𝑑𝛿
=𝑛𝑘
𝑛 𝛿 𝑛 −1
so that the stress can
be written 𝜏 ~𝑘 𝑛 𝛿 𝑛 /𝑑 2
. For 0.3g 6mm particles, 𝑛 =1.15,𝑘 𝑛 =𝑘 1.15
, while for the rest three kinds of
particles 𝑛 =1.5,𝑘 𝑛 =𝑘 𝐻 . In order to take the nonlinearity at contacts into account, stress and shear
64

rate are scaled as (𝜏 𝑑 2−𝑛 /𝑘 𝑛 )
1/𝑛 and (𝑘 𝑛 /𝜌 𝑑 4−𝑛 𝛾 2
)
1/𝑛 , where (𝜏 𝑑 2−𝑛 /𝑘 𝑛 )
1/𝑛 ~𝛿 /𝑑 and
(𝑘 𝑛 /𝜌 𝑑 4−𝑛 𝛾 2
)
1/𝑛 ~𝛿 𝐼 /𝑑 , following the prescription of Campbell (2002). Here 𝛿 𝐼 is the particle
deformation generated by particle inertia. Thus dimensional analysis dictates:
𝜏 𝑑 2−𝑛 𝑘 𝑛 =𝑓 (𝜈 ,
𝑘 𝑛 𝜌 𝑑 4−𝑛 𝛾 2
,𝘀 ,𝜇 ,𝑜𝑡 ℎ𝑒𝑟 𝑚𝑎𝑡𝑒𝑟𝑖𝑎 𝑙 𝑝𝑎𝑟𝑎𝑚𝑒𝑡𝑒𝑟𝑠 )
(26)
Where ε is the coefficient of restitution (or some other dimensionless measure of the energy dissipation
at contacts), µ is the particle surface friction coefficient and “other material parameters” can include a
myriad of dimensionless possibilities such as a difference between dynamic and static friction, rolling
friction, the parameters needed to model the velocity dependence of the coefficient of restitution
(Raman, 1918), work modification of the particle surface, (Mullier et al 1991) etc.
Figure 31 shows the relation between elastically scaled normal stress (𝜏 𝑑 2−𝑛 /𝑘 𝑛 )
1/𝑛 and
(𝑘 𝑛 /𝜌 𝑑 4−𝑛 𝛾 2
)
1/𝑛 at different solid concentrations on logarithmic plots. When scaled in this way, quasi-
static behavior appeared as horizontal lines, Elastic-inertial behavior corresponded to lines with a slope
of −1/2, and lines with slope of −1 indicates inertial behavior.
Note that one cannot expect these scalings to completely collapse the data for all four particle
types. At a minimum, the differences in the surface friction prevent that. Furthermore, as pointed out in
Campbell (2002), one should not expect collapse for different powers n, except in the Elastic-Quasistatic
regime. Still it is clear that scaling the data in this fashion, collapses the curves together.

65

(a)

(b)

66

(c)

(d)

(e)
Figure 31 Elastically scaled normal stress(𝜏 𝑑 2−𝑛 /𝑘 𝑛 )
1/𝑛 as a function of (𝑘 𝑛 /𝜌 𝑑 4−𝑛 𝛾 2
)
1/𝑛 at various concentrations.(a) 𝜈 =
0.48 (b) 𝜈 =0.54 (c) 𝜈 =0.56 (d) 𝜈 =0.58 (e) 𝜈 =0.60. Notice that for 0.3g 6mm particles, 𝑛 =1.15,𝑘 𝑛 =𝑘 1.15
, while for
the rest three kinds of particles 𝑛 =1.5,𝑘 𝑛 =𝑘 𝐻 .

4.2 Controlled stress results

As discussed in Campbell (2005), it is difficult to identify the regimes in controlled stress flows.
This is because the applied stress pretty much controls the measured stresses so that the measured
stress changes little with the shear rate. Campbell (2005) used several indirect means to determine the
67

flow regime. One was that the stress ratio 𝜏 𝑥𝑦
/𝜏 𝑦𝑦
is a constant in the Elastic-quasistatic regime and
rises in the Elastic-inertial regime, making the knee in the diagram an indication of the regime transition.
However, for reasons that are unclear the same behavior is not seen in the shear cell results (or for that
matter in simulations of the shear cell). Thus, we are limited in what can be obtained from controlled
stress measurement to giving the change in concentration as a function of applied load and shear rate.
Those results for 0.12g (6mm) particles are shown with inertial scaling in Figure 32. Here 𝜏 is the
applied stress. Note also that each line in the figure corresponds to different values of 𝑘 ∗
. These
experiments were performed by applying a load to the cell, and then changing the shear rate until it
corresponded to the desired value of 𝑘 ∗
, at which point the concentration was recorded. Then the top
loading was changed and the process repeated.
Figure 32 shows that for much of the stress range, the solid fraction does not change but stays
around 𝜈 ~0.6, This behavior corresponds to a state called critical state which was observed in earlier
studies (Hvorslev, 1937; Schofield & Wroth, 1968). Note that at the higher stress levels, there is a slight
increase in the concentration due to the compressibility of the particles – also a characteristic of the
critical state. At lower stress/high-shear levels the concentration drops indicating a transition from the
Elastic-Quasistatic regime to Elastic-inertial regime and eventually to the full inertial regime occur.
Similar behavior was observed in the simulations of Aharonov & Sparks (1999) and this is generally the
behavior expected from Inertial number models (GDR behavior was Midi, 2004) except here 𝜏 is the
68

Figure 32 Controlled stress experiment data for 0.12g (6mm) particles. Solid concentration 𝜈 is plotted against inertially
scaled normal stress 𝜏 /𝜌 𝑑 2
𝛾 2
.

Figure 33 Comparison of controlled stress data. Solid concentration 𝜈 as a function of inertially scaled normal stress
𝜏 /𝜌 𝑑 2
𝛾 2
for different kinds of particles. For each type of particles, 𝑘 ∗
is varied between 7E+06 and 4E+08. The lines
represent least squares fits to the data across the many k*.

69

applied stress which is different from the pressure P in the Inertial Number due to normal stress
difference effects.
The same data for all 4 particle types are shown in Figure 33. For each type of particles, 𝑘 ∗
is
varied between 7E+06 and 4E+08. As with the 0.12g particles in Figure 32, the inertially scaled stress
does a good job of collapsing the data. The largest deviation is seen for the 0.3g particles that, perhaps
because of its larger surface friction, deviates from the critical state at smaller dimensionless stress than
the other three types. This is perhaps because their larger surface friction makes these particles more
dissipative than the others.

70

5. Simulation results

The shear cell results provide a data set that can be used to calibrate Discrete Element Model
computer simulations. Thus, a goal of this project is to see how close the simulations can match the
experimental data using the best measurements of the material properties and then to test the
sensitivity of the results to those various properties. (However, as shall be shown, this was a somewhat
simplistic view of the problem.) To that end, a simulation was constructed that built upon the model
used to make the preliminary design of the experiment. A side view with the outer wall removed is
shown in Figure 34. The simulation was built to be as close to the experiment as possible with top and
bottom walls roughened by gluing particles following the pattern shown in Figure 28 and the side walls
kept smooth. For computational efficiency. only about an eighth of the circumference was modeled with
the ends closed by periodic boundaries. In theory, this acts as if the simulation was periodically repeated
about the circumference. In practice, it means that as a particle leaves one end of the simulation, it
reenters the other with the same relative velocity and position at which it left.

Figure 34 A snapshot of the computer simulation seen through the outer side wall. The side wall has been removed, so that
the particles may be observed. The blue particles are the roughening particles glued to the driving walls. As there is twice the
density of roughening particles on the upper surface, the roughening particles are more easily seen than on the lower surface.

The simulations used the spring-dashpot contact model shown in Figure 26, assuming a Hertzian
spring and a standard frictional slider based on the linearized Mindlin-Deresiewicz model of Tsuji et al
(1992). (Various frictional models will be tested in section A.1 Review on the Inertial Number model.)
71

The simulations were done against the 0.12g (6mm) particle for which we had the most complete
material data (both our own measurements in Table 1 and those in Maladen et al 2011). In addition,
were assumed a coefficient of restitution 𝜖 =0.96. Table 2 shows the material properties measured by
Maladen et al 2011. Note that Maladen’s friction coefficient measurement of 0.073 is significantly
smaller than the 0.17 from our measurements. They measured the friction coefficient by sticking the
particles on the surfaces of two plates and making one plate slide on the other with an angle. The angle
at with the plates can barely slide continuously is the frictional angle. The restitution coefficient was
measured by fixing one particle on a hard plate and letting another particle fall from the top. A tube was
used to guide the falling particle before the collisions and high-speed camera was used to find out the
speed of the particle. (Note that friction with the tube walls may result in an underestimation of the
coefficient of restitution.) In the following simulations, we will use Maladen’s friction coefficient which,
as we shall see is already too large to match much of the data.
Table 2 Parameters measured by Maladen et al 2011
Measured values from Maladen et al 2011
Hardness 1.7x10
8
kg s
-2
m
-1/2

Restitution coefficient 0.96
Friction coefficient 0.073
Density 1.03 g/cm
3

The results using the strict Maladen parameters are shown in Figure 35.
72

(a)

(b)
73

(c)

(d)
74

(e)
Figure 35 Simulation data in comparison with experiment data for 0.12g(6mm) particles at various concentrations. 𝜇 =
0.07,𝜖 =0.96. (a) ν =0.48 (b) ν=0.52 (c) ν=0.56 (d) ν=0.58 (e) ν=0.60.

At lower concentrations which exhibit inertial behavior, (Figure 35a-c) the experimentally
measured stress values were higher than the stresses from the computer simulation. However, at higher
concentrations in the elastic regimes as in Figure 35 d,e, it was the opposite. This disagreement could be
due to the assumed parameters were not accurate. In order to verify that, material property parameters
are adjusted to see if a match could be achieved between computer simulation and experimental data.
The deviations at the lower concentrations can be understood as the result of too much energy
dissipation. For the purpose of argument, these may be imagined as rapid granular flows where the
stresses are produced by a granular temperature. Too much dissipation means a small temperature and
a correspondingly smaller stress. There are two parameters that dissipate energy, the coefficient of
restitution and the friction coefficient. The coefficient of restitution is already large 𝜖 =0.96. But to
check its effects, Figure 36 shows the comparison between experimental results and simulation results
using a somewhat extreme value, 𝜖 =0.99. At 𝜈 =0.48,0.52, and 0.54 shown in Figure 36a-c, one can
75

see that the stresses increased as a result of a larger coefficient of restitution and the simulation results
matched with experimental results. However, at higher concentrations in the Elastic regimes, particles
within the flow are mostly locked into force chains and have little relative velocities with respect to their
neighbors. This is also the explanation for the small effect of coefficient of restitution in the Elastic
regimes observed by Campbell (2002). As a result, increasing the coefficient of restitution had little
effect on flow stresses at 𝜈 =0.58,𝜈 =0.60 seen in Figure 36d,e.

(a)
76

(b)

(c)
77

(d)

(e)
Figure 36 Simulation data in comparison with experiment data for 0.12g(6mm) particles at various concentrations. 𝜇 =
0.07,𝜖 =0.99 were used for all the cases. (a) ν =0.48 (b) ν=0.54 (c) ν=0.56 (d) ν=0.58 (e) ν=0.60.

At 𝜈 = 0.58 and 0.60, the computer simulations demonstrated much larger stresses than
observed in experiments. Since in the elastic regimes, flow behavior is frictional, a larger friction helps
stabilize the structures of force chains and build up stresses within the flows. Reducing the friction
coefficient would result in a decrease in stresses. This may imply that, the value of 𝜇 =0.07 when
78

particles are subject to no normal compression might have changed when particles are under normal
compression. This is consistent with the results of Mullier, T üz ün, & Walton, (1991); Archard, (1957),
Singer et al (1990) who showed that the friction on a Hertzian contace decreases with the the normal
force on the contact. Furthermore, reducing the friction will reduce the dissipation and raise the
stresses in the inertial regimes.
By adjusting the value coefficient of friction in the elastic regimes, we were able to nearly match
the simulation with experiment. However to get a close match required that different friction
coefficients be used at different solid concentrations. In the elastic regimes, a decreasing value of
coefficient of friction was found with increasing solid concentration thus stress level. For 𝜈 =0.48, 𝜈 =
0.52, 𝜈 =0.54, 𝜇 =0.07 was used; for 𝜈 =0.58, 𝜇 =0.06 was used; for 𝜈 =0.60, 𝜇 =0.02 was used.

(a)
79

(b)

(c)
80

(d)

(e)
Figure 37 Simulation data in comparison with experiment data for 0.12g(6mm) particles at various concentrations. A value of
0.99 was used for restitution coefficient. Different values for coefficient of friction were used at different solid concentrations.
(a) ν =0.48, 𝜇 =0.07 (b) ν=0.54, 𝜇 =0.07 (c) ν=0.54, 𝜇 =0.07 (d) ν=0.58, 𝜇 =0.06 (e) ν=0.60, 𝜇 =0.02.

One would wonder if there is any other way to match simulation and experimental data, after
all, ϵ=0.99 is an extreme for restitution coefficient. Plus, from the results above, it looks like by using a
constant friction coefficient at smaller value, μ=0.02 for example, and adjusting the restitution
81

coefficient might be able to match the data at both low and high concentrations, since it will cause the
stresses to decrease at higher concentrations (ν=0.58,ν=0.60) and the stress to increase in the
lower concentrations (ν=0.48,ν=0.52,ν=0.56). And hopefully this may result in a more realistic
restitution coefficient. The results are shown in Figure 38.
At μ=0.02, simulation and experiment were matched by using varying restitution coefficient
for ν=0.48,0.52,0.54. Coefficient of restitution was set to ϵ=0.95 for ν=0.48; ϵ=0.97 for ν=
0.52,0.54. However, one exception prevented this set of parameters to work out eventually. At ν=
0.58, μ=0.02 is too small that the experimental result was in the elastic-inertial regime, while the
computer simulation predicted the flow in the pure inertial regime. Since the data at ν=0.58 is not
sensitive to changing coefficient of restitution, the only way to match experiment and simulation is by
setting friction coefficient to μ=0.06 for ν=0.58. It shows the data could not be matched without
changing the friction coefficient.

(a)
82

(b)

(c)
83

(d)

(e)
Figure 38 Simulation data in comparison with experiment data for 0.12g(6mm) particles at various concentrations. A
constant value of 0.02 was used for friction coefficient. Different values for restitution coefficient were used at different solid
concentrations. (a) ν =0.48, 𝜖 =0.95 (b) ν=0.54, 𝜖 =0.97 (c) ν=0.54, 𝜖 =0.97 (d) ν=0.58, 𝜖 =0.97 (e) ν=0.60, 𝜖 =0.97.

Test cases were also run for a fixed restitution coefficient at ϵ=0.96. A varying friction
coefficient was used across the concentration range to fit the computer simulation to the experiment.
At ϵ=0.96, for ν=0.48,μ=0.03 was used; for ν=0.52 ,0,54 ,μ=0.01 was used; for ν=0.58,μ=
84

0.06 was used; for ν=0.60,μ=0.02 was used. However, because the required friction coefficient
increases between ν=0.54 and ν=0.58, this set of tested parameters didn’t agree with the
predictions and findings of Archard (1957) and Mullier et al (1991) that showed that the friction
coefficient shold decrease with the applied force. To be consistent with Archard and Mullier et al, i.e.
that μ
ν=0.58
≤μ
ν=0.52
, the minimum restitution coeffcient has to be ϵ=0.99. A large restitution
coefficient has to be used in order to make up for the extra energy dissipated by increasing the friction
coefficient in the inertial regimes.

(a)
85

(b)

(c)
86

(d)

(e)
Figure 39 Simulation data in comparison with experiment data for 0.12g(6mm) particles at various concentrations. A
constant value of 0.96 was used for restitution coefficient. Different values for friction coefficient were used at different solid
concentrations. (a) ν =0.48, 𝜇 =0.02 (b) ν=0.54, 𝜇 =0.01 (c) ν=0.54, 𝜇 =0.01 (d) ν=0.58, 𝜇 =0.06 (e) ν=0.60, 𝜇 =0.02.

The results above suggest that in order to match the simulation and experimental data, different
values of the friction coefficient have to be used in the Elastic regimes. The friction coefficient was found
to have to decrease with increased solid concentration. This result agreed with previous studies by
87

Archard (1957), Singer et al (1990) and Mullier et al (1991) that friction coeffcient of Hertizan Contacts
decreased with increased loads. In the Elastic regimes, particles are in constant contacts. The higher the
solid concentration, the harder the particles are pushed together and the higher the stresses.
Furthermore, In the Elastic regimes at ν=0.58,ν=0.60, flow stresses varied by one order of
magnitude with the shear rate. Thus, even experimental data and simulation data matched at low shear
rate end, by using a constant friction coefficient at each concentration, they would eventually deviate as
stresses increase with shear rate. This is again consistent with Archard (1957), Singer et al (1990) and
Mullier et al (1991) as it indicates that the friction coefficient decreases with stress and therefore the
presure on the contact. An attempt was then made to better fit the simulation and experiment using a
friction coefficient that decreases with the flow stress. The results were shown in Figure 40. In Figure
40a, at 𝜈 =0.58, 𝜇 =0.06 was used for the low shear rate end data, 𝜇 =0.02 for moderate shear rate
range data and 𝜇 =0.01 for high shear rate end data. In Figure 40b, at 𝜈 =0.60, 𝜇 =0.02 was used for
the low shear rate end data, 𝜇 =0.01, for the moderate shear rate range data and 𝜇 =0.005 for the
high shear rate end data. By comparing Figure 40a,b to Figure 37d,e, it’s obvious that using a friction
coefficient that decreases with stress would better match the simulation results to the experiment.
88

(a)

(b)
Figure 40 To better fit the simulation and experiment, decreasing friction coefficients with stresses were used. (a) 𝜈 =0.58,
𝜇 =0.06 would fit the low shear rate end data; 𝜇 =0.02 would fit the moderate shear rate range data; 𝜇 =0.01 would fit
the high shear rate end data. (a) 𝜈 =0.60, 𝜇 =0.02 would fit the low shear rate end data; 𝜇 =0.01 would fit the moderate
shear rate range data; 𝜇 =0.005 would fit the high shear rate end data.

In Figure 41, friction coefficient values found by trials in the elastic regimes for ν=0.58,ν=
0.60 was then showed against the measured normal stress on a log-log plot. However the results of
89

Archard (1957), Singer et al (1990) and Mullier et al (1991) indicate that it is the normal force on each
individual contact and not the flow stress that causes the reduction in the friction coeficient, In the
experiments, there is no access to the indivual contact forces, thus the normal stress τ, which is of
course, a force per unit area, was scaled by ν
2
3
to roughly account for the different number of contacts,
per unit area, expected at different solid concentrations. Although only 4 different values of friction
coefficient were used, it was still quite obvious that in the elastic regimes, the friction coefficient
decreases with increasing normal stress. A straight line was fitted through the data and the slope of the
fitted line was found to be -0.75. This value is out of the range predicted by Archard (1957), which is
between -0.32 and 0. However, considering the friction coefficient was found by trials and only 4
different values were used for the fitting, such a crude fit was the best we could get. The fitted
relationship between 𝜇 and 𝜈 2
3
is:

𝜇 =26(
𝜏 𝜈 2
3
)
−0.75

(27)

90

Figure 41 A graph shows that in the elastic regimes friction coefficient μ decreased with increased normal flow stress τ. The
normal flow stress τ is scaled by ν
2
3
which accounts for the different number of contacts at different solid concentrations.

Figure 42 shows the results of using the fitted friction coefficient. The fit is quite good although
not as good as in Figure 40. The latter is to be expected as Figure 40 fitted the friction coefficient on a
point by point basis. But the results do indicate that the friction coefficient decrease with increasing
normal force as in Archard (1957), Singer et al (1990) and Mullier et al (1991), and that including this
effect is important for an accurate match between simulation and experiment. But the results show that
it is possible to fit a simulation to experimental data although one cannot use measured particle
properties. Instead one must go through the laborious process of running many simulations and trying
to fit to the data. Furthermore, it is not clear that these fitted properties are in anyway universal in the
sense that they can be used for anything other than an annular shear cell. It does show that an accurate
match requires complicated material modelling.

91

(a)

(b)
Figure 42. Simulation cases using 𝜇 , whose values are calculated from the fitted curve in Figure 41 are compared with
experimental results. (a) 𝜈 =0.58 (b) 𝜈 =0.60.

92

6. Conclusions and Summary

The goal of this work was to experimentally confirm, as much as possible, the Elastic theory of
granular flows (Campbell, 2002, 2005). And indeed, the data show the major flow regimes, Elastic-
Quasistatic, Elastic-inertial and pure Inertial, and the transitions between them. In particular, they show
the transition between pure Inertial and Elastic-inertial at high shear rates and low concentrations when
the shear forces particles together faster than the elastic forces pull them apart so that force chains
form at low density. Also at large concentrations, we saw the expected transition between Elastic-
Quasistatic and Elastic-inertial. All of this was in concert with the flowmap drawn out by Campbell
(2002) which is shown in Figure 17. Furthermore, the non-linear elastic scaling proposed in that paper
did much to collapse the data even though the contact forces of the various materials exhibited
different power law behaviors. (A complete collapse would not be expected due to differences in
material properties such as the surface friction.) Attempts to study controlled stress flows were
somewhat less successful as peculiarities of the shear cell experiment prohibited distinguishing between
Elastic-inertial and pure inertial flows. However, the results did show the expected critical state behavior
and generally scaled well with an inertially scaled stress.
Computer simulations using Discrete Element Method(DEM) were performed and compared
with experimental results. Through computer simulations, it was found that in the inertial regimes, by
using material properties measured by Maladen et al (2011) and a restitution coefficient of 0.99, the
experiment agreed with simulation’s prediction. While in the elastic regimes, when particles were under
constant compression, it was found in order to match experimental results and computer simulations,
the friction coefficient used in the simulation had to be lowered at higher stress levels. A rough
relationship between the friction coefficient in the elastic regimes and flow stress was shown in
93

equation (27). Previous studies by Archard (1957), Singer et al (1990) and Mullier et al (1991) found that
the coeffcient of friction for a single particle under Hertzian contact decreased with increased load. Our
study indicated that at least in an annular shear cell, the collective flow bahavior of elastic spherical
particles exhibiting the same trend. A more advanced material model is needed to further increase the
accuracy of predictions made by computer simulations.

94

Acknowledgement

The work was supported by National Science Foundation under grant CBET-0828514 for which
the authors are extremely grateful. We would like to express special thanks to Yunpeng Zhang for his
help on the uniaxial compression tests.

95

Appendix

A.1 Review on the Inertial Number model

A reviewer of the paper based on the experimental work, claimed that all this data falls into the
paradigm of the Inertial number models. Clearly the inertial data must. But the real issue is whether the
𝜇 (𝐼 ) models can explain the Elastic regimes for which they do not even contain the necessary
dimensional parameters. So, in response here are the Elastic regime data, expressed in terms of the 𝜇 (𝐼 )
model. There are two parts of the following data that are important. If the 𝜇 (𝐼 ) models are correct, then
(obviously) 𝜇 should be a solely valued function of 𝐼 . Secondly according to the theory, the solid
concentration should also be solely valued function of 𝐼 . As the following data shows, neither is true.
Figure 43 shows the data solely from the Elastic regimes (Elastic-Quasistatic or Elastic-Intertial)
as a function of 𝐼 . As you can see. 𝜇 is independent of 𝐼 and that the 𝜇 (𝐼 ) models fail in the elastic
regimes, even though the stresses are varying inertially in most of these cases.

96

(a)

(b)
0
0.05
0.1
0.15
0.2
0.25
0.3
0.35
0.4
0.45
0.5
0.1 0.12 0.14 0.16 0.18 0.2
Effective Friction µ
inertia number I
ν=0.54 Elastic
ν=0.54-0.3g
0
0.05
0.1
0.15
0.2
0.25
0.3
0.35
0.4
0.45
0.5
0.06 0.08 0.1 0.12 0.14 0.16 0.18 0.2
Effective Friction µ
Inertial number I
ν=0.56 Elastic
ν=0.56-0.2g
ν=0.56-0.3g
97

(c)

(d)
Figure 43 Plots of the effective friction 𝜇 =𝜏 𝑥𝑦
/𝜏 𝑦𝑦
as a function of the Inertial number 𝐼 = 𝛾𝑑 /
√
𝜏 𝑦𝑦
. All of these data are
taken solely from the Elastic regimes (Elastic-Inertial or Elastic Quasistatic). (a) ν=0.54, (b) ν=0.56, (c) ν=0.58, (d) ν=0.60.

Furthermore, the 𝜇 (𝐼 ) models predict that the solid concentration ν is also a function of 𝐼 ,
which means that each concentration, should collapse to a single inertial number 𝐼 . In Figure 44, the
solid concentration 𝜈 was plotted against Inertial number 𝐼 . Yet 𝐼 is seen to vary by an order of
magnitude or more at constant ν, again indicating that ν is not a function of 𝐼 . Note also that the ranges
of 𝐼 overlap, so that there at least four concentrations over a range of 0.04 (a huge variation in the
0
0.05
0.1
0.15
0.2
0.25
0.3
0.35
0.4
0.45
0.5
0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16
Effective Friction µ
Inertial number I
n=0.58 Elastic
ν=0.58-0.12g
ν=0.58-0.2g
ν=0.58-0.3g
ν=0.58-0.27g
0
0.05
0.1
0.15
0.2
0.25
0.3
0.35
0.4
0.45
0.5
0 0.02 0.04 0.06 0.08 0.1
Effective Friction µ
Inertial number I
v=0.6 Elastic
ν=0.6-0.12g
ν=0.6-0.2g
ν=0.6-0.3g
ν=0.6-0.27g
98

granular world) corresponding to the same value of 𝐼 . Once again, the 𝜇 (𝐼 ) models fail in the Elastic
zones.

Figure 44 A graph shows solid concentration 𝜈 against Inertial number 𝐼 .

A.2 Comparison of three different frictional models

All the simulations in 5.Simulation results used the linearized Mindlin-Deresiewicz frictional
model of Tsuji et al (1992). There are more complicated models, some which take hundreds of code
lines to implement, and, in the middle of the above work, it was decided to see if a different frictional
model could account for the behavior that we ultimately ascribed to a varying friction coefficient.
For this we tested three different models from the open literature. In the direction normal to
the particle contacts, all cases assumed Hertzian behavior and they differ only by their behavior in the
tangential direction. The first or “linear” model was the linearized Mindlin-Deresiewicz model of Tsuji et
al (1992).
0.52
0.53
0.54
0.55
0.56
0.57
0.58
0.59
0.6
0.61
0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0.18
solid concentration
Inertia number I
0.12g partiles
0.2 gparticles
0.3g particles
0.27g particles
99

The second model “incslip,” which was first introduced by Walton (1993) adopted a simplified
incremental Mindlin-Deresiewicz. It assumes the stiffness in the direction of the existing friction force is
calculated as:

𝐾 𝑡 =
{

𝐾 𝑡 0
(1−
𝑇 −𝑇 ∗
𝜇 𝐹 𝑁 −𝑇 ∗
)
𝛾 ,𝑓𝑜𝑟 𝑖𝑛𝑐𝑟𝑒𝑎𝑠𝑖𝑛𝑔 𝑇 𝐾 𝑡 0
(1−
𝑇 −𝑇 ∗
𝜇 𝐹 𝑁 +𝑇 ∗
)
𝛾 ,𝑓𝑜𝑟 𝑑𝑒𝑐𝑟𝑒𝑎𝑠𝑖𝑛𝑔 𝑇
(28)
Where 𝐾 𝑡 0
is the initial tangential stiffness, 𝑇 is the current tangential force and 𝑇 ∗
is the tangential
force of the last “turning point.” (A turning point in when the tangential motion has halted and is about
to reverse direction.) For the first loading curve, 𝑇 ∗
is zero. 𝛾 is a constant that is often set to
1
3
for
Mindlin-Deresiewicz contact Law.
For the third model “ScaLUR” (Walton O. R., 2011)(unpublished Walton model, ScaLUR) it
rewrites the original Mindlin-Deresiewicz tangential contact model in incremental form using an explicit
loading-unloading-reloading curve (Tüzün & Walton, 1992; Tomas, 2007). For the first loading, the
tangential force 𝑇 𝐿 and tangential displacement 𝛿 𝐿 are defined as:

𝑇 𝐿 =𝜇 𝐹 𝑁 [1−(1−
𝐾 𝑡 0
𝛿 𝐿 𝛽 𝜇 𝐹 𝑁 )
𝛽 ],𝛿 𝐿 =
𝛽𝜇 𝐹 𝑁 𝐾 𝑡 0
[1−(1−
𝑇 𝐿 𝜇 𝐹 𝑁 )
1
𝛽 ]
(29)
Where 𝛽 =
2
3
. When unloading and reloading happen, the tangential force 𝑇 𝑈 is:

𝑇 𝑈 =2𝜇 𝐹 𝑁 [1−(1−
𝐾 𝑡 0
𝛿 𝑈 2𝛽𝜇 𝐹 𝑁 )
𝛽 ],𝛿 𝑈 =
2𝛽𝜇 𝐹 𝑁 𝐾 𝑡 0
[1−(1−
𝑇 𝑈 2𝜇 𝐹 𝑁 )
1
𝛽 ]
and after the next direction reversal, the reloading force 𝑇 𝑅 is :
(30)

𝑇 𝑅 =2𝜇 𝐹 𝑁 [1−(1−
𝐾 𝑡 0
𝛿 𝑅 2𝛽𝜇 𝐹 𝑁 )
𝛽 ],𝛿 𝑅 =
2𝛽𝜇 𝐹 𝑁 𝐾 𝑡 0
[1−(1−
𝑇 𝑅 𝜇 𝐹 𝑁 )
1
𝛽 ]
(31)
The total tangential force 𝑓 𝑇 and displacement 𝛿 𝑇 are
100

𝑓 𝑇 =𝑇 𝐿 −𝑇 𝑈 +𝑇 𝑅 ,𝛿 𝑇 =𝛿 𝐿 −𝛿 𝑈 +𝛿 𝑅 (32)
When the normal forces changes, the model first uses the explicit force-displacement relations
described above to scale the old tangential force or displacement. There are two basic scenarios: (1) if
the normal force has increased, then keep the tangential force unchanged and find a new displacement
vector producing the same tangential force (this scales the old tangential displacement vector); (2) if the
normal force has decreased, then keep the old tangential displacement vector and find the new
tangential force (this scales down the tangential force). In other words, the scaling is done differently
depending on the whether the normal force is loading or unloading. After dealing with the change in
normal force and scaling the old tangential force or displacement vector, the new tangential
displacement vector is added to the old one and, from them, the new tangential force is calculated. Of
the three models, this is the only one that avoids the concerns of Elata et al (1996) and properly
conserves energy in all situations.
Simulation results from the three different models described above were compared in Figure 43.
All the cases were run for the 0.12g (6mm) particles at different solid concentration with three
coefficients of friction, 𝜇 =0.01,𝜇 =0.1 and 𝜇 =0.5. There are two conclusions to be drawn from
these figures. (a) Even though three models used different approaches to the friction behavior on the
contact of particles, there was little effect on the stresses; (b) Flows at different solid concentrations
reacted differently to an increasing friction. Figure 45a,b show that, at 𝜈 =0.48, 0.52, the normal
stresses decreased with increasing friction. At 𝜈 =0.48, stresses for 𝜇 =0.1 are 45% lower than
stresses for 𝜇 =0.01 and stresses with 𝜇 =0.5 are 10% lower than stresses with 𝜇 =0.1. At 𝜈 =0.52,
stresses for 𝜇 =0.1 are 40% lower than stresses for 𝜇 =0.01 and stresses with 𝜇 =0.5 are 10% lower
than stresses for 𝜇 =0.1.
101

(a)

(b)
102

(c)

(d)*
103

But then at 𝜈 =0.56,0.58 and 0.60, Figure 45c,d,e, increasing the coefficient of friction, not
only caused an increase in the stresses, it also changed the flow regime from inertial to quasi-static. At
𝜈 =0.56, the flows are in the inertial regime for 𝜇 =0.01 and 𝜇 =0.1, while for 𝜇 =0.5, the flow is in
the quasi-static regime with stresses increased by almost two orders of magnitude. At 𝜈 =0.58, for 𝜇 =
0.01, the flow was in the inertial regime. For 𝜇 =0.1, with an increasing shear rate, the flow entered
the elastic-inertial regime from the elastic-quasi-static regime. For 𝜇 =0.5, the flow was in the elastic-
quasi-static regime. At 𝜈 =0.60, for 𝜇 =0.01, flow was in the elastic-inertial regime, while for 𝜇 =
0.1 𝑎𝑛𝑑 𝜇 =0.5, the flows were in the elastic-quasi-static regime. A comparison for the corresponding
shear stresses are shown in Figure 46.

(e)
Figure 45 Normal stress comparison of results from different simulation models with various coefficients of friction. Elastically
scaled normal stress 𝜏 𝑑 1
2/𝑘 𝐻 as a function of 𝑘 𝐻 /𝜌 𝑑 5
2𝛾 2
(a) 𝜈 =0.48 (b) 𝜈 =0.52 (c) 𝜈 =0.56 (d) 𝜈 =0.58 (e) 𝜈 =0.60.
104

(a)

(b)
105

(c)

(d)
106

(e)
Figure 46 Shear stress comparison of results from different simulation models with various coefficients of friction. Elastically
scaled shear stress 𝜏 𝑑 1
2/𝑘 𝐻 as a function of 𝑘 𝐻 /𝜌 𝑑 5
2𝛾 2
(a) 𝜈 =0.48 (b) 𝜈 =0.52 (c) 𝜈 =0.56 (d) 𝜈 =0.58 (e) 𝜈 =0.60.

These results show that no improvement can be attained using complex frictional models. The
spring-dashpot-slider model performs equally well as the models trying to follow the hysteresis behavior
described by Mindlin and Deresiewicz (1953). However, it’s obvious that friction plays a significant role
in determining the bulk behavior of granular flows. At small solid concentrations, when particles are
relatively free to move with respect to each other, friction works as an energy dissipation mechanism so
that larger friction reducing the stresses. At larger solid concentrations, when particles are packed
closely, friction helps develop force chains within the flows so that increasing the friction coefficient not
only increases the stresses but also forces transitions between flow regimes.

107

References

Archard, J. F. (1957). Elastic Deformation and the Laws of Friction. Proc. R. Soc. Lond. A, 243 190-205.
Atkinson, J. H., & Bransby, P. L. (1982). The mechanics of soils. Maidenhead: McGraw-Hill.
Bagnold, R. A. (1954). experiments on a gravity-free dispersion of large solid spheres in a newtonian
fluid under shear. Proceedings of the Royal Society of London, Vol. 225, No. 1160 (Aug. 6, 1954),
pp. 49-63.
Batchelor, G. K. (1967). An introduction to fluid mechanics. cambridge: cambridge university press.
Bowden, F. B., & Tabor, D. (2001). Friction and Lubrication of Solids. Oxford: Oxford University Press.
Brennen, C., & Pearce, J. C. (1978). Granular material flow in two-dimensional hoppers. Trans. ASME E J.
Appl. Mech, 45, 43-50.
Brown, R. L., & Richards, J. C. (1970). Principles of Powder Mechanics. oxford: pergamon.
Campbell, C. (1993). Boundary interactions for two-dimensional granular flows:Part I. Flat boundaries,
asymmetric stresses and couple stresses. J. Fluid Mech., 247,111-136.
Campbell, C. (1997). Self diffusion in granular flows . J. Fluid. Mech., 348, 85-101.
Campbell, C. (2003). A problem related to the stability of force chains. Granular matter , 5, 129-134.
Campbell, C. (2011). Clusters in dense-inertial granular flows. J. Fluids Mech., 687, 341-359.
Campbell, C. S. (1989). The stress tensor for simple shear flows of a granular. J. Fluid Mech, 203, 449-
473.
Campbell, C. S. (2002). Granular shear flows at the elastic limit. J. Fluid Mech., 465, 261-291.
Campbell, C. S. (2005). Stress controlled elastic granular shear flows. J. Fluid Mech., 539, 273-297.
Campbell, C. S., & Brennen, C. E. (1985). Computer simulation of granular shear. J. Fluid Mech, 151, 167-
188.
Campbell, C. S., Cleary, P., & Hopkins, M. A. (1995). Large landslide simulations: global deformation,
velocities and basal friction. J. Geophys. Res., 100, 8267-8283.
Chapman, S., & Chowling, T. G. (1964). the mathematical theory of nonuniform gases. cambridge:
cambrige university press.
Chung, T. J. (1988). Continuum Mechanics. New York: prentice-hall.
108

Cundall, P. (1974). A computer model for rock-mass behavior using interactive graphics for input and
output of geometrical data. U.S. Army Corps of Engineers (Missouri River Division), Tech. Rep.
No. MRD-2074.
da Cruz, F., Emam, S., Prochnow, M., Roux, J.-N., & Chevoir, F. (2005). Rheophysics of dense granular
materials: discrete simulation of plane shear flows. Phys. Rev.E, 72(2) 021309.
Daniels, K., & Behringer, R. (2005). Hysteresis and competition between disorder and crystallization in
sheared and vibrated granular flow. Phys. Rev. Lett., 94, 168001.
Daniels, K., & Behringer, R. (2006). Characterization of a freezing/melting transition in a vibrated and
sheared granular medium. J. Stat. Mech., 7, P07018.
Di Renzo, A., & Di Maio, F. P. (2004). Comparison of contact-force models for the simulation of collisions
in DEM-based granular flow codes. Chemical Engineering Science, 59(3), 525-541.
Drucker, D. C., & Prager, W. (1952). Soil mechanics and plastic analysis or limit design. Q. Appl. Math.,
10, 157-165.
Duran, J. (2010). Sands, Powders and Grains: An introduction to the physics of granular materials. new
york: springer.
Elata, D., & Berryman, J. G. (1996). Contact force-displacement laws and the mechanical behavior of
random packs of identical spheres. Mechanics of Materials, 24(3), 229-240.
Ertas, D., & Halsey, T. (2002). Granular gravitational collapse and chute flow. Europhys., Lett. 60, 931-
935.
Feda, J. (1982). Mechanics of Particulate Materials. Amsterdam: Elsevier.
Forterre, Y., & Pouliquen, O. (2008). Flows of dense granular media. Annu. Rev. Fluid Mech., 40(1) 1–24.
Fung, Y. C. (1977). A first Course in Continuum Mechanics. englewood cliffs: prentican-hall.
Goldhirsch, I., Noskowicz, S. H., & Bar-Lev, O. (2005). Nearly smooth granular gases. Phys, Rev, Lett., 95,
068002-1-4.
Goldsmith, W. (1960). Impact: The Theory and Physical Behavior of Colliding Solids. London: Edward
Arnold Publ.
Haff, P. K. (1983). Grain flow as a fluid-mechanical phenomenon. J. Fluid Mech., 134, 401-420.
Howell, W., Behringer, R., & Veje, C. (1999). Fluctuations in granular media. Chaos 9, 559-572.
Hunter, S. C. (1983). Mechanics of Continuous Media. Chichester: Ellis Horwood.
Hvorslev, M. J. (1937). Über die Festigkeitseigenschaften gestörter bindiger Böden ( on the physical
properties of disturbed cohesive soils). Ingeniorvidenskabelige Skrifter, A45.
109

Jackson, R. (1983). Some mathematical and physical aspects of continuum models. In R. Meyer, Theory
of Dispersed Multiphase Flow (pp. 291-337). New York: Academic Press.
Jenike, A. (1965). Gravity Flow of Frictional-Cohesive Solids: Convergence to Radial Stress Fields. ASME J.
Appl. Mech. , 32, pp. 205–207.
Jenike, A. W. (1961). GRAVITY FLOW OF BULK SOLIDS. Bulletin 108, University of Utah Engineering
experment station.
Jenike, A. W. (1964). Steady Gravity Flow of Frictional-Cohesive Solids in Converging Channels. ASME J.
Appl. Mech., 31, pp.5-11.
Jenkins, J. T., & Askari, E. (1991). Boundary conditions for granular flows: phase interfaces. J. Fluid Mech,
223, 497-508.
Jenkins, J. T., & Richman, M. W. (1985). Grad's 13-moment system for a dense gas of inelastic particles.
Arch. Rat. Mech. Anal., 87, 355-377.
Jenkins, J. T., & Savage, S. B. (1983). A theory for the rapid flow of identical, smooth, nearly elastic,
spherical partucles. J. Fluid Mech., 130, 187-202.
Johnson, K. L. (1985). Contact Mechanics. Cambridge: Cambridge University Press.
Johnson, P. C., Nott, P., & Jackson, R. (1990). Frictional-Collisional equations of motion for particulate
flows and their application to chutes. J. Fluid Mech., 210, 501-535.
Kaza, K. R. (1982). the mechanics of flowing granular materials. Ph.D. thesis, University of Houseton.
Kumaran, V. (2004). Constitutive relations and linear stability of a sheared granular flow. J. Fluid Mech.,
1-42.
Kumaran, V. (2006). The constitutive relation for the granular flow of rough particles, and its application
to the flow down an inclined plane. J. Fluid Mech., 1-42.
Kumaran, V. (2008). Dense granular flow down an inclined plane: from kinetic theory to granular
dynamics. J. Fluid Mech., 599, 121-168.
Kumaran, V. (2009). Dynamics of dense sheared granular flows, PartI: structure and diffusion. J. Fluid
Mech., 632, 109-144.
Kumaran, V. (2009). Dynamics of dense sheared granular flows, PartII: the relative velocity distribution.
J. Fluid Mech., 632, 109-144.
Lois, G., Lemaître, A., & Carlson, J. (2007). Spatial force correlations in granular shear flow. I. Numberical
evidence. Phys. Rev. E, 76(2) 021302.
110

Lun, C. K. (1991). Kinetic theory for granular flow of dense, slightly inelastic, slightly rough spheres. J.
Fluid Mech, 223, 539-559.
Lun, C. K., Savage, S. B., Jeffery, D. J., & Chepurniy, N. (1984). Kinetic theories for granular flow: inelastic
particles in Couette flow and slightly inelastic particles in a general flow field. J. Fluid Mech, 140,
223-256.
Lun, C., & Savege, S. (1987). A simple kinetic theory for granular flow of rough, inelastic, spherical
particles. Trans. ASME J. Appl. Mech, 54:47-53.
Maladen, R., Ding, Y., Umbanhowar, P., Kamor, A., & Goldman, D. (2011). Mechanical models of sandfish
locomotion reveal principles of high performance subsurface sand-swimming. J. R. Soc.
Interface, 8, 1332-1345.
Malvern, L. E. (1969). Introduction to the Mechanics of a Continuous Medium. englewood cliffs: prentice-
Hall.
McCoy, B. J., Sandler, S. I., & Dahler, J. S. (1966). Transport properties of polyatomic fluids. IV. The
kinetic theory of a dense gas of perfectly rough spheres. J. Chem. Phys., 45, 3485-3512.
MiDi, G. (2004). On dense granular flows. Eur. Phys. J. E, 14(4) 341-365.
Mindlin, R. D., & Deresiewicz, H. (1953). Elastic sphere in contact under varying oblique forces. Journal of
Applied Mechanics(20), 327-344.
Mullier, M., Tüzün, U., & Walton, O. R. (1991). A single-particle friction cell for measuring contact
frictional properties of granular materials. Powder Technology , 65(1-3), 61-74.
Ogawa, S. (1978). Multitemperature theory of granular materials. In S. C. Cowin, & M. Satake, Proc. US-
Japan Seminar on continuum Mechanical and Statistical Approaches in the Mechanics of
Granular Materials (pp. 208-217). gakujutsu bunken fukyu-kai.
Pidduck, F. B. (1922). The kinetic theory of a special type of rigid molecule. Proc. Roy. Soc. A, 101, 101-
112.
Poptapov, A. V., & Campbell, C. S. (1996). Computer simulation of hopper flows. Phys, Fluids, 8, 2884-
2894.
Potapov, A., & Campbell, C. (1996). Computer simulation of hopper flows. Phys. Fluids, A Fluid Dyn., 8,
pp 2884–2894.
Pouliquen, O. (1999). Scaling laws in granular flows down rough inclined planes. Phys. Fluids, 542-548.
Pouliquen, O., & Forterre, Y. (2009). A non-local rheology for dense granular flows. Philos. Trans. A.
Math. Phys. Eng. Sci., 367, 5091–5107.
111

Pouliquen, O., Cassar, C., Jop, P., Forterre, Y., & Nicolas, M. (2006). Flow of dense granular material:
towards simple constitutive laws. J. Stat. Mech. Theory Exp., (07) P07020.
Roscoe, K. H. (1970). the influence of strains in soil mechanics. geotech, 20, 129-170.
Roscoe, K. H., Schofield, A. N., & Wroth, C. P. (1958). on the yield of soils. Geotech, 8, 22-53.
Savage, S. B., & Jeffery, D. J. (1981). The stress tensor in a granular flow at high shear rates. J. Fluid
Mech. , 199,177-215.
Savage, S., & Sayed, M. (1984). stresses developed by dry cohesionless granular materials sheared in an
annular shear cell. J. Fluid Mech., 199,177-215.
Schofield, A., & Wroth, C. P. (1968). Critical state soil mechanics. London: McGraw-Hill.
Serrin, J. (1959). Mathematical principles of classical fluid mechanics. . In S. Flugge, Handbuch der Physik
(p. Vol. VIII/1). Berlin: Springer-Verlag.
Silbert, L., Ertas, D., Grest, G., Halsey, T., Levine, D., & Plimpton, S. (2001). Granular flow down an
inclined plane: Bagnold scaling and rheology. Phys. Rev. E, Vol.64 051302.
Silbert, L., Landry, J., & Grest, G. (2003). Granular flow down a rough inclined plane: Transition between
thin and thick piles. Phys. Fluids , 15, 1 .
Singer, L., Bolster, R. N., Wegand, J., & Fayeulle, S. (1990). Hertzian stress contribution to low friction
behavior of thin MoS2 coatings. Appl. Phys. Lett. 57, 995-997.
Sokolovskii, V. V. (1965). statics of granular media. oxford: pergamon.
Takagi, D., McElwaine, J., & Huppert, H. (2011). Shallow granular flows. Phys. Rev. E, 83(3), 031306.
Tomas, J. (2007). Adhesion of ultrafine particles - a micromechanical approach. Chemical Engineering
Science, 62(7), 1997-2010.
Torquato, S. (1995). Nearest-neighbor statistics for packings of hard spheres and disks. Phys. Rev. E,
51,3170-3182.
Tsuji, Y., Tanaka, T., & Ishida, T. (1992). Lagrangian numerical simulation of plug flow of cohesionless
particles in a horizontal pipe. Powder Technology, 71(3), 239-250.
Tüzün, U., & Walton, O. R. (1992). Micromechanical modeling of load-dependent friction in contacts of
elastic spheres. J. Phys. D: Appl. Phys., 25, A44.
Tüzün, U., Adams, M. J., & Briscoe, B. J. (1988). An interface dilation model for the prediction of wall
friction in a particulate bed. Chemical Engineering Science, 43(5), 1083-1098.
Vardoulakis, I., & Aifantis, E. C. (1991). A gradient flow theory of plasticity for granular materials. Acta
Mech, 87, 197-217.
112

Walton, O. R. (1993). Numerical simulation of inclined chute flows of monodisperse, inelastic, frictional
spheres. Mechanics of Materials, 16(1-2), 239-247.
Walton, O. R. (2011). The Scale-Load-Unload-Reload (ScaLUR) algorithm for friction between spheres in
contact . Personal Communication.
Walton, O. R., & Braun, R. L. (1986). Viscosity, granular-temperature, and stress calculations for shearing
assemblies of inelastic, frictional disks. Journal of Rheology, 30, 949-980.
Wang, D., & Campbell, C. (1992). Reynolds' analogy for a shearing granular. J. Fluid Mech., 244, 527-546.
Zhang, X., & Vu-Quoc, L. (1999). An accurate and efficient tangential force–displacement model for
elastic frictional contact in particle-flow simulations. Mechanics of Materials, 31(4), 235-269.
Zhang, X., & Vu-Quoc, L. (2007). An accurate elasto-plastic frictional tangential force-displacement
model for granular-flow simulations: Displacement-driven formulation. Journal of
Computational Physics, 225, 730-752.
Zhang, Y., & Campbell, C. S. (1992). The interface between fluid-like and solid-like behaviour in two-
dimensional granular flows. J. Fluid Mech., 237, 541-568.

Abstract (if available)

Abstract This paper reports annular shear cell measurements of granular flows with an eye towards experimentally confirming the flow regimes laid out in the Elastic theory of granular flows. Tests were carried out on four different kinds of plastic spherical particles for both constant volume flows and constant applied stress flows. In particular, observations were made of the new regime in that model, the Elastic-inertial regime, and the predicted transitions between the Elastic-Inertial and both the Elastic-Quasistatic and pure Inertial regimes. Results from Discrete Element Method(DEM) computer simulation using spring-dashpot-slider particles contact model were compared to experiment results evaluating its performance. The results indicated that the particles’ friction coefficient decrease with increasing normal force, and that including this effect is important for an accurate match between simulation and experiment. It’s possible to fit a simulation to experimental data, although one cannot use measured particle properties. Instead one must go through the laborious process of running many simulations and trying to fit to the data. It shows that an accurate match requires complicated material modelling.

Linked assets

University of Southern California Dissertations and Theses

Conceptually similar

PDF

A Boltzmann model for tracking aerosols in turbulent flows

PDF

Reconstruction and estimation of shear flows using physics-based and data-driven models

PDF

Pulsatile and steady flow in central veins and its impact on right heart-brain hemodynamic coupling in health and disease

PDF

Numerical study of flow characteristics of controlled vortex induced vibrations in cylinders

PDF

Supersonic flow study of a capsule geometry using large-eddy simulation

PDF

Flow and thermal transport at porous interfaces

PDF

Chemical and mechanical deformation of porous media and materials during adsorption and fluid flow

PDF

Particle simulations of colloid thruster beam

PDF

Large eddy simulations of laminar separation bubble flows

PDF

An experimental study of shock wave attenuation

PDF

Pattern generation in stratified wakes

PDF

Continuum modeling of reservoir permeability enhancement and rock degradation during pressurized injection

PDF

Aerodynamics at low Re: separation, reattachment, and control

PDF

Laser manipulation of atomic and molecular flows

PDF

Passive rolling and flapping dynamics of a heaving Λ ﬂyer

PDF

The effect of lattice structure and porosity on thermal conductivity of additively-manufactured porous materials

PDF

Acquisition of human tissue elasticity properties using pressure sensors

PDF

Studies of combustion characteristics of heavy hydrocarbons in simple and complex flows

PDF

Analysis of the existence of long range internal stresses in deformed metals

PDF

Human visual perception of centers of optic flows

Asset Metadata

Creator Guo, Tongtong (author)

Core Title An experimental study of the elastic theory of granular flows

School Viterbi School of Engineering

Degree Doctor of Philosophy

Degree Program Mechanical Engineering

Publication Date 09/25/2017

Defense Date 09/15/2017

Publisher University of Southern California (original), University of Southern California. Libraries (digital)

Tag computer simulation,fluid mechanics,granular flow,granular materials,OAI-PMH Harvest,particles

Language English

Contributor Electronically uploaded by the author (provenance)

Advisor Campbell, Charles (committee chair), Kassner, Michael (committee member), Luhar, Mitul (committee member)

Creator Email guotongtong2003@gmail.com,tongtong@usc.edu

Permanent Link (DOI) https://doi.org/10.25549/usctheses-c40-435676

Unique identifier UC11264117

Identifier etd-GuoTongton-5768.pdf (filename),usctheses-c40-435676 (legacy record id)

Legacy Identifier etd-GuoTongton-5768.pdf

Dmrecord 435676

Document Type Dissertation

Rights Guo, Tongtong

Type texts

Source University of Southern California (contributing entity), University of Southern California Dissertations and Theses (collection)

Access Conditions The author retains rights to his/her dissertation, thesis or other graduate work according to U.S. copyright law. Electronic access is being provided by the USC Libraries in agreement with the a...

Repository Name University of Southern California Digital Library

Repository Location USC Digital Library, University of Southern California, University Park Campus MC 2810, 3434 South Grand Avenue, 2nd Floor, Los Angeles, California 90089-2810, USA