Granular material

Static behaviors

Coulomb friction law

Coulomb regarded internal forces between granular particles as a friction process, and proposed the friction law, that the force of friction of solid particles is proportional to the normal pressure between them and the static friction coefficient is greater than the kinetic friction coefficient. He studied the collapse of piles of sand and found empirically two critical angles: the maximal stable angle ${\displaystyle {\theta }_m}$ and the minimum angle of repose ${\displaystyle {\theta }_r}$. When the sandpile slope reaches the maximum stable angle, the sand particles on the surface of the pile begin to fall. The process stops when the surface inclination angle is equal to the angle of repose. The difference between these two angles, ${\displaystyle \mathrm{\Delta }\theta ={\theta }_m-{\theta }_r}$, is the Bagnold angle, which is a measure of the hysteresis of granular materials. This phenomenon is due to the force chains: stress in a granular solid is not distributed uniformly but is conducted away along so-called force chains which are networks of grains resting on one another. Between these chains are regions of low stress whose grains are shielded for the effects of the grains above by vaulting and arching. When the shear stress reaches a certain value, the force chains can break and the particles at the end of the chains on the surface begin to slide. Then, new force chains form until the shear stress is less than the critical value, and so the sandpile maintains a constant angle of repose.

Janssen effect

In 1895, H. A. Janssen discovered that in a vertical cylinder filled with particles, the pressure measured at the base of the cylinder does not depend on the height of the filling, unlike Newtonian fluids at rest which follow Stevin's law for hydrostatic pressure. Janssen suggested a simplified model with the following assumptions:

1. The vertical pressure, ${\displaystyle {\sigma }_{zz}}$, is constant in the horizontal plane;
2. The horizontal pressure, ${\displaystyle {\sigma }_{rr}}$, is proportional to the vertical pressure ${\displaystyle {\sigma }_{zz}}$, where ${\displaystyle K=\frac{{\sigma }_{rr}}{{\sigma }_{zz}}}$ is constant in space;
3. The wall friction static coefficient ${\displaystyle \mu =\frac{{\sigma }_{rz}}{{\sigma }_{rr}}}$ sustains the vertical load at the contact with the wall;
4. The density of the material is constant over all depths.

The pressure in the granular material is then described in a different law, which accounts for saturation:

${\displaystyle p(z)=p_{\mathrm{\infty }}[1-\exp (-z/\lambda )]}$,

where ${\displaystyle \lambda =\frac{R}{2\mu K}}$ and ${\displaystyle R}$ is the radius of the cylinder, and at the top of the silo ${\displaystyle z=0}$.

The given pressure equation does not account for boundary conditions, such as the ratio between the particle size to the radius of the silo. Since the internal stress of the material cannot be measured, Janssen's speculations have not been verified by any direct experiment.

Rowe stress and dilatancy relation

In the early 1960s, engineer Peter Walter Rowe studied dilatancy effect on shear strength in shear tests and proposed a relation between them.

The mechanical properties of assembly of mono-dispersed particles in 2D can be analyzed based on the representative elementary volume, with typical lengths, ${\displaystyle {\ell }_1,{\ell }_2}$, in vertical and horizontal directions respectively. The geometric characteristics of the system is described by ${\displaystyle \alpha =\arctan (\frac{{\ell }_1}{{\ell }_2})}$ and the variable ${\displaystyle \beta }$, which describes the angle when the contact points begin the process of sliding. Denote by ${\displaystyle {\sigma }_{11}}$ the vertical direction, which is the direction of the major principal stress, and by ${\displaystyle {\sigma }_{22}}$ the horizontal direction, which is the direction of the minor principal stress.

Then stress on the boundary can be expressed as the concentrated force borne by individual particles. Under biaxial loading with uniform stress ${\displaystyle {\sigma }_{12}={\sigma }_{21}=0}$ and therefore ${\displaystyle F_{12}=F_{21}=0}$.

At equilibrium state:

${\displaystyle \frac{F_{11}}{F_{22}}=\frac{{\sigma }_{11}{\ell }_2}{{\sigma }_{22}{\ell }_1}=\tan (\theta +\beta )}$,

where ${\displaystyle \theta }$, the friction angle, is the angle between the contact force and the contact normal direction.

