List of coefficients as introduced in the constitutive relations of SKT (standard kinetic theory).
Particulate systems and granular matter display dynamic or static, fluid‐ or solid‐like states, respectively, or both at the same time. The mystery of bridging the gap between the particulate, microscopic state and the macroscopic, continuum description is one of the challenges of modern research. This book chapter gives an overview of recent progress and some new insights about the collective mechanical behavior of granular, deformable particles.
- solid‐fluid granular behavior
- micro‐macro transition
- numerical simulations
Dune migration, landslides, avalanches, and silo instability are a few examples of systems where granular materials play an important role. Furthermore, handling and transport of these materials are central to many industries such as pharmaceutical, agricultural, mining, and construction and pose many open questions to the researchers. In spite of their ubiquity, understanding and predicting the flow behavior of granular materials is still a major challenge for science and industry. Even in a seemingly simple system such as dry sand, the presence of large numbers of internal degrees of freedom leads to highly nonlinear effects making it difficult to relate the microscopic grain‐level properties to the macroscopic bulk behavior.
Granular systems can show properties commonly associated with either solid or liquid. They can behave like a fluid, that is, yielding under an applied shear stress. On the other hand, they can also behave like solids, being able to resist applied stresses without deforming, showing also interesting anisotropic structure (contact‐and force‐networks) [1, 2]. Lucretius (ca. 98–55 B.C.) was among the first ones to recognize this interesting behavior of soil‐like materials, when he wrote “
Due to their microscopic, discrete nature and their interesting macroscopic, bulk behavior response, granular materials are studied using both discrete and continuum mechanics frameworks. In the realm of the discrete approach, several numerical techniques that are able to reproduce the single particle motions with the given micromechanical properties of the grains have been developed. In such an approach, the dynamic behavior is studied by integrating the Newton’s equations of motion for each grain using micromechanical properties and specific interaction law. Following the pioneer work by Goldhirsch [5, 6], several numerical techniques have been developed to obtain continuum fields from discrete particle data.
Using these numerical methods, one can study the flow behavior of the idealized grains, characterized by some specific micromechanical properties, which might not exist in the nature, but is helpful in understanding the underlying physics of their global behavior. In spite of their versatile applicability and benefits, these numerical methods have limitations such as excessive computational requirements, round off or truncation errors, and an intrinsic dynamic that is sometimes not reflecting the experimental reality. On the other hand, continuum models give a macroscopic view to investigate granular material behaviors. Continuum mechanics theories solve the conservation equations for the whole medium, that is, the balance of mass, momentum, and when necessary, energy. Although the balance laws are easily deducible, defining the constitutive relations poses the bigger challenge. The latter relate stresses and strains taking into account the physics of the grain‐grain interaction.
The goal of the present book chapter is to study the constitutive behavior of granular systems using particle, numerical simulations, and micro‐macro transition. In particular, we focus on the different mechanical responses of a granular material in dense and dilute conditions, corresponding to the fluid and solid behaviors, respectively. In order to systematically analyze the influence of some crucial material parameters, which affect the flow behavior, we focus on an idealized material composed of frictionless, spherical particles, in the absence of any interstitial fluids. Moreover, in order to concentrate on the rheology of particulate systems, disregarding boundary effects, we have considered two system setups which allow simulating steady and homogeneous flows.
This chapter is organized as follows. Section 2 introduces the general rheological framework to describe the flow behavior of granular materials. In the same section, we also briefly review some existing granular rheological models. The particle simulations along with micro‐macro transition are introduced in Section 3, where different system setups that are used to study the steady and homogeneous granular flows are shown. Finally, in Section 4, we present a comprehensive comparison of the existing simulation data with frictionless particles in dilute and dense regimes. In the same section, we highlight the effect of various micromechanical properties (coefficient of restitution, polydispersity, and particle stiffness) on the macroscopic fields (stresses and volume fraction). We present a comparison of these results with the theoretical models in two regimes: the kinetic theory in the dilute regime, and a recently proposed generalized rheological model in the dense regime.
2. Granular rheology
2.1. A micromechanical based continuum approach
Despite the fact that granular materials are discontinuous media, their behavior is commonly described by a continuum approach. Continuum mechanics theories solve the conservation equations of the whole medium, that is, the balance of mass, momentum, and when necessary, energy. Although the balance laws are easily deducible, the big challenge is the definition of the constitutive relations, that is, the rheology. The latter captures the macroscopic behavior of the system, incorporating the microscale grain‐grain interaction dynamics.
A granular flow can undergo different behaviors depending on both properties at the particle level and the macroscopic characteristic of the flow (i.e., velocity and concentration). At the microscopic level, each particle is characterized by its shape, dimension, material, and contact properties. For the sake of simplicity, in this chapter an assembly of identical spheres, of diameter
In the framework of continuum mechanics, dimensionless numbers are often introduced in order to describe the material behavior. These dimensionless numbers are defined as the ratio of different time scales or forces, thus signifying the relative dominance of one phenomenon over another.
In the case of granular flows, the macroscopic time scale associated with the shear rate parallel to the flow plays an important role. Then, it is convenient to scale all the quantities using the particle diameter, particle density, and shear rate , so that the dimensionless pressure and stiffness are given as and , respectively. On the other hand, when particle deformability becomes relevant, quantities are usually made dimensionless using the particle stiffness; pressure and shear rate are then expressed as and . In the following sections, we will see how these dimensionless numbers are used to characterize granular flows in their different regimes, namely fluid‐like and solid‐like.
2.2. Continuum models
In the early modeling attempts, granular flow is envisaged as existing in either dense solid‐like or loose gas‐like regimes. Early works using shear cell experiments observed these regimes by varying the shear rate and allowing the bed to dilate or compact. Granular materials exhibit solid‐like behavior if the particles are packed densely enough and a network of persistent contacts develops within the medium, resulting in a jammed mechanically stable structure of the particles. On the other hand, when the grains are widely spaced and free to move in any direction, interacting only through collisions, the medium is unjammed and behaves like a fluid .
In the fluid‐like limit, the system is very dilute and the grains interact mainly through binary, instantaneous, uncorrelated collisions. One of the first rheological models for granular flows in this regime was proposed in 1954 by Bagnold . This empirical model, derived from experiments in two‐dimensional plane shear flows, basically states that the stresses are proportional to the square of the strain rate. This simple law, now known as “Bagnold scaling,” has been the first to understand the physics of granular dynamics at large deformations and has been verified for dry grains in a number of experimental and numerical studies [9–12]. In the fluid‐like regime, the generalization of kinetic theory of granular gases provides a meaningful hydrodynamic description.
On the other hand, when the system is very dense, its response is governed by the enduring contacts among grains, which are involved in force chains; the deformations are extremely slow because the entire network of contacts has to be continuously rearranged (jammed structure). In these conditions, the granular material behaves like a solid, showing an elastic response in which stresses are rate independent. The corresponding flow regime is usually referred to as quasi‐static. Slowly deforming quasi‐static dense granular material has been mainly investigated in the framework of geo‐mechanics. There, the majority of the constitutive models are based on the theories of elasto‐plasticity and visco‐plasticity [13–16], and many of them have been conceived by starting from the well‐known critical state theory [17, 18].
In the transition phase, the grains interact via both force chains and collisions. None of the models cited above is able to deal with this phase‐transition of granular materials from a solid‐like to a fluid‐like state and vice‐versa. Intensive studies of the granular rheology at the phase transition have been conducted in the last decades, for example, by Campbell , Ji and Shen [20, 21], and Chialvo et al.  using 3D simulations of soft frictional spheres at imposed volume fractions. In these works, the authors derived a flow‐map of the various flow regimes and analyzed the transition areas. In particular, they found that, for a collection of particles, the solid‐fluid transition occurs in the limit of zero confining pressure at the critical volume fraction
More recently, new theories have been developed to model the phase transition. The French research group GDR‐MiDi  has suggested that dense granular materials obey a local, phenomenological rheology, known as μ(I)‐rheology, that can be expressed in terms of relations between three nondimensional quantities: volume fraction, shear to normal stress ratio, usually called
Below we present a summary of the two continuum theories that well describe the flow behavior in the limits and their extension to the intermediate regime. Kinetic theory in its standard form (SKT) provides a meaningful hydrodynamic description for frictionless particles in the very dilute regime, while μ(I)‐rheology holds for both frictionless and frictional particles for dense flows. It is important to mention that both theories work only for ideal systems, made of rigid, perfectly elastic, monodisperse particles. Finally, the extension of μ(I)‐rheology to deal with soft and deformable particles is also introduced.
2.2.1. Standard kinetic theory (SKT)
The term “granular gas” is used in analogy with a (classical) molecular gas, where the molecules are widely separated and are free to move in all directions, interacting only through instantaneous, uncorrelated collisions. The main differences between molecular and granular gases are that in the latter case part of the energy is irreversibly lost whenever particles interact and the absence of strong scale separation. These facts have numerous consequences on the rheology of granular gases, one of which being the sizeable normal stress differences .
Analogous to the molecular gases (or liquids), the macroscopic fields velocity and mass density are defined for granular systems . An additional variable of the system, the granular temperature,
Following the statistical mechanics approach, the kinetic theory of granular gases rigorously derives the set of partial differential equations given by the conservation laws of mass, momentum, and energy (the latter describing the time development of the granular temperature) for the dilute gas of inelastically colliding particles.
In this section, we summarize the standard kinetic theory (SKT) for the case of steady and homogeneous flows for a collection of ideal particles, that is, they are rigid, monodisperse, frictionless with diameter,
Further, by substituting the constitutive relations for
SKT was rigorously derived under very restrictive assumptions. In particular, the granular system is assumed to be monodisperse and composed of spherical, frictionless, and rigid particles, interacting only through binary and uncorrelated collisions [7, 40, 41]. Several modifications to the SKT have been introduced in the literature accounting for different effects: interparticle friction [4, 7, 42–44], nonsphericity , or polydispersity . As one example, Jenkins [47, 48] extended the kinetic theory to account for the existence of correlated motion among particles at high concentration.
2.2.2. Traditional µ(I) rheology
A convincing, yet simple phenomenological model that predicts the flow behavior in moderate‐to‐dense regime is the µ(I) rheology. Once again, this rheological law is based on the assumption of homogeneous flow of idealized rigid, monodisperse particles, though the extra constraint of frictionless particles can be dropped. According to this empirical model, only three dimensionless variables are relevant for steady shear flows of granular materials: the volume fraction
To account for the polydispersity of particles, the generalized inertial number taking into account the average diameters of the particles was introduced by . Traditional µ(I) rheology had been successful in describing the flow behavior of homogeneous flows (both dense and fast). But it has failed to capture the slow and nonhomogeneous flow, where a shear rate gradient is present. Researchers have made significant efforts into developing nonlocal models for granular flows .
2.2.3. Soft µ(I) rheology
When particles are not perfectly rigid, instead they have a finite stiffness (or softness), the binary collision time is nonzero and hence presents an additional timescale, which is ignored in the standard inertial number phenomenology. A dimensionless number signifying the finite softness of the particles is the dimensionless pressure , which is needed to describe the flow behavior, as proposed recently in Refs. [30–32].
with the dimensionless pressure
The other dimensionless number needed for the full flow characterization is the volume fraction
3. Numerical simulations
Since a few decades, dynamic particle simulations have been a strong tool to tackle many challenging issues related to understanding the flow behavior of particulate systems.
The molecular dynamics or discrete element methods (DEM) is the term given to the numerical procedure, which is used to simulate assemblies of discrete particles. Molecular dynamics (MD) was originally introduced to simulate the motion of molecules [53–55]. It is essentially the simultaneous numerical solution of Newton’s equation for the motion of individual particles, for which the position, velocity, and acceleration are computed at each time step. Through averaging of positions, velocities, and forces of the particles, the macroscopic fields of the whole system, such as the density, mean velocity, and stresses can be obtained in terms of the micromechanical properties. This helps in revealing insights of the behavior of granular materials, which cannot be captured by experiments. In particular, with MD methods, the role of micromechanical properties of the grains on the macroscopic collective behavior of the system can be analyzed.
Particle simulation methods include three different techniques: The discrete element method (DEM), the event‐driven (ED), and the contact dynamics method (CD). All these methods simulate the inelastic and frictional nature of the contacts among grains through microscopic coefficients (i.e., the coefficients of restitutions and the interparticle friction coefficient). In DEM, deformations of particles during contacts are modeled allowing a finite overlap between grains, whereas in the other two methods, the particles are assumed to be infinitely rigid. Since the results presented in this chapter are obtained by using DEM simulations, below we briefly present an overview of DEM. Readers interested in the latter two methods are referred to Refs. [56–58].
3.1. Discrete element method (DEM)
The discrete element method (DEM) is a family of numerical methods for simulating the motion of large numbers of particles. In DEM, the material is modeled as consisting of finite number of discrete particles, with given micromechanical properties. The interactions between particles are treated as dynamic processes with states of equilibrium developing when the internal forces balance. As previously stated, the granular material is considered as a collection of discrete particles interacting through contact forces. Since the realistic modeling of the deformations of the particles is extremely complicated, the grains are assumed to be nondeformable spheres which are allowed to overlap . The general DEM approach involves three stages: (i) detecting the contacts between elements; (ii) calculating the interaction forces among grains; and (iii) computing the acceleration of each particle by numerical integrating the Newton’s equations of motion while combining all interaction forces. This three‐stage process is repeated until the entire simulation is complete. Based on the fundamental simulation flow, a large variety of modified codes exist and often differ only in terms of the contact model and some techniques used in the interaction force calculations or the contact detection.
In this chapter, we focus on the standard linear spring‐dashpot (LSD) model. Considering two particles,
where and are the normal and tangential component of the overlap at the contact among particle
Collisions may be described using the coefficients of normal and tangential restitution,
3.2. Micro‐macro transition
A research goal in the granular community is to derive macroscopic continuum models based on relevant micromechanical properties. This means to bridge the gap between the microscopic properties and the macroscopic mechanical behavior. The methods and tools for this so‐called micro‐macro transition are often applied to small so‐called representative volume elements (RVEs), where all particles can be assumed to behave similarly. Note that both time‐ and space‐averaging are required to obtain reasonable statistics, the latter being appropriate in the case of steady states.
As previously introduced in Section 2.1, the average velocity of
The stress tensor is of particular interest for the description of any continuum medium. In the case of granular assemblies, previous studies have proposed stress‐force relationships for idealized granular systems that relate average stress in the assembly to fundamental parameters that are explicitly related to statistical averages of inter‐particle load transmission and geometrical arrangement. When referring to a homogeneous volume element
where are the eigenvalues of the stress tensor in Eq. (10). With the development of computational power, nowadays one can simulate reasonable number of particles in a granular system and retrieve good statistical information by micro‐macro procedure. The simulations and coarse‐graining presented in this section were undertaken using the discrete element method (DEM) open‐source code Mercury‐DPM (www.mercurydpm.org).
3.3. Simulation setups
There are two popular ways to extract continuum quantities relevant for flow description such as stress, density, and shear rate from the discrete particle data. The traditional one is ensemble averaging of “microscopic” simulations of homogeneous small samples, a set of independent RVEs. A recently developed alternative is to simulate a nonhomogeneous geometry where dynamic, flowing zones and static, high‐density zones coexist. By using adequate local averaging over equivalent volume (inside which all particles can be assumed to behave similarly), continuum descriptions in a certain parameter range can be obtained from a single simulation.
In Section 4 we will combine results from (a) simple shear RVE and (b) split‐bottom shear cell. The setups are briefly introduced and shown in Figure 1 (see Refs. [30, 49] for more details) and relevant numerical parameters are reported in Table 1. When dimensionless quantities (see Section 2.1) are matched and averaging zones are properly selected, the behaviors from different setups are comparable and a wide flow range can be explored.
3.3.1. Simple shear RVE
The collection of spheres of mean diameter
Variables governing the problem are the volume fraction
Constant‐volume steady simple shear samples are placed in a cuboid box (Figure 1a). The height of the computational domain as
In the case of RVE under constant normal stress condition (Figure 1b), granular systems with polydispersity
3.3.2. Split‐bottom ring shear cell
A common feature of natural slow granular flow is the localization of strain in shear bands, which are typically of few particle diameters width. A specialized geometry proposed recently which allows one to impose an external deformation at constant rate is so‐called split‐bottom geometry (Figure 1c). In this geometry, stable shear bands of arbitrary width can be achieved allowing for a detailed study of microstructure associated with the flow of granular materials in the steady state. Unlike the previous setups, in the split‐bottom geometry, the granular material is not sheared directly from the walls, but from the bottom. The bottom of the setup that supports the weight of material above it is split in two parts, the two parts move relative to each other and creates a wide shear band away from sidewalls. The resulting shear band is robust, as its location exhibits simple and mostly grain independent properties.
In this geometry, due to inhomogeneous flow, granular packings with contrasting properties and behavior coexist, that is, high‐density static to quasi‐static areas and dilated dynamic flowing zones are found in the same system. A superimposed grid meshes the granular bed and averaging is performed within each grid volume. Inside a grid volume all particles are assumed to behave similarly and information for a wide parameter range can be obtained using a single numerical experiment, for example, at increasing pressure levels along the depth of the cell. In the following sections, when presenting data from split‐bottom cell simulations, only grid‐points in the center of the shear band will be considered, where the shear rate
|Parameter/Setup||Symbol||Constant Volume||Constant Pressure||Split‐Bottom Shear Cell|
|Simple Shear||Simple Shear|
|Boundary conditions||Periodic||Periodic||Periodic in azimuthal direction|
|Number of particles||2000||4096||37,000|
|Polydispersity||1||2 and 3||2|
|Coeff. of restitution||0.7||0.8||0.8|
|Dimensionless shear rate||3||10−5–100||10−2–10−5|
4. Rheological flow behavior
In this section, we compare the results from various flow setups discussed above for low‐to‐high volume fractions. We vary various particle and contact properties to understand how the particle micromechanical properties influence the macroscopic flow behavior. We have compared different datasets from different setups and/or authors, and numbered as follows: [A] Peyneau et al. ; [B] Chialvo and Sundaresan ; [C] Shi et al. (unpublished); [D] Singh et al. [30, 63], and [E] Vescovi and Luding . Unless specified, we will only use the data labels in the following discussion for the sake of brevity.
4.1. Influence of coefficient of restitution
Figure 2 presents a data collection from two different setups and plots the dimensionless pressure against volume fraction. It shows data with constant pressure simulations from data [A] together with the constant volume simulation results of data [B], for frictionless monodisperse rigid particles. As expected, the data from the two setups are in good agreement. We observe that the restitution coefficient
For the dilute case, a granular gas with high restitution coefficient, for example,
4.2. Influence of polydispersity
Figure 3 shows the variation of the nondimensional pressure with volume fraction for different polydispersity for constant pressure (data [A] and [C]), constant volume (data [B]) homogeneous shear flow simulations, together with the local shear band data from nonhomogeneous shear flows (data [D]). We observe that for low‐to‐moderate volume fractions, pressure is weakly increasing with volume fraction. The data from different shear setups and different polydispersity collapse and agree with the predictions of SKT. However, for higher volume fractions (
4.3. Effect of particle stiffness
In Figure 4, we show the dimensionless pressure as a function of volume fraction for various values of dimensionless particle stiffness, ranging from 103 to 107. The vertical dashed line shows the monodispersed critical volume fraction as in Figure 3. For the sake of comparison, rigid cases (data [A] and [B]) are also plotted. As expected, for the rigid case, pressure diverges close to the critical volume fraction. For soft particles, the deviation from the rigid case is a function of particle stiffness and depending on the system volume fractions (even for the softest particles the deviation from the rigid limit is small for volume fractions smaller than 0.55). When decreasing the volume fraction below 0.5, all different stiffness data tend to collapse. The solid line is the same standard kinetic theory as in Figure 3 where the assumption of rigid particle breaks down for volume fractions
4.4. Combining both particle stiffness and polydispersity in the dense regime
Figure 5 displays dimensionless pressure plotted against volume fraction for both constant volume (data [E]) and normal stress (data [C]) setups with three polydispersities and dimensionless contact stiffnesses, in the moderate to dense volume fraction regime. Diamonds represent constant volume simulation for monodisperse particles while stars and triangles refer to the constant pressure simulation data for polydispersity 2 and 3, respectively, and different color represent different particle stiffness. For
4.5. From dilute to dense, from “liquid” to “solid,” universal scaling
Figure 6 shows the pressure nondimensionalized in two possible ways (a) using shear rate and (b) using particle stiffness (as introduce in Section 2.1) plotted against the distance from the critical volume fraction for the data from different simulations using frictionless particles. Figure 6a shows a good data collapse for the volume fractions below the critical volume fraction (unjammed regime), or the so‐called fluid regime. In the special case of nearly rigid particles or small confining stress, the scaled pressure diverges at the critical volume fraction, which indicates that the granular fluid composed of rigid particles under shear cannot reach a denser shear jammed state. For the data with softer particles, flow is possible even above the critical volume fraction. For low to moderate volume fractions, the agreement of our data with the rigid case is excellent, while for high volume fractions (especially close to the critical volume fraction) deviations are considerable. The data collapse in the low volume fraction regime shows that the Bagnold scaling relationship between pressure and volume fraction is not strongly affected by particle stiffness, polydispersity, and shear setups, but was influenced by the restitution coefficient (see Figure 2). The “fluid” experiences the energy loss more prominent due to collisions.
For larger volume fractions, the scaling does not collapse the data. Note the deviation between constant volume (data [E]) and constant pressure (data [C]) due to the small difference in the dimensionless stiffness as shown in the legend.
Figure 6b shows the same data but only the soft particle simulations ([C] and [E]) with pressure nondimensionalized by the particle stiffness. In this way, we observe a data collapse for high volume fractions,
4.6. So much for the granular rheology
While up to now, the focus was on understanding the relation between pressure and volume fraction, a granular rheology also must consider the shear stress.
Figure 7 shows the steady state shear stress ratio,
This chapter gives an overview of recent progress in understanding and theoretically describing the collective mechanical behavior of dissipative, deformable particles in different states, both fluid‐like and solid‐like. Particulate systems and granular matter display collisional, dilute and solid, mechanically stable states, either switching forth and back, or both at the same time. In which state the system resides depends not only on material properties like, for example, their discrete nature (elastic stiffness), the dissipation (restitution coefficient) or the size distribution (polydispersity) of the particles, but also on the density of the system and balance between the energy input by (shear) stress or strain‐rate and the energy dissipation by collisions or plastic deformations. Realistic material properties like friction and cohesion as well as nonsphericals particles go beyond the scope of this chapter.
One extreme case of low and moderate density collisional flows (for weak to moderate dissipation and arbitrary polydispersity) is well described by standard kinetic theory (SKT) up to system volume fractions about 0.5, beyond which the elastic behavior of longer‐lasting contacts becomes dominant. Open challenges involve very soft particles for which basic theoretical assumptions of kinetic theory fail, for example, due to multiple contacting particles.
The other extreme case of quasi‐static flow of elastic, mechanically stable solid‐like structures are approximately described by the classical μ(I)‐rheology in the limit of rigid particles, but require a softness correction for comparatively large confining stresses. Remarkably, dissipation, as quantified by the coefficient of restitution, dominates the collisional flows in the dilute regime, while the particle stiffness, the polydispersity, and the friction (data not shown here) are the controlling microparameters for denser quasi‐static and jammed flows.
The mystery of bridging the gap between the collisional, dilute, and the denser quasi‐static, elastic solid‐like regimes is not completely solved yet. The particulate, microscopic states are well understood by particle simulations that via so‐called micro‐macro transition can guide the development of macroscopic, continuum constitutive relations that allow to predict the state and characteristics where a granular system resides in. A unified description that ranges from dilute to dense, from rapid to slow, from soft to rigid, etc., is still one of the great challenges of today’s research.
This chapter provided a few methods and some phenomenology, as well as an overview of recent literature in this field, with theories that can describe the extremes. Various recent works attempted to combine those limit‐cases and provide first combined, generalized theories that go beyond the classical states. However, due to dissipation, friction, cohesion, and nonsphericity of realistic materials, this poses still plenty of challenges for today’s research. Our own ongoing research focuses on providing simple unified/generalized theories, also for systems with attractive forces and with anisotropic microstructures, which were not addressed in this chapter.