${\displaystyle {\theta }_{\mu }}$, which describes the angle that if the tangential force falls within the friction cone the particles would still remain steady. It is determined by the coefficient of friction ${\displaystyle \mu =tg{\varphi }_u}$, so ${\displaystyle \theta \le {\theta }_{\mu }}$. Once stress is applied to the system then ${\displaystyle \theta }$ gradually increases while ${\displaystyle \alpha ,\beta }$ remains unchanged. When ${\displaystyle \theta \ge {\theta }_{\mu }}$ then the particles will begin sliding, resulting in changing the structure of the system and creating new force chains. ${\displaystyle {\mathrm{\Delta }}_1,{\mathrm{\Delta }}_2}$, the horizontal and vertical displacements respectively satisfies

${\displaystyle \frac{\dot{{\mathrm{\Delta }}_2}}{\dot{{\mathrm{\Delta }}_1}}=\frac{\dot{{\epsilon }_{22}}{\ell }_2}{\dot{{\epsilon }_{11}}{\ell }_1}=-\tan \beta }$.

Granular gases

If the granular material is driven harder such that contacts between the grains become highly infrequent, the material enters a gaseous state. Correspondingly, one can define a granular temperature equal to the root mean square of grain velocity fluctuations that is analogous to thermodynamic temperature. Unlike conventional gases, granular materials will tend to cluster and clump due to the dissipative nature of the collisions between grains. This clustering has some interesting consequences. For example, if a partially partitioned box of granular materials is vigorously shaken then grains will over time tend to collect in one of the partitions rather than spread evenly into both partitions as would happen in a conventional gas. This effect, known as the granular Maxwell's demon, does not violate any thermodynamics principles since energy is constantly being lost from the system in the process.

Ulam model

Consider ${\displaystyle N}$ particles, particle ${\displaystyle i}$ having energy ${\displaystyle {\epsilon }_i}$. At some constant rate per unit time, randomly choose two particles ${\displaystyle i,j}$ with energies ${\displaystyle {\epsilon }_i,{\epsilon }_j}$ and compute the sum ${\displaystyle {\epsilon }_i+{\epsilon }_j}$. Now, randomly distribute the total energy between the two particles: choose randomly ${\displaystyle z\in \left[0,1\right]}$ so that the first particle, after the collision, has energy ${\displaystyle z\left({\epsilon }_i+{\epsilon }_j\right)}$, and the second ${\displaystyle \left(1-z\right)\left({\epsilon }_i+{\epsilon }_j\right)}$.

The stochastic evolution equation:

${\displaystyle {\epsilon }_i(t+dt)=\{\begin{array}{ll}{\epsilon }_i(t)& probability:1-\mathrm{\Gamma }dt\\ z\left({\epsilon }_i(t)+{\epsilon }_j(t)\right)& probability:\mathrm{\Gamma }dt\end{array}}$,

where ${\displaystyle \mathrm{\Gamma }}$ is the collision rate, ${\displaystyle z}$ is randomly picked from ${\displaystyle \left[0,1\right]}$ (uniform distribution) and j is an index also randomly chosen from a uniform distribution. The average energy per particle: ${\displaystyle \begin{array}{cc}⟨\epsilon (t+dt)⟩& =\left(1-\mathrm{\Gamma }dt\right)⟨\epsilon (t)⟩+\mathrm{\Gamma }dt\cdot ⟨z⟩\left(⟨{\epsilon }_i⟩+⟨{\epsilon }_j⟩\right)\\ & =\left(1-\mathrm{\Gamma }dt\right)⟨\epsilon (t)⟩+\mathrm{\Gamma }dt\cdot {\displaystyle \frac{1}{2}}\left(⟨\epsilon (t)⟩+⟨\epsilon (t)⟩\right)\\ & =⟨\epsilon (t)⟩\end{array}}$.

The second moment:

${\displaystyle \begin{array}{cc}⟨{\epsilon }^2(t+dt)⟩& =\left(1-\mathrm{\Gamma }dt\right)⟨{\epsilon }^2(t)⟩+\mathrm{\Gamma }dt\cdot ⟨z^2⟩⟨{\epsilon }_i^2+2{\epsilon }_i{\epsilon }_j+{\epsilon }_j^2⟩\\ & =\left(1-\mathrm{\Gamma }dt\right)⟨{\epsilon }^2(t)⟩+\mathrm{\Gamma }dt\cdot {\displaystyle \frac{1}{3}}\left(2⟨{\epsilon }^2(t)⟩+2{⟨\epsilon (t)⟩}^2\right)\end{array}}$.

Now the time derivative of the second moment:

${\displaystyle {\displaystyle \frac{d⟨{\epsilon }^2⟩}{dt}}=li{m}_{dt\to 0}{\displaystyle \frac{⟨{\epsilon }^2(t+dt)⟩-⟨{\epsilon }^2(t)⟩}{dt}}=-{\displaystyle \frac{\mathrm{\Gamma }}{3}}⟨{\epsilon }^2⟩+{\displaystyle \frac{2\mathrm{\Gamma }}{3}}{⟨\epsilon ⟩}^2}$.

In steady state:

${\displaystyle {\displaystyle \frac{d⟨{\epsilon }^2⟩}{dt}}=0\Rightarrow ⟨{\epsilon }^2⟩=2{⟨\epsilon ⟩}^2}$.

Solving the differential equation for the second moment:

${\displaystyle ⟨{\epsilon }^2⟩-2{⟨\epsilon ⟩}^2=\left(⟨{\epsilon }^2(0)⟩-2{⟨\epsilon (0)⟩}^2\right)e^{-\frac{\mathrm{\Gamma }}{3}t}}$.

However, instead of characterizing the moments, we can analytically solve the energy distribution, from the moment generating function. Consider the Laplace transform: ${\displaystyle g(\lambda )=⟨e^{-\lambda \epsilon }⟩={\int }_0^{\mathrm{\infty }}{e}^{-\lambda \epsilon }\rho (\epsilon )d\epsilon }$,

where ${\displaystyle g(0)=1}$, and ${\displaystyle {\displaystyle \frac{dg}{d\lambda }}=-{\int }_0^{\mathrm{\infty }}\epsilon e^{-\lambda \epsilon }\rho (\epsilon )d\epsilon =-⟨\epsilon ⟩}$.

the n derivative:

${\displaystyle {\displaystyle \frac{d^{n}g}{d{\lambda }^n}}={\left(-1\right)}^n{\int }_0^{\mathrm{\infty }}{\epsilon }^n{e}^{-\lambda \epsilon }\rho (\epsilon )d\epsilon =⟨{\epsilon }^n⟩}$,

now:

${\displaystyle e^{-\lambda {\epsilon }_i(t+dt)}=\{\begin{array}{ll}{e}^{-\lambda {\epsilon }_i(t)}& 1-\mathrm{\Gamma }t\\ e^{-\lambda z\left({\epsilon }_i(t)+{\epsilon }_j(t)\right)}& \mathrm{\Gamma }t\end{array}}$

${\displaystyle ⟨e^{-\lambda \epsilon \left(t+dt\right)}⟩=\left(1-\mathrm{\Gamma }dt\right)⟨e^{-\lambda {\epsilon }_i(t)}⟩+\mathrm{\Gamma }dt⟨e^{-\lambda z\left({\epsilon }_i(t)+{\epsilon }_j(t)\right)}⟩}$

${\displaystyle g\left(\lambda ,t+dt\right)=\left(1-\mathrm{\Gamma }dt\right)g\left(\lambda ,t\right)+\mathrm{\Gamma }dt{\int }_0^1\underset{=g^2(\lambda z,t)}{\underbrace{⟨e^{-\lambda z{\epsilon }_i(t)}⟩⟨e^{-\lambda z{\epsilon }_j(t)}⟩}}dz}$.

Solving for ${\displaystyle g(\lambda )}$ with change of variables ${\displaystyle \delta =\lambda z}$:

${\displaystyle \lambda g(\lambda )={\int }_0^{\lambda }{g}^2(\delta )d\delta \Rightarrow \lambda g^{\prime }(\lambda )+g(\lambda )=g^2(\lambda )\Rightarrow g(\lambda )={\displaystyle \frac{1}{\lambda T+1}}}$.

We will show that ${\displaystyle \rho (\epsilon )={\displaystyle \frac{1}{T}}{e}^{-\frac{\epsilon }{T}}}$ (Boltzmann Distribution) by taking its Laplace transform and calculate the generating function:

${\displaystyle {\int }_0^{\mathrm{\infty }}{\displaystyle \frac{1}{T}}{e}^{-\frac{\epsilon }{T}}\cdot e^{-\lambda \epsilon }d\epsilon ={\displaystyle \frac{1}{T}}{\int }_0^{\mathrm{\infty }}{e}^{-\left(\lambda +\frac{1}{T}\right)\epsilon }d\epsilon =-{\displaystyle \frac{1}{T\left(\lambda +\frac{1}{T}\right)}}{e}^{-\left(\lambda +\frac{1}{T}\right)\epsilon }{|}_0^{\mathrm{\infty }}={\displaystyle \frac{1}{\lambda T+1}}=g(\lambda )}$.

Jamming transition

Granular systems are known to exhibit jamming and undergo a jamming transition which is thought of as a thermodynamic phase transition to a jammed state. The transition is from fluid-like phase to a solid-like phase and it is controlled by temperature, ${\displaystyle T}$, volume fraction, ${\displaystyle \varphi }$, and shear stress, ${\displaystyle \mathrm{\Sigma }}$. The normal phase diagram of glass transition is in the ${\displaystyle {\varphi }^{-1}-T}$ plane and it is divided into a jammed state region and unjammed liquid state by a transition line. The phase diagram for granular matter lies in the ${\displaystyle {\varphi }^{-1}-\mathrm{\Sigma }}$ plane, and the critical stress curve ${\displaystyle \mathrm{\Sigma }(\varphi )}$ divides the state phase to jammed\unjammed region, which corresponds to granular solids\liquids respectively. For isotropically jammed granular system, when ${\displaystyle \varphi }$ is reduced around a certain point, ${\displaystyle J}$, the bulk and shear moduli approach 0. The ${\displaystyle J}$ point corresponds to the critical volume fraction ${\displaystyle {\varphi }_c}$. Define the distance to point ${\displaystyle J}$, the critical volume fraction, ${\displaystyle \mathrm{\Delta }\varphi \equiv \varphi -{\varphi }_c}$. The behavior of granular systems near the ${\displaystyle J}$ point was empirically found to resemble second-order transition: the bulk modulus shows a power law scaling with ${\displaystyle \mathrm{\Delta }\varphi }$ and there are some divergent characteristics lengths when ${\displaystyle \mathrm{\Delta }\varphi }$ approaches zero. While ${\displaystyle {\varphi }_c}$ is constant for an infinite system, for a finite system boundary effects result in a distribution of ${\displaystyle {\varphi }_c}$ over some range.

The Lubachevsky-Stillinger algorithm of jamming allows one to produce simulated jammed granular configurations.

Pattern formation

Excited granular matter is a rich pattern-forming system. Some of the pattern-forming behaviours seen in granular materials are:

• The un-mixing or segregation of unlike grains under vibration and flow. An example of this is the so-called Brazil nut effect where Brazil nuts rise to the top of a packet of mixed nuts when shaken. The cause of this effect is that when shaken, granular (and some other) materials move in a circular pattern. Some larger materials (Brazil nuts) get stuck while going down the circle and therefore stay on the top.
• The formation of structured surface or bulk patterns in vibrated granular layers. These patterns include but are not limited to stripes, squares and hexagons. These patterns are thought to be formed by fundamental excitations of the surface known as oscillons. The formation of ordered volumetric structures in granular materials is known as Granular Crystallisation, and involves a transition from a random packing of particles to an ordered packing such as hexagonal close-packed or body-centred cubic. This is most commonly observed in granular materials with narrow size distributions and uniform grain morphology.
• The formation of sand ripples, dunes, and sandsheets

Some of the pattern-forming behaviours have been possible to reproduce in computer simulations. There are two main computational approaches to such simulations, time-stepped and event-driven, the former being the most efficient for a higher density of the material and the motions of a lower intensity, and the latter for a lower density of the material and the motions of a higher intensity.

Acoustic effects

Some beach sands, such as those of the aptly named Squeaky Beach, exhibit squeaking when walked upon. Some desert dunes are known to exhibit booming during avalanching or when their surface is otherwise disturbed. Granular materials discharged from silos produce loud acoustic emissions in a process known as silo honking.

Granulation

Granulation is the act or process in which primary powder particles are made to adhere to form larger, multiparticle entities called granules.

Crystallization

When water or other liquids are cooled sufficiently slowly, randomly positioned molecules rearrange and solid crystals emerge and grow. A similar crystallisation process may occur in randomly packed granular materials. Unlike removing energy by cooling, crystallization in granular material is achieved by external driving. Ordering, or crystallization of granular materials has been observed to occur in periodically sheared as well as vibrated granular matter. In contrast to molecular systems, the positions of the individual particles can be tracked in the experiment. Computer simulations for a system of spherical grains reveal that homogeneous crystallization emerges at a volume fraction ${\displaystyle \varphi =0.646\pm 0.001}$. The computer simulations identify the minimal ingredients necessary for granular crystallization. In particular, gravity and friction are not necessary.

Computational modeling of granular materials

Several methods are available for modeling of granular materials. Most of these methods consist of the statistical methods by which various statistical properties, derived from either point data or an image, are extracted and used to generate stochastic models of the granular medium. A recent and comprehensive review of such methods is available in Tahmasebi and other (2017). Another alternative for building a pack of granular particles that recently has been presented is based on the level-set algorithm by which the real shape of the particle can be captured and reproduced through the extracted statistics for particles' morphology.