Granular and particle-laden flows: from laboratory experiments to field observations

2.1.1. Particle/particle forces.

2.1.2. Fluid/particle forces.

In the presence of an interstitial fluid, particles are subjected to additional forces. First, a particle undergoes a fluid resistance force which is opposite to its relative motion. This drag force, F_d, can be expressed formally as:

$$\begin{eqnarray}&&{{F}_{\text{d}}}=C\left({{u}_{\text{f}}}-{{u}_{\text{p}}}\right)\end{eqnarray} \tag{ 2.1 }$$

where u_p and u_f are the velocities of the particle and fluid, respectively. The factor C depends on various parameters such as the particle Reynolds number, the local solid volume fraction ϕ, etc. For dilute flows of smooth spherical particles, at small particle Reynolds number, the fluid drag reduces to the Stokes force ($C=3\pi d{{\eta}_{\text{f}}}$, where η_f is the fluid dynamic viscosity and d the grain diameter), while at high Reynolds particle number, the drag is independent of the fluid viscosity and scales quadratically with the relative particle velocity ($C\approx 0.5(\pi /8){{d}^{2}}{{\rho}_{\text{f}}}\left({{u}_{\text{f}}}-{{u}_{\text{p}}}\right)$, where ρ_f is the fluid mass density). Note that (2.1) holds in principle only for steady particle motion. For unsteady motion, additional forces should be considered, such as the virtual mass effect and the Basset force (Johnson 1998).

Second, when two particles come close to each other, they may undergo a repulsive force known as lubrication force. This force arises from hydrodynamic pressure in the interstitial fluid being squeezed out from the space between two particles and is expressed as:

$$\begin{eqnarray}&&{{F}_{\text{L}}}\approx (3\pi /2){{\eta}_{\text{f}}}\delta {{u}_{\text{p}}}/h\end{eqnarray} \tag{ 2.2 }$$

where h is the gap between the particles and δu_p their relative approaching velocity. As a consequence, the interparticle contact forces may be drastically altered in the presence of an interstitial fluid. Gondret et al (2002) showed experimentally that the dissipation in a binary collision of particles immersed in a fluid is increased in comparison with the dry case and can be expressed in terms of an effective restitution coefficient e_eff which follows a master curve depending solely of the Stokes number St  =  ρ_pδu_p d/η_f where ρ_p is the particle mass density. At small Stokes number (i.e. St  <  10) the ratio e_eff/e_dry drops to zero, it approaches 1 at high Stokes number (St  >  10⁴) and, on a relatively large interval, e_eff can be approximated by (Yang and Hunt 2006):

$$\begin{eqnarray}&&{{e}_{\text{eff}}}\approx {{e}_{\text{dry}}}-62\left(1+{{e}_{\text{dry}}}\right)/\text{St}\text{.}\end{eqnarray} \tag{ 2.3 }$$

If the interstitial fluid is the air, lubrication force can be safely neglected for millimeter particle size, whereas in a liquid, the latter may play a significant role.

Third, particles may be sensitive to the turbulent fluctuations of the fluid flow. This is the case for fine particles in a turbulent fluid flow. The associated transport is called 'turbulent suspension'. The dimensionless number for determining whether turbulent suspension is effective or not is the Rouse number: $\text{Rouse}={{w}_{\text{p}}}/\kappa {{u}_{\ast}}$. It is a ratio between the particle fall velocity w_p and the vertical fluctuating grain velocity as a product of the von Karman constant κ and the shear velocity u_*. Rouse number smaller than unity indicates that turbulent suspension is effective.

Among granular flows and particle liquid mixtures, three major categories of flows can be identified according to their solid fraction (the typical values given here correspond to spherical monosized particles): very dilute (ϕ  <  0.01), dilute (0.1  <  ϕ  <  0.5) and dense flows (0.5  <  ϕ  <  0.59). For very dilute flows, particle/particle interactions play a minor role. These flows include dilute turbidity currents, powder snow avalanches and pyroclastic currents. The motion of particles is driven by the turbulent liquid flow and gravity forces. On the contrary, for dilute and dense flows, particle/particle interactions play an important role. These flows include aerial rock and debris avalanches, landslides, submarine avalanches, debris flows and dense pyroclastic currents. Dense flows can be divided into three subcategories depending on the nature of the particle/particle interactions: (i) particle inertia-dominated regime, (ii) viscous resistance-dominated regime and (iii) fluid inertial resistance-dominated regime. This classification is based on the analysis of the time scale of the particle displacement (see Courrech du Pont et al (2003a) and Cassar et al (2005) for further details).

For these different regimes, we can introduce a dimensionless ratio I defined as the ratio of the characteristic particle time scale to the characteristic time equal to the inverse of the flow shear rate $\overset{\centerdot}{{\gamma}}\,$ (for sake of simplicity we consider here a simple shear flow):

• Particle inertia-dominated regime (i): $I={{I}_{\text{pi}}}=\sqrt{{{\rho}_{\text{p}}}{{\overset{\centerdot}{{\gamma}}\,}^{2}}{{d}^{2}}/{{P}_{\text{p}}}}$
• Viscous resistance-dominated regime (ii): $I={{I}_{\text{v}}}=\frac{{{\eta}_{\text{f}}}\overset{\centerdot}{{\gamma}}\,}{{{P}_{\text{p}}}}$
• Fluid inertia-dominated regime (iii): $I={{I}_{\text{fi}}}=\sqrt{{{\rho}_{\text{f}}}{{\overset{\centerdot}{{\gamma}}\,}^{2}}{{d}^{2}}/{{P}_{\text{p}}}}$.

P_p is the confining pressure, ρ_p and ρ_f are the particle and fluid density respectively. Note that we drop the numerical factor in the definition of I to make it consistent with the usual definition of inertial number (see equation (2.10)). The inertial regime (or Bagnoldian regime) has raised a great interest in the granular physics community and various theoretical approaches have been developed to model this flow regime. Two emblematic approaches are the granular kinetic theory for collisional granular flows and the so called 'μ(I)' rheology for dense granular flows.

Alternative dimensionless numbers have been introduced in the literature (Iverson 1997) on the basis of the comparison of particle and fluid stresses. One can mention the Savage and Bagnold number defined respectively as the ratio of the inertial grain shear stress to the weight of flowing layer per unit surface

$$\begin{eqnarray}&&{{N}_{\text{Sav}}}=\frac{{{\rho}_{\text{p}}}{{\overset{\centerdot}{{\gamma}}\,}^{2}}{{d}^{2}}}{\left({{\rho}_{\text{p}}}-{{\rho}_{\text{f}}}\right)gH}\end{eqnarray} \tag{ 2.4 }$$

and the ratio of inertial grain shear stress to viscous shear stress,

$$\begin{eqnarray}&&{{N}_{\text{Bag}}}=\frac{\phi {{\rho}_{\text{p}}}{{\overset{\centerdot}{{\gamma}}\,}^{2}}{{d}^{2}}}{(1-\phi ){{\eta}_{\text{f}}}\overset{\centerdot}{{\gamma}}\,}\end{eqnarray} \tag{ 2.5 }$$

where H is the thickness of granular layer. The Savage number N_Sav is nothing but the square of the inertial number I_pi in which the confining pressure is driven by the hydrostatic force (P_p ~ (ρ_p  −  ρ_f)gH). The Bagnold number is a simple combination of I_pi and I_v: N_Bag ~ $I_{\text{pi}}^{2}$/I_v. High Bagnold number thus refers to particle inertia-dominated regime.

We focus here on the inertial regime, where the stress due to the grains is far larger than that exerted by the fluid, so that one may ignore the effects of the ambient fluid and consider inter-particle interactions only.

2.3.1. Collisional granular flows: kinetic theory.

In rapid granular flows, the grains interact mainly through binary collisions in a way that is reminiscent of molecules interacting in a gas. There are nevertheless significant differences between a 'granular gas' and a molecular gas. The disparity of sizes between grains and molecules is not the direct cause of these differences. The fact that grain collisions are not elastic, and thus dissipative, has major implications concerning the behaviour of grain assemblies: granular gases are always in non-equilibrium states.

It is common to define the granular temperature in a slightly different way than for molecular gases: T is defined as the average (over volumes and time ranges) of the square of the grain velocity fluctuations. This statistical definition does not need to refer to a state of equilibrium. One of the main applications of the concept of granular temperature is the construction of kinetic theories. The term 'kinetic theory' is often used in different contexts. The more general one is a mesoscopic approach where the state of the system is described using the distribution function specifying the distribution of particle positions and velocities. In most cases, it cannot be solved exactly. The other one, called Chapman–Enskog expansion, is an approximation where the two particle distribution function is expressed as the product of the single particle velocity distribution functions, and the pair distribution at contact. This approach neglects the correlations. The hydrodynamic equations of motion resulting from the application of these theories resemble the compressible Navier–Stokes equations:

$$\begin{eqnarray}&&\partial \left({{\rho}_{\text{p}}}\phi \right)/\partial t+\nabla \centerdot \left({{\rho}_{\text{p}}}\phi \mathbf{u}\right)=0\end{eqnarray} \tag{ 2.6 }$$

$$\begin{eqnarray}&&{{\rho}_{\text{p}}}\phi \left(\partial \mathbf{u}/\partial t+\mathbf{u}\centerdot \nabla \mathbf{u}\right)=\nabla \centerdot \boldsymbol{\sigma }+\,\mathbf{f}.\end{eqnarray} \tag{ 2.7 }$$

Here, ϕ is the volume fraction, u the mean velocity, σ the stress tensor and f the body force density. The equation for the energy density, e, contains a sink term Γ_d that represents the loss of energy due to the inelasticity of the collisions (q is the energy flux vector):

$$\begin{eqnarray}&&{{\rho}_{\text{p}}}\phi \left(\partial e/\partial t+\mathbf{u}\centerdot \nabla e\right)=-\nabla \centerdot \mathbf{q}+\boldsymbol{\sigma }:\nabla \mathbf{u}-{{ \Gamma }_{d}}.\end{eqnarray} \tag{ 2.8 }$$

This term is responsible for the existence of many phenomena that characterize granular gases. The hydrodynamic description is completed by constitutive relations expressing σ, q and Γ_d in terms of the hydrodynamic fields. Kinetic theory provides local constitutive relations through the application of the Chapman–Enskog expansion. The constitutive relation for τ is expressed as in a Newtonian fluid:

$$\begin{eqnarray}&&\boldsymbol{\sigma }=\left[-p+(\lambda -2\eta /3)\nabla \centerdot \mathbf{u}\right]1+\eta \left[\left(\nabla \mathbf{u}+{{\left(\nabla \mathbf{u}\right)}^{T}}\right)\right].\end{eqnarray} \tag{ 2.9 }$$

The energy flux vector $\mathbf{q}=-\kappa \nabla \centerdot T$, the pressure p  =  ρ_pϕ(1  +  4ϕg₁₂)T, the viscosities λ and η, the conductivity κ and the sink term Γ_d are given in terms of the volume fraction, the granular temperature T  =  2e/3 and the pair distribution at contact g₁₂(ϕ) (Jenkins and Savage (1983), see also Delannay et al (2007) and references herein).

The solutions of these equations have had some success in describing relatively dilute flows or moderately dense flows in the absence of gravity (Azanza et al 1999, Forterre and Pouliquen 2001, Xu et al 2003).

In the dense limit, this approximation generally fails. An improvement of the models to incorporate the effect of correlations is necessary to describe the dense flow regime. The challenges for applying such ab initio theories for dense granular flows are discussed in Delannay et al (2007) and Kumaran (2015).

2.3.2. Dense granular flows.

In the dense flow regime, collisions still exist, but they can involve more than two bodies and become so frequent that very strong energy dissipation occurs (inelastic collapse). Contacts between particles tend to be enduring and long-lasting; the contact network can percolate through particles. The enduring frictional contacts induce correlations between the velocities and/or locations of the particles.

There have been several attempts to extend the simple kinetic theory to include higher volume fractions, greater dissipation in collisions, and multiple or enduring contacts. The existing constitutive relations of the kinetic theory need to be modified by the incorporation of the correlations associated with the enduring and/or repeated contacts. However, the appropriate form of this modification has yet to be determined. Up to now no constitutive equations are available in the dense regime.

Several studies (GDR-MiDi 2004, Da Cruz et al 2005) have pointed out that simple dimensional arguments could provide an interesting framework for constitutive laws for dense granular flows. These may be phrased as empirical constitutive relations for the shear stress and normal stress in steady homogeneous shearing. An assembly of rigid frictional spheres of diameter d and mass density ρ_p, confined under a normal stress P in between two bumpy planes, is sheared at a given shear rate $\overset{\centerdot}{{\gamma}}\,$ by applying a shear stress τ. If one neglects the finite size of the sample, the non-dimensional inertial number I defined above (section 2.2) is the only dimensionless parameter of the problem. As a consequence, the shear stress has to be proportional to P times a function of I. The ratio of tangential to normal stresses is the effective coefficient of friction μ:

$$\begin{eqnarray}&&\tau /P=\mu (I)\,\,\,\,\,\text{with} ~ \,\,\,\,I=\overset{\centerdot}{{\gamma}}\,d/\sqrt{P/{{\rho}_{\text{p}}}}.\end{eqnarray} \tag{ 2.10 }$$

The form of the function μ(I) can be obtained from the numerical simulations of plane shear (Da Cruz et al 2005) or indirectly by experimental measurements for flow down inclined planes (see section 3.2). The function μ(I) goes from a minimum value μ_s for very low I (quasi-static regime) up to an asymptotical finite value μ₂ when I increases. It can be expressed by the following expression (Jop et al 2005):

$$\begin{eqnarray}&&\mu (I)={{\mu}_{\text{s}}}+\left({{\mu}_{2}}-{{\mu}_{\text{s}}}\right)/\left(1+{{I}_{0}}/I\right)\end{eqnarray} \tag{ 2.11 }$$

where I₀ is a constant. The values of coefficients are material-dependent, for instance typical values for glass bead are μ_s  =  tan(21°), μ₂  =  tan(33°) and I₀  =  0.28 (GDR-Midi 2004, Jop et al 2005).

Another relation, also derived from dimensional analysis, gives the variation of the volume fraction; ϕ is, in general, an affine decreasing function of I (Da Cruz et al 2005). This relation completes the constitutive laws:

$$\begin{eqnarray}&&\phi ={{\phi}_{\max}}-\left({{\phi}_{\max}}-{{\phi}_{\min}}\right)I.\end{eqnarray} \tag{ 2.12 }$$

Typical values are ϕ_max  =  0.6 and ϕ_min  =  0.5 (Pouliquen et al 2005). The latter equation is only for I lower that 1, i.e. for volume fraction larger than ϕ_min.

A generalization of this approach to tensorial constitutive relations for incompressible 3D flows has been proposed (Jop et al 2006):

$$\begin{eqnarray}&&\boldsymbol{\sigma }=-P\mathbf{I}+\boldsymbol{\tau },\,\text{with}\,\boldsymbol{\tau }=\eta \overset{\centerdot}{{\gamma}}\,,\,\text{and}\,\eta =\frac{\mu (I)P}{\|\overset{\centerdot}{{\gamma}}\,\|},I=\frac{\|\overset{\centerdot}{{\gamma}}\,\|d}{\sqrt{P/{{\rho}_{\text{p}}}}}\end{eqnarray} \tag{ 2.13 }$$

P represents an isotropic pressure, by definition it is the negative one-third the trace of the stress tensor, and τ is the deviatoric component. The effective viscosity η is a function of the second invariant $\|\overset{\centerdot}{{\gamma}}\,\|$ of the strain rate tensor $\overset{\centerdot}{{\gamma}}\,$. Its definition is related to the friction coefficient law μ(I) (2.11).

This relation can be extended for the compressible case by replacing the strain rate tensor by its deviatoric component in (2.13) (Börzsönyi et al 2009, Cortet et al 2009). It implies, in addition, a univocal dependence of the volume fraction on the inertial number.

Considering (2.13), in the limit of a shear rate going to zero, one can see that the material flows only if the following Drucker–Prager-like yield criterion is satisfied: $\|\tau \|$  >  μ_s P. Below the threshold, the medium behaves locally as a rigid body. Note that the constitutive relation (2.13) can be separated into two terms: a Drucker–Prager yield stress term involving the constant friction coefficient μ_s, and a viscous term with a viscosity depending on both the pressure and the norm of the strain rate tensor (Ionescu et al 2015). The specificity compared to classical Bingham or Herschel–Bulkley fluids is that the yield stress depends on the local pressure and the effective viscosity depends both on the shear rate and on the local pressure.

The μ(I) rheology had some success in reproducing quantitatively experimental observations in relatively simple situations such as granular column collapses (see section 3.4.2). Nevertheless, (2.13) contains very strong assertions: stress and strain rate tensors have to be aligned, and the second invariants should satisfy the relation:

$$\begin{eqnarray}&&\frac{\|\boldsymbol{\tau }\|}{P}=\mu (I)\end{eqnarray} \tag{ 2.14 }$$

where the function μ(I) follows (2.11). The tensorial relation (2.13) has been tested by numerical simulations using DEM (see section 2.7) in different situations (Börzsönyi et al 2009, Cortet et al 2009). It turns out that the alignment between stress and strain rate tensors fails. Interestingly, the relation (2.14) works better. But it does not seem to give a friction coefficient law similar to (2.11). A logarithmic decay is observed (Cortet et al 2009) for vanishing values of I (leading to very low values of the friction coefficient) and a non-monotonic variation of μ(I) is reported at larger values of I (with high value for the friction coefficient) (Börzsönyi et al 2009). A granular rheology that generalises the μ(I) model and incorporates first and second normal stress differences has been recently proposed (McElwaine et al 2012).

Another questionable point is that such a rheology does not use the notion of granular temperature which is at the base of the kinetic theory (Jenkins and Richman 1985) even in the case of dense flows. The friction function is purely phenomenological without explicit links between the coefficients used and the grain properties. Moreover, the variation of the effective friction μ with the inertial number I—which has been obtained empirically and not theoretically—has been recently questioned (Holyoake and McElwaine 2012). Clearly, a link between the μ(I)-rheology and the grain-scale physics is still missing. Note also that it has been recently shown (Barker et al 2015) that, in the incompressible limit, the μ(I) rheology is only mathematically well-posed for intermediate values of the inertial number. At high and low inertial numbers the equations are ill-posed, i.e. the growth rate of small perturbations increases infinitely fast in the high wavenumber limit. The practical significance of this is that 2D time-dependent simulations in the ill-posed regime will exhibit short wavelength instabilities that will get progressively stronger as the resolution is enhanced, i.e. the results will be grid-dependent. This points out that there is some important missing physics (force chains, effect of grain stiffness) for both high and low inertial numbers where the problem is ill-posed.

The μ(I)-rheology is a local rheology in its current form, it is thus not valid close to the jamming transition (Deboeuf et al 2005, Staron et al 2010), or when non-local effects are not negligible (Nichol et al 2010, Reddy et al 2011). As already mentioned, the enduring frictional contacts induce correlations between the velocities and locations of the particles. These spatial correlations become very important, both in the force network and velocity fluctuations, near the jamming transition. To incorporate these non-local effects in the constitutive relations is not an easy task. It is a generic question emerging in the study of other complex fluids near jamming: foams, emulsions, colloidal suspensions. Some of the attempts to develop non local theories are described in GDR-MiDi (2004), Delannay et al (2007), Henann and Kamrin (2013) and Jop (2015). For example, Jenkins (2006, 2007) modified the kinetic theory by introducing a length associated with the size of particle clusters into the expression of the rate of collisional dissipation. This theory has been extended to very dissipative, frictional spheres (Jenkins and Berzi 2010). The latter theory has the capacity to reproduce, at least qualitatively, experimental results of dense, unconfined, inclined granular flows over a rigid bumpy base and dense, confined flows over an erodible base (see section 3.2).

2.4.1. Introduction.

2.4.2. Single-phase model.

2.4.3. Two-phase model.

2.4.4. Steady and fully developed regimes.

We consider a unidirectional, steady and fully developed flow down an incline (see figure 3), consisting of a mixture of granular material and interstitial fluid, each of constant specific mass density ρ_p and ρ_f, respectively.

Zoom In Zoom Out Reset image size

In this simple flow configuration, the momentum equations take the following forms,

$$\begin{eqnarray}\begin{array}{*{35}{l}} \frac{\text{d}{{\tau}_{xz,\text{p}}}}{\text{d}z} & +\,C\phi \left({{u}_{\text{f}}}-{{u}_{\text{p}}}\right)\,+{{\rho}_{\text{p}}}\phi g\sin \theta =0 \\ {} & -\,\frac{\text{d}{{P}_{\text{p}}}}{\text{d}z}-\phi \frac{\text{d}{{P}_{\text{f}}}}{\text{d}z}+{{\rho}_{\text{p}}}\phi g\cos \theta =0 \end{array}\end{eqnarray} \tag{ 2.15 }$$

for the particle phase, and

$$\begin{eqnarray}\begin{array}{*{35}{l}} \frac{\text{d}{{\tau}_{xz,\text{f}}}}{\text{d}z} & -\,C\phi \left({{u}_{\text{f}}}-{{u}_{\text{p}}}\right)\,+{{\rho}_{\text{f}}}(1-\phi )g\sin \theta =0 \\ {} & -\,(1-\phi )\frac{\text{d}{{P}_{\text{f}}}}{\text{d}z}+{{\rho}_{\text{f}}}(1-\phi )g\cos \theta =0 \end{array}\end{eqnarray} \tag{ 2.16 }$$

for the fluid phase. u_p and u_f stand for the velocities of the solid and fluid constituents, respectively, τ_xz,p and τ_xz,f for the particle and fluid shear stress, P_p and P_f for the particle and fluid pressure, and ϕ for the solid volume fraction. The factor C which quantifies drag, depends on various parameters (Reynolds number, solid fraction ...). The simplest form of the latter factor that incorporates viscous and form drag and concentration dependence reads as (Dallavalle 1943, Richardson and Zacki 1954):

$$\begin{eqnarray}&&C=\frac{{{\rho}_{\text{f}}}}{d{{(1-\phi )}^{3.1}}}\left(\frac{3}{10}|{{u}_{\text{f}}}-{{u}_{\text{p}}}|\,+18.3\frac{{{\eta}_{\text{f}}}}{{{\rho}_{\text{f}}}d}\right).\end{eqnarray} \tag{ 2.17 }$$

To complete the system equations, closure models for both the fluid and particle shear stress should be supplemented. For turbulent liquid flows, the fluid shear stress can be described by a simple classical turbulent closure based on the mixing length theory (Berzi and Jenkins 2008) as:

$$\begin{eqnarray}&&{{\tau}_{xz,\text{f}}}=-(1-\phi ){{\rho}_{\text{f}}}\left({{\nu}_{\text{f}}}+{{\nu}_{\text{t}}}\right)\frac{\partial {{u}_{\text{f}}}}{\partial z}\end{eqnarray} \tag{ 2.18 }$$

where ν_t is the turbulent viscosity modeled by ${{\nu}_{\text{t}}}={{l}^{2}}{{\partial}_{z}}{{u}_{\text{f}}}$. The spatial length l is the mixing length and is proportional to the distance from the flow base (l  =  κz, where κ is the von Karman constant) for a turbulent liquid flow free of particles. The presence of particles is however expected to influence the turbulence and thus to modify the mixing length. In absence of any experimental data, the simplest approximation is to use the law (2.18). However, more elaborate models have been developed to account for the influence of particles on fluid turbulence (Jenkins and Hanes 1998, Hsu et al 2003, 2004, Cantero et al 2012).

The closure for the particle phase is another crucial issue which is not yet solved. Depending of the flow regime, different closure models have been proposed.

• (i)  
  For collisional flows of heavy grains (called sheet flows), where collisional dissipation dominate viscous dissipation (Jenkins and Hanes 1998), the particle shear stress can be derived from the kinetic theory of granular gas (see section 2.3.1). The particle pressure and the particle viscosity are thus given in terms of the granular temperature. This simple model has the capacity to predict profiles of particle concentration and particle mean velocity in steady and fully developed sheet flows.
• (ii)  
  For dense liquid particle mixture, several recent approaches employing the 'μ(I)' rheology for the particle shear stress closure (see section 2.3.2) have been developed by Cassar et al (2005). For dense viscous granular suspension, a more refined relation has been proposed by Boyer et al (2011). According to the flow regime (inertial or viscous), the pertinent inertial parameter to be considered is respectively I_pi and I_v. In this context, Berzi and Jenkins (2008) have considered the problem of steady and fully developed flows of particles and fluid down a slope. They made the assumption that the fluid rheology is described by an eddy viscosity model (see (2.18)) and the granular phase by a shear-rate-dependent friction law μ(I). They investigated the case where the particle inertia dominates the fluid drag (i.e. particle inertia-dominated regime). In other words, they assumed that I  =  I_pi. This approach was shown to be successful for predicting quantitatively the velocity and density profiles observed in the experiments of Armanini et al (2005) in which a mixture of fluid and particles is recirculated through an inclined channel. This approach is restricted to steady and fully developed flows and therefore disregard dilatancy effects.

For describing bed load transport in laminar flows, Aussilous et al (2013) have implemented a two-phase continuum model using a shear-rate-dependent friction law μ(I_v) for the granular phase and a volume-fraction-dependent viscosity for the fluid phase. This approach provides reasonable predictions for particle velocity and particle density profiles.

2.4.5. Unsteady flows.

We shall now discuss the capability of the two-phase approaches to model the transient dynamics during which the mixture can contract or dilate.

During the transient, part of the relative vertical displacements between the fluid and the granular phase are caused by the dilation and compaction of the granular material.

These displacements affect the pressure field. The classical approach is to split the fluid pressure into two contributions:${{P}_{\text{f}}}=P_{\text{f}}^{\text{h}}+P_{\text{f}}^{\text{e}}$, where $P_{\text{f}}^{\text{h}}$ is the hydrostatic fluid pressure, verifying $\nabla P_{\text{f}}^{\text{h}}={{\rho}_{\text{f}}}\mathbf{g}$ and $P_{\text{f}}^{\text{e}}$ is the excess pore-fluid pressure.

The key issue is thus to provide a closure relation between the particle volume fraction and the excess pore pressure. Iverson (1997) argued that it is reasonable to assume that the particle volume fraction is a function of only the excess pore pressure and that the change in particle volume fraction due to changes in excess pore pressure can be expressed in terms of the compressibility coefficient $K=\partial \phi /\partial P_{\text{f}}^{\text{e}}$ of the sediment-water mixture. Within this hypothesis, the equation for the particle volume fraction reduces to a standard diffusion equation that can readily be solved.

Other closure equations can be found in the literature. For example, Pailha and Pouliquen (2009) consider an equation relating the relative variation of the particle volume fraction to the shear-induced dilatancy angle. More precisely, their model is inspired from the critical state theory proposed by Roux and Radjai (1998) and is based on the concept of the dilatancy angle ψ, which gives the rate of dilation (compaction) of the material under shear. Based on this theory, the closure equations for the particle shear stress and pressure are modified and rewritten as:

$$\begin{eqnarray}&&{{\tau}_{xy,\text{p}}}-\tau _{xy,\text{p}}^{\text{eq}}=\tan \psi {{P}_{\text{p}}}\end{eqnarray} \tag{ 2.19 }$$

$$\begin{eqnarray}&&\phi -{{\phi}^{\text{eq}}}=\frac{\tan \psi}{K}\end{eqnarray} \tag{ 2.20 }$$

supplemented by

$$\begin{eqnarray}&&\frac{1}{\phi}\frac{\text{d}\phi}{\text{d}t}=\tan \psi \overset{\centerdot}{{\gamma}}\,\end{eqnarray} \tag{ 2.21 }$$

where $\overset{\centerdot}{{\gamma}}\,$ is the local particle shear rate. $\tau _{xy,\text{p}}^{\text{eq}}$ and ϕ^eq are the stress and the volume fraction in the steady regime (Pailha and Pouliquen 2009). This model represents the simplest shear-rate dependent critical state theory for granular/fluid mixtures and has been used successfully to model submarine avalanches (see next section).

The physical principles on which the equations of continuum models are based are the balance laws of mass and momentum (and energy) and constitutive closure equations describing the mechanical behaviour of the material involved. These equations are then used together with boundary conditions to predict the motion of the system.

The behaviour of granular flows contrasts with the one of single phase newtonian fluids as the former depend crucially on the nature of the boundaries (i.e. the properties of the base atop which the flow occurs and the possible presence of lateral confinement). Most experimental and numerical studies generally consider bumpy, flat frictional or erodible basal boundaries. Bumpy boundaries usually consist of particles glued on a flat substrate. In most studies, glued particles are the same as those of the flowing material. However, the bumpiness can be increased or decreased by using grains larger or smaller than the flowing grains. The spatial repartition of the grains can also be used to modify the bumpiness of the boundary.

The conditions at the boundaries, for granular flows, differ from the one used for Newtonian fluids. For the latter, no-slip conditions at solid boundaries is often achieved. In contrast, a granular flow slips relative to a boundary; so, in addition to collisional dissipation at the boundaries, fluctuation energy can be generated by slip (Richman 1988). Consequently, boundaries can either provide or remove fluctuation energy from the flow.

Establishing ab initio boundary conditions for granular gases is complex: in addition to the difficulty of solving the Boltzmann equation near a solid boundary, one has to know the properties of the collision between the particle and the boundary. The feed-back effect of the boundary on the distribution function of the incoming particles must be taken into account. Approximate boundary conditions at solid surfaces describing the way energy fluctuation is provided to the flow can be applied to flows down bumpy inclined (Richman and Marciniec 1990, Jenkins 1994, Delannay et al 2007, Kumaran 2015).

Within the framework of the μ(I) rheology, boundary conditions also cause difficulties. For example, the effect of changing the degree of bumpiness of a bumpy base is not clear: does it only modify the shear and velocity in a boundary layer? Or does it change the whole bulk rheology (Goujon et al 2007). It has been shown that relatively small changes at the boundary can induce transitions from disordered to ordered flow states (Silbert et al 2002, Kumaran and Maheshwari 2012, Kumaran and Bharathraj 2013). The slip velocity at the boundary (flat or, to a lesser extent, bumpy) should be related to the physical characteristics of the grains and the boundary (Artoni et al 2009, 2012).

Free surface flows induce other questions concerning boundary conditions. Phase transition from liquid to gas may occur at the free surface. As the pressure decreases, the number I becomes very large and the kinetic regime appears. This transition is neither predicted nor described by the μ(I) rheology. Another transition also appears in free conditions: as will be seen in section 3.3, unconfined flows can create lateral levées, which correspond to erodible boundaries. It can be seen as a jamming/unjamming transition. The same behaviour is also observed on very dissipative bases (Louge et al 2015), and when there is the spontaneous formation of a heap, called SSH (see section 3.2). The flowing 'liquid' zone is in contact with a quasi-static part where the displacement of the grains decreases exponentially. The number I vanishes when approaching this part. The μ(I) rheology does not describe this quasi-static regime, which probably requires a non-local rheology, as already mentioned (section 2.3.2).

In the case of a flow over an erodible base composed of the same material, different approaches have been developed to describe the erosion/deposition process. One of them is the BCRE model (Bouchaud et al 1994, Boutreux et al 1998). Assuming a thin flow on a pile slope close to the angle of repose, this model provides a continuum phenomenological description of the surface dynamics by two variables: b(x,t), the local height of the pile (profile of immobile grains) and h(x,t), the thickness of flowing grains. These variables are dynamically coupled by an interaction term Γ(b,h) allowing conversion from flowing to static grains and vice versa (deposition/erosion). The mass conservation is then expressed as:

$$\begin{eqnarray}&&\begin{array}{*{35}{l}} \frac{\partial b}{\partial t}=- \Gamma \\ \frac{\partial h}{\partial t}+U\frac{\partial h}{\partial x}=D\frac{{{\partial}^{2}}h}{\partial {{x}^{2}}}+ \Gamma \end{array}\end{eqnarray} \tag{ 2.22 }$$

where U is the constant drift downhill velocity of the rolling grains, and D is a constant diffusion coefficient. To close the model, a phenomenological relation between the interaction term Γ and the local slope θ has to be written. A first approach is to set the amplitude of exchange rate with the distance to a neutral inclination angle θ_n: $\Gamma ={{\gamma}_{\text{c}}}R\left(\theta -{{\theta}_{\text{n}}}\right)$ so that if θ  >  θ_n erosion is dominant (else it is deposition), γ_c is the collision frequency of a mobile grain with the static phase. When h is thicker than a few grain diameters (h  >  λ) the previous exchange term is not realistic, grains from the upper layers do not interact with the static phase (Boutreux et al 1998): $\Gamma ={{v}_{\text{up}}}\left(\theta -{{\theta}_{\text{n}}}\right)$, where v_up is a constant corresponding to γ_c λ. In their recent review on granular surface flows, Iverson and Ouyang (2015) have listed the alternative formulations of the exchange rate term that have been proposed in the literature. We shall point out that this simple phenomenological model can be formally derived from hydrodynamic models based on Saint-Venant equations (e.g. Gray 2001, Aradian et al 2002, Bouchut et al 2008). However, as highlighted in Bouchut et al (2008), the energy balance of these models is not always consistent. More refined descriptions of erosion and deposition processes can be developed from the knowledge of the constitutive law of the granular fluid coupled with a yield criterion without additional phenomenological law for the erosion rate (Lusso et al 2016a, Bouchut et al 2016b). The exchange rate cannot be derived explicitly except in very simple configurations.

For fluid particle mixture continuum models, the same difficulties arise for accounting properly for what happens at the boundaries and for deriving appropriate boundary conditions. Recently, an interesting approach has been introduced by Ancey and Heyman (2014) as a microstructural model of bedload transport, which describes the advection and scattering of coarse particles carried by a turbulent water stream down a sloping granular bed. Fluctuations of the particle flux are generated by particle exchanges with the bed consisting of particle entrainment and deposition. The evolution of the number of moving particles is described probabilistically using a coupled set of reaction-diffusion master equations.

Another key issue in fluid particle mixture continuum models is to describe what happens at the free surface. In particular it is important to describe different free surfaces for the solid and fluid phases to make it possible for the fluid to be sucked or expelled from the granular material when dilatancy effects are accounted for (Bouchut et al 2015).

Practically, solving the mass and momentum conservation full set of equations requires a prohibitive computational cost when applied to natural flows over real topography and simplifications are necessary.

A first approximation is to assume that the flow is thin compared to its extent along the slope (i.e. small aspect ratio). This is the case for most geophysical flows that are a few meters thick and travel distances of several hundred meters to several kilometers. Asymptotic development of the equations based on this thin-layer approximation makes it possible to neglect some terms such as the vertical acceleration, generally leading to the hydrostatic pressure assumption. Thin layer models solving the resulting simplified 2D or 3D equations have been proposed based on the partial fluidization method that takes into account an order parameter describing the state of the granular matter from static to flowing (Aranson and Tsimring 2002, Mangeney et al 2007b, Aranson et al 2008) or by using visco-plastic models with a yield stress separating static and flowing zones (Lusso et al 2016a, Bouchut et al 2016b).

Assuming thin-layer flows together with averaging the equations over the depth of the flow makes it possible to further reduce the computational cost as proposed in the pioneering work of Savage and Hutter (1989) partially based on former german (Voellmy 1955) and russian (Grigorian et al 1967) works on snow avalanches (figure 4). In that case, the whole granular material is assumed to be flowing. The resulting system (so-called Saint-Venant or shallow-water) involves two unknowns, the thickness and the mean (depth-averaged) velocity of the flowing mass. This hydrodynamic model is widely used to simulate natural flows. The first order asymptotic development leads essentially to a hydrostatic pressure (including centrifugal forces) even though attempts have been made to include the vertical acceleration (Denlinger and Iverson 2004). The mass and momentum conservation equations in the simplest 1D case (i.e. flow over a 2D topography b(x)) then read

$$\begin{eqnarray}&&\begin{array}{*{35}{l}} \frac{\partial h}{\partial t}+\frac{\partial}{\partial X}\left(h\bar{u}\right)\,=0 \\ \frac{\partial}{\partial t}\left(h\bar{u}\right)\,+\frac{\partial}{\partial X}\left(h{{{\bar{u}}}^{2}}\right)\,=gh\sin \theta +\,\frac{\partial}{\partial X}\left(g\cos \theta \frac{{{h}^{2}}}{2}\right) \\ \,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,-\,\mu h\left(g\cos \theta +\frac{{{{\bar{u}}}^{2}}}{R}\right)\frac{{\bar{u}}}{|\bar{u}|} \end{array}\end{eqnarray} \tag{ 2.23 }$$

where $\bar{u}(X,t)$ is the depth-averaged velocity in the downslope direction X (see figure 4) and $|\bar{u}|$ it's absolute value, $h(X,t)$ is the flow thickness, $\theta (X,t)$ the slope angle, and $R(X,t)$ the curvature radius of the topography. In that case, the momentum equation reduces to the balance between the acceleration, the gravity force, the force related to the pressure gradient and the friction force including centrifugal acceleration effects.

Zoom In Zoom Out Reset image size

Formally, the friction terms derived from the Coulomb friction law based on this depth-averaged asymptotic method appear both in the basal shear stress and in the pressure gradient term when anisotropy of normal stresses is taken into account (e.g. Gray et al 1999, Denlinger and Iverson 2004). However, because accurate development of the thin-layer approximation for a Coulomb friction law has never been done over a 3D arbitrary topography and because anisotropy effects have been questioned for granular flows (e.g. Ertas et al 2001), the models generally only consider the term related to the basal shear stress to describe the friction dissipation (e.g. Mangeney-Castelnau et al 2005). As a result, these models loose the relative importance of the sliding and shearing contributions to the flow. Consequently, one of the key issues of the method is the modelling of the averaged friction dissipation at the base. The basal friction is generally described phenomenologically in the framework of Coulomb friction using friction coefficient which may be constant or depend on the inertial number (see for example Pouliquen and Forterre 2002, Mangeney-Castelnau et al 2003, Mangeney et al 2007a, Pirulli and Sorbino 2008). Note, however, that fitting the basal friction dissipation (the only dissipation term in the model) to reproduce deposit extent of laboratory or field scale events result implicitly in an empirical quantification of the whole effective friction within the flow (Lucas et al 2014).

One of the main difficulties in deriving thin-layer depth-averaged models is to appropriately take into account topography effects, which is a key issue for application to natural landslides. Indeed, these flows are thin in the direction perpendicular to the flow and not in the vertical direction (figure 4). For these reasons, a series of asymptotic developments have been performed in a variable reference frame linked to the topography for complicated but specific topographies (see Pudasaini et al (2007) for a review). By considering a fixed reference frame but applying the thin-layer approximation in the direction perpendicular to the slope, Bouchut et al (2003) and Bouchut and Westdickenberg (2004) derived a model that takes into account accurately all the curvature terms for flows over an arbitrary topography, including centrifugal forces ignored in other models (see appendix A of Mangeney et al (2007a)). These terms can be important in the dynamics of natural landslides (Favreau et al 2010, Moretti et al 2015).

By assuming that the whole material flows, the Savage–Hutter model fails to capture the features of flows for which interaction with an erodible base prevails. This is the case for example for dam-break granular flows (see section 3.4). Indeed the transition between static and flowing regions in a depth-averaged model requires the determination of an equation for the evolution of the moving interface between the flowing and non-flowing parts, an interface that is not always well defined. As discussed in section 2.5, this point is currently a matter of debate, notably because a quantitative study of the erosion/deposition rate at the grain scale is still missing. Moreover, the proposed approaches are generally overly schematic to be extended to natural flows over real topography (Bouchut et al 2008). A key new approach using the depth-averaged μ(I)-rheology and prescribed Bagnold velocity profiles (Gray and Edwards 2014) explicitly model the formation of erosion-deposition waves (Edwards and Gray 2015), where completely stopped regions can form between mobile wave crests. A depth-integrated model taking into account a linearization of the μ(I)-rheology and prescribed S-shaped velocity profile was also proposed by Capart et al (2015). This model reproduces granular flows over erodible beds. However, by prescribing a specific shape of the velocity profiles, these models cannot reproduce the different profiles observed in highly transient flows. A new alternative approach based on mutlilayer shallow models has been recently proposed that makes it possible to calculate the different velocity profiles within the shallow granular layer and thus the static/flowing transition (Fernandez-Nieto et al 2016).

The depth-averaged approach has been also used for modeling particles/liquid mixtures. In particular, depth-averaged mixture models and two-phase models have been developed by Iverson and Denlinger (2001), Denlinger and Iverson (2001) and Pitman and Lee (2005), Pelanti et al (2008) respectively. In these versions of the two-phase flow equations, the flowing layer is assumed to be thin. This model allows describing, within a tractable approach, unsteady and non-uniform configurations encountered in real geophysical flows. In these works, the granular phase is described as a Coulomb material but dilatancy is disregarded. Bouchut et al (2015) showed that the constitutive equation required to close the system of two-phase flows has been replaced by two free surface boundary conditions in Pitman and Lee (2005), and Pelanti et al (2008), leading to a non-dissipative energy equation for the resulting model. This problem can be solved by prescribing the incompressibility of the solid phase.

The coupling between the dilatancy of the granular layer and the pore pressure was taken into account in depth-averaged mixture models by Iverson (1997, 2005). In this approach, Iverson explicitly incorporates an equation of the pore fluid pressure as discussed earlier. The landslide is described as a rigid block sliding down a slope, the source of dilatancy being localized at the base. The pore pressure builds up at the base and diffuses through the granular phase which is assumed to behave as a poro-elastic medium. The limit of this approach is that the long term evolution and the fully developed regime cannot be described because no shear rate dependence is incorporated in the rheology of the granular phase. Later on Iverson and George (2014) incorporated the shear rate dependence.

Pailha and Pouliquen (2009) combine the depth-averaged two-phase approach with the rheogical model for fluid/granular mixtures complemented by the critical state theory by Roux and Radjai (1998). This model is able to describe both the dilatancy effects experienced by the granular skeleton during the initial deformations and the rheology of wet granular media when the flow is fully developed. The predictions of the model agree quantitatively with the experiments by Pailha et al (2008) who investigated the role of the initial volume fraction on the dynamics of underwater avalanches. However, this model assumes uniform flow along the slope aligned direction. The extension of this model to non-uniform flow in the framework of thin-layer depth-averaged flows is proposed in Bouchut et al (2016a). Following Iverson et al (2000), these models assume that, when dilation occurs the fluid is sucked into the granular material, the pore pressure decreases and the friction force on the granular phase increases. On the contrary, in the case of contraction, the fluid is expelled from the mixture, the pore pressure increases, and the friction force diminishes. In order to enable the fluid to enter or exit the mixture, Bouchut et al (2016a) proposed a two-layer model with a fluid layer on top of the two-phase mixture layer.

The description of submarine granular flows and generated water waves should combine the description of the two-phase granular flows and its interaction with the water column. Up to now, these processes are described by oversimplified approaches that decouple the submarine flow from the water motion (Le Friant et al 2003) or by two-layer depth-averaged models that account explicitly from the avalanche/water interaction but without describing the dilatation/compression processes within the granular flow (Fernandez-Nieto et al 2008). Furthermore, in that case, the thin-layer approximation should be imposed in the downslope/normal direction for the flowing avalanche while it should be imposed in the horizontal/vertical direction for the description of generated water waves.

Other physical processes such as polydispersity and related granular segregation have been introduced in thin-layer depth-averaged models, successfully reproducing key experimental observations (e.g. Gray and Thornton 2005, Gray and Kokelaar 2010, Gray and Ancey 2011, Wiederseiner et al 2011, van der Vaart et al 2015). As landslide characteristics change during time even for a single event, one of the challenges is to incorporate a fluid phase as well as erosion/deposition, fragmentation and segregation processes within a single model.

Developed in the late 70s (Cundall and Strack 1979), discrete element modeling (DEM) may be seen as an outgrowth of molecular dynamics simulations used in statistical physics. In the past three decades, it has been widely used to study flowing or static granular media. It considers individual grains (usually, for simplicity's sake spheres or polyhedra) interacting through frictional collisions and/or enduring contacts. Interactions between grains are approximated by idealized force models that dissipate energy (see below). The position and velocity of each grain are obtained by integrating Newton's equations of motion. Relevant transport quantities and bulk properties as the velocity field, stress tensor, or local packing fraction can then be computed. Contrary to molecular systems, dissipation is an important characteristic of granular systems. It is thus necessary to employ realistic approximations to model energy loss in interacting grains. To that end, there are essentially two basic approaches: 'Hard' Particle Models—contact dynamics (CD) and event-driven (ED)—which consider grains to be infinitely stiff, and 'soft' particle models (SPM) for which overlapping grains interact at compliant elastic, frictional contacts with dissipative elements providing inelasticity. In the latter case, the interaction is a function of an allowed overlap between colliding grains that model the contact deformation. As an example, the normal force F_n between contacting grains can be approximated by the association of a spring of normal stiffness k_n and a dashpot of normal damping coefficient γ_n:

$$\begin{eqnarray}&&{{F}_{\text{n}}}={{k}_{\text{n}}}\delta +{{\gamma}_{\text{n}}}\frac{\text{d}\delta}{\text{d}t}\end{eqnarray} \tag{ 2.24 }$$

where δ is the grain overlap. In the presence of friction and contrary to the CD method (see below), the Coulomb law has to be regularized to avoid the discontinuity at zero sliding velocity. As an example, the tangential force can be modeled by a tangential spring

$$\begin{eqnarray}&&{{F}_{\text{t}}}={{k}_{\text{t}}}{{u}_{\text{t}}}\end{eqnarray} \tag{ 2.25 }$$

where u_t is the tangential displacement at the contact whose magnitude is truncated to satisfy the Coulomb law |F_t|  ⩽  μ|F_n|, where μ is the friction coefficient (Cundall and Strack 1979).

Contact dynamics (CD), by considering the grains as perfectly rigid bodies, neglects grain deformation. Thus, deformation of the granular medium is exclusively due to grain rearrangements. The volume exclusions of the perfectly rigid grains are formulated by unilateral constraints: the gap between grain surfaces must be larger or equal to zero. Consequently, the number of degrees of freedom depends on the number of aforementioned constraints and thus on the number of contacts that may change with time. For frictional grains, a second type of constraint has to be taken into account: the non-sliding condition if the Coulomb yield threshold is not satisfied. The contact laws are non-smooth, as they express no functional dependence between relative velocities and contact forces. At each time step, iterative methods are used to determine both forces and velocities. The time step should be chosen much smaller (typically 100 times smaller) than the characteristic time of the system (e.g. the inverse of the shear rate). By contrast, in SPM simulations, the number of degrees of freedom does not change in time and integrating the equations of motion requires a time step that is small enough to resolve the grain deformations during collisions. Thus, the time steps used in SPM are generally smaller than those used in CD, leading to larger computational times.

The event-driven method (ED) consists in advancing the system from event to event by repeating the following steps: (i) predict the next collision event between two grains, (ii) bring forward the system to that instant, and (iii) compute the velocities after the collision. This method is very precise because the grain dynamics can be solved analytically between collisions. It is also efficient as the effective time step matches the mean free time. Unfortunately it assumes binary collisions and is thus not suitable to simulate enduring contacts. Therefore ED is relevant to model collision-dominated flows but not systems where enduring contacts prevail. In the latter case, SPM and CD are more suitable.

DEM simulations mimic laboratory experiments. There are advantages of these simulations over experiments since they permit the measurement of any physical quantity without perturbing the system. Moreover, all parameters such as size and mass distribution, friction properties, grain restitution coefficient and eventually rigidity (i.e. Young modulus) can be tuned independently. For example, one can freeze the grain rotation or use perfectly elastic grains, which is difficult, if not impossible in the lab. The attractiveness of these simulations is thus obvious but their ability to reproduce laboratory experiments (i.e. a large number of grains over a long period of time) is recent. This is mainly due to the large amount of computational time needed to carry out such simulations. Simulating a natural avalanche in reasonable times is still impossible because of (i) the very large number of grains involved (several millions or billions) and (ii) their wide and generally unknown disparity of properties, which reduces the computation efficiency. However, it is an inescapable tool at laboratory scale.

DEM can be coupled with a continuum fluid model to describe fluid/particles flows. These hybrid methods have been developed recently for modelling bedload transport in rivers or aeolian sand transport in deserts (Carneiro et al 2011, Duran et al 2012) and can be in principle used for describing debris flows.

The existence of the unjamming transition is fundamental in controlling the behaviour of granular materials. The latter are jammed at rest and can sustain some load, but if a threshold shear stress is exceeded, part of the material starts to flow. In nature, geological processes like landslides or rock avalanches involve an unjamming transition of granular media.

For free surface flows of granular systems driven by gravity, the unjamming transition above a critical shear stress is evidenced by the existence of an angle of maximum stability of a pile: θ_a. It is the angle at which the flow starts; the angle of the pile relaxes then towards the smaller angle of repose θ_r. Parameters like humidity, system dimensions, friction between grains, bottom roughness, and packing fraction can influence the value of θ_a. Many studies have been devoted to these angles, most of them focusing on the succession of avalanches in a rotating drum (Courrech du Pont et al 2003a, 2003b, GDR-MiDi 2004, Fischer et al 2008), or on a continuously fed pile (Jaeger et al 1989, Frette et al 1996). An alternative method for investigating the unjamming transition is to incline gently in the gravity field an initially horizontal box filled with grains (Bolthenhagen 1999, Aguirre et al 2000, 2001, Kiesgen de Richter et al 2012), or a plane on which has been prepared a granular layer of uniform thickness h (Pouliquen 1999, GDR-MiDi 2004). The grains start flowing when the plane inclination reaches a critical angle θ depending on the layer thickness h. The thickness then slowly decreases; the material slows down, stops, and leaves a static layer of material on the bed of thickness h_stop(θ). By starting with different thickness h of the initial layer the curve h_stop(θ) can be obtained. It divides the phase diagram (h, θ) in a region where no flow is possible (h  <  h_stop) and a region where flow may occur (h  >  h_stop). That critical thickness h_stop, which is influenced both by the bulk material and by the bumpiness of the incline (Goujon et al 2003), decreases with the inclination angle (see figure 5): θ₁ corresponds to the angle where h_stop(θ) diverges, its value is generally close to the angle of repose of the corresponding static pile; θ₂ corresponds to the angle where h_stop(θ) vanishes. A similar phase diagram was obtained in numerical simulations (DEM) where periodic boundary conditions are imposed both along the flow direction and in the transverse direction. The bumpy bottom was made of fixed particles (Silbert et al 2001).

Zoom In Zoom Out Reset image size

We will first focus on steady and fully developed flows, referred in the remaining of the paper as SFD flows, i.e. flows whose quantities, such as the flow height, mass flow rate Q, etc, change neither with time (steady) nor with the position in the main flow direction (fully developed). Such systems have yielded a wealth of information useful to understand geophysical flows although the latter are generally unsteady and not fully developed. For sake of simplicity, physicists first focused on dry and nearly monosized grain flows which have two major advantages: (i) segregation is weak or absent and (ii) the effect of the surrounding and interstitial fluid can be neglected as long as grains are large enough. Despite this simple configuration, the observed regimes are numerous and the nature of the boundaries (flat, bumpy or erodible base and sidewalls, etc) is of crucial importance.

Most of the studies of granular gravity-driven flows on inclines reported in the literature deal with flows of grains on bumpy surfaces, both in experiments and in DEM numerical simulations (GDR-MiDi 2004, Delannay et al 2007 and references therein). Such a configuration tends to avoid the inherent discontinuity and sliding observed at a flat base by gluing on the surface of the incline grains of the same nature as those involved in the flow. In these conditions a no-slip boundary condition at the base of the flow is observed, at least for moderate inclination angles. The use of very wide channels in most experiments is an attempt to avoid the influence of confinement due to the lateral boundaries (e.g. sidewalls). For the same purpose, most numerical simulations use periodic boundary conditions in the transverse direction.

One of the main observations for steady and fully developed flows down 'unconfined' bumpy inclines is their limited domain of existence. They have indeed been observed only in a finite range of inclinations and thicknesses. As seen above, the flowing region is bounded by the h_stop(θ) curve in the phase diagram (θ–h). Moreover, a too thick layer and/or a too inclined plane are supposed to lead to accelerated flows (Pouliquen 1999, Silbert et al 2001, GDR-MiDi 2004). However, it is worth noting that the experimental attainment of SFD flows is restricted by the physical length of the chute: a flow may require a length longer to that of the incline to reach its SFD state. Thus the measurement of the maximal angle above which an 'accelerated regime' is observed may depend on the setup used (Brodu 2013). Numerical simulation are subject to the same kind of restriction as the runs are limited in time and the above mentioned maximal angle may depend of this time.

Unidirectional SFD flows possess some noticeable properties (Pouliquen 1999, Silbert et al 2001, GDR-MiDi 2004):

• (a)  
  the volume fraction is constant across the depth of the flow -and independent of the value of this depth- except for the boundary layers near the bottom and the surface;
• (b)  
  for thick flows, the velocity follows the so-called Bagnold profile i.e. it varies with the depth to the power 3/2;
• (c)  
  for thin flows, when h is slightly larger than h_stop (and thus θ slightly larger than θ₁), the velocity profile is close to linear (Silbert et al 2001, Rajchenbach 2003);
• (d)  
  the Froude number defined by $\bar{u}$/(gh)^1/2, where $\bar{u}$ is the depth averaged velocity, is empirically found to be proportional to the ratio h/h_stop. This is consistent with the claim that when h_stop  =  0 (i.e. θ  >  θ₂), the velocity diverges (accelerating flows).

The inclined-plane geometry has an interesting property: for SFD flows, if sidewall friction is negligible, the tangential and normal forces exerted on the base correspond exactly to the components of the flow weight. More generally, if we consider a material slab of arbitrary thickness and bounded by the free flow surface, the tangential and normal forces on the base of the slab correspond also exactly to the components of the slab weight. Their ratio, the effective friction coefficient, is thus always equal to the tangent of the angle of inclination: tanθ. In the frame of the μ(I) rheology described previously (section 2.3), it means that, for flows which are both SFD and invariant in the direction transverse to the flow, the inertial number is uniform in the flow and solely function of the inclination angle. It is thus easy to show that, for those flows, the Bagnold profile results from the μ(I) rheology (GDR-MiDi 2004). The number I can then be related to h_stop(θ) through the aforementioned linear relationship between the Froude number and the ratio h/h_stop. The measure of the number I for an inclination angle θ, gives the function μ(I), as μ  =  tanθ. Note that there is no reason why a SFD flow should not vary in the direction transverse to the flow. If it varies, neither the local ratio of shear to normal stress, nor the inertial number, is expected to remain uniform.

Up to this point, the overall situation appears to be simple, with a limited domain of the parameter space (θ–h) in which unidirectional Bagnold type dense SFD flows occur. As expected, when increasing the inclination angle θ, these flows become unstable. Non-transverse invariant stripe state with spanwise vortices emerges naturally from their destabilization (Börzsönyi et al 2009). Nevertheless, this new state is also SFD, and, as explained before (section 2.3.2), the 3D extension of the μ(I) rheology is not able to describe it.

Besides these studies on flows in wide channels (where side-wall friction is neglected), extensive experimental works have shed light on the crucial influence of sidewalls (Taberlet et al 2003, 2004, 2008, Bi et al 2005, 2006, Jop et al 2005, Richard et al 2008, Brodu et al 2013, 2015) on the flow behaviour. The so-called Sidewall-Stabilized Heap (SSH) is a spectacular effect: at any given inclination angle, above a critical flow rate, a static pile forms along the incline. The pile, whose angle increases with the flow rate, is stabilized by the flow atop it (Taberlet et al 2003). If the flow rate is constant in time, the flow atop this sidewall-stabilized heap is SFD. The influence of sidewalls on SFD flows on top of a static pile has been studied by carrying out experiments in setups of different widths, up to 600 particle diameters (Jop et al 2005).

The properties of flows occurring atop a SSH are different from that observed for SFD flows over bumpy inclines. For instance (i) the main-flow direction velocity profile is no more concave but convex, (ii) the volume fraction varies from zero to the random close packing with a characteristic length that scales with the tangent of the angle of the SSH and thus with the flow rate. Consequently, the gas, liquid and solid behaviours of granular systems are present at the same time in such systems and the corresponding rheology is much more complex than that observed for negligible confinement, especially at high flow rate and/or important confinement. Therefore, although the application of the μ(I) rheology to such flows could remain valid for regions that are neither quasi-static nor too rapid (Jop et al 2006) or for moderate confinement and flow rates, it is unable to capture the full behaviour of the flow in the general case (Jop 2015). Moreover, in a similar context of confined surface flows, but not SFD, we have seen (section 2.3.2) that the 3D extension of the μ(I) rheology failed.

Surprisingly relatively few experimental studies (Savage 1979, Johnson et al 1990, Ahn et al 1991, Louge and Keast 2001) have considered the case of flat frictional inclines. In contrast with a bumpy plane, the mean-flow direction velocity profile involves a considerable slip velocity at the base (Ahn et al 1991, Artoni et al 2012). In their experiments, Louge and Keast (2001) observed SFD flows only for a given range of inclination angles (θ_min,θ_max), independently of the flow height. Here also, the experimental attainment of SFD flows is restricted by the physical length of the chute. Nevertheless, if sidewall friction is negligible, it can be deduced from the force balance that an upper bound exists for θ; its tangent is lower than the friction coefficient between the flat incline and a grain (Brodu et al 2015). The case of θ_min is more complex and its value, when it is defined, depends crucially of the confinement (Brodu et al 2013).

Compared to the bumpy incline configuration, flows on flat frictional surfaces involve a much faster overall velocity, due to the presence of a basal layer of rolling grains, upon which slides the bulk of the unidirectional flow. The bottommost layer of grains is then an effective base for the flow bulk, which follows the aforementioned Bagnold scaling (Delannay et al 2007, Artoni et al 2012, Brodu et al 2013). The analogy with an effective bumpy base extends to the presence, at larger inclination angle, of a convective regime with similar velocity and density profiles. The μ(I) rheology is valid in the case of unidirectional flows but its use in the case of convective regime is more problematic (Brodu et al 2013).

Recent extensive numerical works (Brodu et al 2015) studied flows down flat frictional inclines at high angle of inclination and/or flow rate, in the presence of sidewalls. They show that the system can achieve numerous steady and fully developed regimes which present non-trivial features including heterogeneous volume fraction, secondary flows, symmetry breaking and dynamically maintained order. Some of these regimes have been observed experimentally (Holyoake and McElwaine 2012). Note that the very nature of the incline or of the sidewalls (flat or bumpy) is not crucial to obtain such regimes as long as grain–wall friction prevails over grain–grain friction.

A polydisperse system is more complex to study since it segregates that is, the smallest grains percolate to the bottom, while the largest are squeezed up toward the free surface (Savage and Lun 1988, Wiederseiner et al 2011). This process takes time to reach the steady state and it introduces a spatial heterogeneity within the system. The former point requires very long inclines to reach SFD flows or the use of numerical simulation with periodic boundary conditions in the main flow directions. For instance, Tripathi and Khakhar (2011) propose an adaption of the μ(I) rheology and test it by means of DEM simulations.

Granular/fluid systems in a steady and fully developed state have been investigated through flume experiments (e.g. Armanini et al 2005, 2009, Larcher et al 2007, Kaitna et al 2014). For example, Armanini et al (2005) examined particle/fluid mixture in steady uniform conditions by setting up a closed loop with an inclined channel and a high–speed conveyor belt. Such conditions make it possible to analyze the internal structure of the flow using high-speed video imaging through the sidewalls of the flume. These flows show a stratification of rheological regimes (Armanini et al 2005, Larcher et al 2007): an external region, near to the free surface, dominated by nearly instantaneous contacts among the particles (collisional regime), an internal region dominated by prolonged contacts among the particles (frictional regime) and possibly a static bed in which the particles are immobile. Depending on the ratio of the solid to liquid discharges and on the channel inclination, different flow regimes can be sustained: either flows of the mixture over a rigid, non-erodible base (moderate slope and large volume mixture) or flows over an erodible base (larger slopes) (see for further details the review of Berzi et al (2010) on debris flows). The experimental data drawn from these steady configurations provide invaluable information for testing two phase-flow models.

Another relevant configuration is truly unconfined flows, i.e. flows whose lateral extension is small compared to the corresponding size of the setup. In practice to obtain such kind of flows, grains are released from a reservoir located at the top of the plane by opening a gate whose width is significantly smaller than that of the incline. Depending on the input flux and on the plane inclination, different flow regimes and associated deposit morphologies are observed when a long lasting supply is imposed upslope of beds with inclination angles close to the friction angles of the granular material involved (figure 5). Low flux and/or low slope angle lead to intermittent avalanches and rounded deposits, intermediate flux and slope promote self-channelized quasi-steady flows, and high flux and slope cause a roll-wave regime (Félix and Thomas 2004, Takagi et al 2011). Note that series of adjacent self-channelized flows form when the gate has the same width as that of the incline, and this phenomenon is called granular fingering (Pouliquen et al 1997, Pouliquen and Vallance 1999). In that case small perturbations at the flow front are amplified, which leads to more or less merged granular fingers.

In the intermediate regime, the flow creates its own channel with an almost constant width bordered by lateral levées along the slope and generates a levées-channel deposit as observed for some pyroclastic flows, snow avalanches, debris or submarine flows (see figures 2(c) and (d) and section 4). The question is as to whether the physical processes at the origin of self-channeling and levées-channel deposit are essentially the same for all these geophysical flows. Important issues are to understand and quantify the role of these parameters in the flow dynamics and deposit morphology. For geophysical flows, the challenge is also to find scaling laws relating the deposit morphology to the flow dynamics.

The detailed and comprehensive laboratory experiments of Félix and Thomas (2004) show that self-channeling and levées-channels deposits can be obtained with dry granular flows. They show that self-channeling results from the building up of static lateral zones that develop behind the front during the flow, creating a steady flow in a natural channel of approximately constant width along the slope. The levées-channel morphology appears when the supply of granular material is stopped and the central channel is drained and its thickness decreases. The levées form from the static lateral zones, even though part of these zones collapse during the channel draining (Félix and Thomas 2004, Kokelaar et al 2014). When the imposed flux increases, both the total width of the flow (and deposit) w and the front velocity increase, although the velocity seems to saturate at high flux (Takagi et al 2011), whereas the thickness of the flow h_fl and of the levées h_l only slightly increase and the thickness of the central channel of the deposit h_c remains almost the same. When the slope angle increases, the width w and the front velocity also increase while the morphology is characterized by a smaller thickness of the levées and of the central channel deposit.

Further laboratory experiments and simulations showed that the total width of the lobe w increases slightly with time (Deboeuf et al 2006; Mangeney et al 2007a, Takagi et al 2011). Deboeuf et al (2006) assumed that flows, in their experiments with glass beads, reach asymptotically a steady regime, in the long time limit. On the contrary, Takagi et al (2011) showed that, at long time (after about one hour), glass beads flows become unstable with fluctuations of their width and thickness. On the other hand, long term flows of sand particles reach a steady regime (Takagi et al 2011). These contradictory assumptions and observations lead to different interpretation of the flow in terms of rheology (see section 4 of Takagi et al 2011). This long term behavior, however, does not correspond to the time scale of natural flows that travel and deposit in much shorter time (maximum of ~10 min at the field scale corresponding to 10s of seconds at the laboratory scale if the characteristic time scale is taken as  √(h/g)).

Félix and Thomas (2004) showed that the morphologies of the flow and of the deposit are strongly sensitive to the degree of polydispersity (glass beads with different sizes) while the width and velocity of the flowing central channel seem to only slightly depend on it. In particular the levées-channel morphology is significantly enhanced when the degree of polydispersity increases, i.e. when the size distribution of the grains is wider (figure 6(B)). As polydispersity increases, small and large particles concentrate in the channel and the levees, respectively. With colored particles it is possible to evidence the segregation process: large particles migrate to the flow surface while small particles percolate downward through kinematic sieving and accumulate in the basal part of the flow. Because large particles at the flow surface have larger velocity than average, they migrate forward to the front zone where they (i) are directed to the quasi-static borders where they get trapped, or (ii) get entrained to the flow base and then are segregated towards the surface and recirculated (Pouliquen and Vallance 1999, Félix and Thomas 2004). Based on the segregation model of Gray and Thornton (2005) and on an empirical 3D velocity field resembling the experimental observations of water saturated granular material, Johnson et al (2012) suggest that the segregation process may be more complicated, involving recirculation of large particles in spiral trajectories within the flow head and that are further advected to the flow edges and deposited to form large-particle-enriched levees. Using angular coarse grains together with more rounded fine particles, as it is frequently the case in nature, Kokelaar et al (2014) showed that the run-out distance of self-channelized flows strongly depend on polydispersity. In that case, the material enriched in coarse irregular shaped particles within the lateral borders experience higher friction than the fines-enriched more rounded material flowing within the central channel (Pouliquen and Vallance 1999, Takagi et al 2011, Kokelaar et al 2014). Run-out distance is enhanced as the coarse-grained borders prevent lateral spreading and as the channel consists mostly of finer grains with lower friction.

Zoom In Zoom Out Reset image size

When granular flows are water-saturated, scaling issues are encountered in small-scale laboratory experiments, in particular due to the relatively strong effect of surface tension, pore-fluid shear resistance and cohesion at this scale compared to real situations and to the faster dissipation of excess pore-fluid pressure (Iverson 2015). The large-scale USGS debris flow facility (95 m long, 2 m wide, and 1.2 m deep channel) in Oregon makes it possible to overcome this problem and to mimic real debris flows made of water-saturated mixtures of sand, gravel, and silt or clay. Typical volumes of ~10 m³ are released from rest on a 31° inclined channel connected to a 2.4° unconfined surface (figure 7). The experimental flows reproduce the architecture of debris flows commonly observed in nature, which consist of a dilated, coarse-grained snout with intense particle collisions and of a nearly liquified, fine-grained body with surface roll waves (Iverson et al 2010). When reaching the runout area, self-channelized flows are observed. They are very similar to the self-channelized dry granular flows described above with levées-channel morphology, despite the very different nature of the material involved and the very different slope and input flux (Johnson et al 2012).

Zoom In Zoom Out Reset image size

The presence of a fluid phase, however, may change the segregation processes and the frictional properties of the flow as shown for example by fingering instabilities of propagating granular fronts (Pouliquen and Vallance 1999). Indeed, the presence of a viscous liquid may decrease the efficiency of segregation processes, which inhibits fingering. On the contrary, size-segregation occurs in granular material fluidized with air or water, thus increasing the permeability of the front and lateral quasi-static borders and reducing it in the inner flow so as to enhance retention of the pore-fluid pressure and significantly increase the flow mobility. This effect is promoted by the presence of fine silt and clay particles in water-saturated mixtures where high pore pressure is maintained up to deposition and even during later consolidation (Iverson 1997, Major and Iverson 1999). This change in permeability also adds to the size-segregation of irregular coarse particles and more rounded fine particles to increase the difference of friction between the inner low-friction flow and the high friction borders. Therefore, the morphology of the deposit results from both the flow dynamics and the flow composition, making it difficult to relate the deposit characteristics to the flow properties. Félix and Thomas (2004), however, proposed a relation to estimate the flow velocity from the thickness of the levées and of the central channel (see section 4).

Numerical simulations based on a thin-layer depth-averaged model using the μ(I) rheology at the leading order makes it possible to reproduce qualitatively most of the laboratory of self-channeling flows without taking into account segregation processes (Mangeney et al 2007a). In these configurations, the μ(I) rheology leads to increase the friction at the border of the flow where the thickness decreases rapidly (Mangeney et al 2007a). The simulations reproduce the self-channeling of the flow with quasi-static borders, the slight increase of the total flow width with time, the draining of the central channel when the supply stops, and the resulting levées-channel morphology. The simulated levées-channel deposit looks very similar to the deposit observed in the laboratory for almost monodisperse mixtures of beads (figure 6). The velocity field is well reproduced with a quasi-constant downslope velocity in the central channel, with velocity decreasing and rotating toward the cross-slope direction when approaching the lateral borders of the granular mass (see e.g. figure 17 of Mangeney et al (2007a) and figure 9 of Johnson et al (2012)). Analysis of the downslope and cross-slope forces shows that the balance between the friction force related to the μ(I) rheology and the thickness gradients plays a key role in the rotation of velocities near the front, in the occurrence of a quasi-static flow, and in the generation of static zones just behind the front (figures 13 and 18 of Mangeney et al (2007a)). By coupling a depth-averaged description of the preferential transport of large particles towards the front (Gray and Kokelaar 2010) with a classical thin-layer depth-averaged model, Woodhouse et al (2012) simulated finger formation in a bidisperse granular flow. Numerical simulations reproduce qualitatively the experimentally observed formation of a front rich in large particles, its instability and the subsequent evolution of elongated fingers bounded by large-rich lateral levees. However, the number of fingers increase with the numerical resolution, highlighting the ill-posedness of the equations in the investigated regime related to the lack of important physical processes in the model. The depth-averaged μ(I)-rheology (Gray and Edwards 2014, Baker et al 2016) provide a plausible means of regularizing the model and hence determining the development of segregation-induced fingers.

Simulations suggest that, contrary to the total width w, the width of the central channel w_c (almost equal to the width of the quasi-steady flow) is almost constant for a given flux and slope. By assuming that the flow thickness is about the same as the levées thickness and that the thickness of the central channel is about equal to h_stop, Mangeney et al (2007a) suggest that the width of the central channel is proportional to the input flux as found in the experiments of Takagi et al (2011), making it possible to estimate the flux (or the velocity) from the measured heights of the central channel and lateral levées of the deposit.

These experimental and numerical studies suggest that the mechanisms governing self-channelling and levée formation could be generic to a wide variety of geophysical granular flows, involving or not polydispersity and/or a fluid phase. In other words, grain size segregation and accumulation of coarse particles at flow surface and margins is not necessarily the cause of self-channelling (though it may contribute to enhance it). Note, however, that grain size segregation seems to be a prerequisite for granular fingering (Pouliquen et al 1997, Pouliquen and Vallance 1999, Woodhouse et al 2012).

However, up to now, simulations of self-channeling flows and levees formation or analysis based on simple models as those described above were compared essentially qualitatively to experimental results. To go further in the understanding and constrain of the physical processes at work, quantitative comparison of simulated and observed thickness and velocity profiles should be performed. Accurate simulation of these processes would require two-phase models describing both the mass and momentum conservation of the solid and fluid constituents and their interaction coupled with segregation models. A key issue is to better describe the behavior at small Froude and inertial numbers corresponding to the regime experienced by the border of the flow (figure 14 of Mangeney et al (2007a)).

3.4.1. Insight from laboratory experiments.

Unsteady flow generated from the quasi-instantaneous release of a granular column (dam-break configuration, figure 8) was considered in many experimental studies over the last decade. Such experiments represent simple analogues of geophysical granular flows and can also serve as tests for numerical simulations. Studies involved various types of granular materials with negligible (granular flows) to high internal pore fluid pressure (gas-particle flows and liquid-particle flows), and were conducted in unidirectional (2D) or axisymmetrical (3D) configurations on horizontal or inclined beds. They focused on flow kinematics and on scaling laws relating the morphological characteristics of the flow deposits (run-out, height) to the geometry of the initial column released.

Zoom In Zoom Out Reset image size

Experiments show that the flow duration is controlled by a timescale proportional to that of particle free-fall, t_f  =  (h₀/g)^1/2 where h₀ is the initial column height. Emplacement occurs in three phases: the flow front accelerates after release, then propagates at nearly constant velocity in some configurations, and finally decelerates until motion ceases (Lajeunesse et al 2004, 2005, Lube et al 2004, 2005, Balmforth and Kerswell 2005, Roche et al 2008). The run-out distance of dry flows is travelled by material originating from the outer surface of the upper third of the initial column (Thompson and Huppert 2007). Only part of the granular mass flows while particle sedimentation occurs as the flow propagates, resulting in a basal static region overridden by an upper and frontal part, both being separated by an interface that propagates upwards until the flow is consumed (Lajeunesse et al 2005, Mériaux 2006, Doyle et al 2007, Lube et al 2007, 2011, Girolami et al 2008, 2010, Farin et al 2014, Ionescu et al 2015, Lusso et al 2016a). Gas-particle flows with high pore pressure notably propagate like inertial fluid gravity currents at velocity U ~ (2gh₀)^1/2 during most their emplacement, despite their high particle concentration, before they enter a frictional granular regime once pore pressure has decreased sufficiently through diffusion (Roche et al 2008).

Experiments on horizontal slope show that the characteristics of the final deposit depend primarily on the column aspect ratio defined as

$$\begin{eqnarray}&&R=\frac{{{h}_{0}}}{{{x}_{0}}}\end{eqnarray} \tag{ 3.1 }$$

where x₀ is the initial column length (see table 3 and references therein). The deposit shape varies from truncated cone, with slope close to the repose angle of the granular material at low R, to cones with slope down to a few degrees at high R (Lajeunesse et al 2004, 2005, Lube et al 2004, 2005, Balmforth and Kerswell 2005, Siavoshi and Kudrolli 2005). The normalized run-out distance

$$\begin{eqnarray}&&{{x}^{\ast}}=\frac{{{x}_{\text{f}}}-{{x}_{0}}}{{{x}_{0}}}\end{eqnarray} \tag{ 3.2 }$$

and deposit height

$$\begin{eqnarray}&&{{h}^{\ast}}=\frac{{{h}_{\text{f}}}}{{{x}_{0}}}\end{eqnarray} \tag{ 3.3 }$$

where x_f and h_f are the deposit length and height (figure 8), are used to define scaling laws of the form

$$\begin{eqnarray}&&{{x}^{\ast}}=\lambda {{R}^{n}}\text{,} ~ {{h}^{\ast}}={{\lambda}^{\prime}}{{R}^{{{n}^{\prime}}}}\end{eqnarray} \tag{ 3.4 }$$

though a clear law for h^* has not been defined yet. Note that h₀ instead of x₀ can be used as well to normalize these parameters, so that x^*  =  λRⁿ⁻¹ (resp. h^*  =  λ'R^n'−1). Data for h^* do not show systematic trends, but those for x^* reveal that the value of the coefficient λ and of the exponent n depend respectively on the material properties and on the initial configuration (table 3). Almost all studies show that n  =  1 (i.e. x^*  =  λR) at low R, typically less than 1.6–3, so that ${{x}_{\text{f}}}-{{x}_{0}}=\lambda {{h}_{0}}$. In other words, the distance travelled by the flow from the edge of the initial column is simply proportional to the initial column height. At higher R, n  <  1 for unidirectional flows but n approaches 1 as the channel width increases and side effects become weaker (Balmforth and Kerswell 2005, Mériaux 2006). Results of Lacaze et al (2008) are a notable exception, but their experiments involve flows only one particle wide whose dynamics might differ fundamentally from that of flows that are typically hundreds of particle diameter wide. In the axisymmetrical configuration n  =  1/2 at high R, which may reflect spreading of the granular mass (Lajeunesse et al 2004). Farin et al (2014) recently showed that these scaling laws are not valid for granular collapse over slopes higher than about one half of the friction angle of the material involved (θ_c). In that case, a last regime of slow, almost steady propagation occurs after the initial acceleration and deceleration phases typical of granular collapses over smaller slopes. In that case the normalized quantities vary also with the volume of the released mass.

Table 3. Scaling laws in experiments on horizontal dam-break granular flows.

Configuration reference	x*  =  (xf  – x0)/x0	h*  =  hf/x0
Unidirectional
Balmforth and Kerswell (2005)	∝R0.55±0.05 (w  =  1 cm, R  >  2)	∝R0.45±0.05 (R  >  1.5)
	∝R0.9±0.1 (w  =  20 cm, R  >  2)
Lajeunesse et al (2005)	∝R (R  <  3)	∝R (R  <  0.7)
	∝R2/3 (R  >  3)	∝R1/3 (R  >  0.7)
Lube et al (2005)	1.6 R (R  <  1.8)	∝R2/5 (R  >  1.15)
	2.2 R2/3 (R  >  2.7)
	1.2 R (R  <  1.8, symmetric)
	1.9 R2/3 (R  >  2.7, symmetric)
Lacaze et al (2008)a	∝R (R  <  2.2)	∝R0.45±0.02 (d  =  2.5 mm, R  >  1)
	∝R0.72±0.03 (d  =  2.5 mm, R  >  2.7)	∝R0.39±0.02 (d  =  5 mm, R  >  1)
	∝R0.83±0.03 (d  =  5 mm, R  >  2.7)
Mériaux (2006)b	0.8 R (sand, R  <  2)	R (R  <2)c
	R (glass beads, R  <  2)	1.1 R0.45±0.05 (sand, R  >  2)
	∝R0.70±0.05 (sand, R  >  2)	R0.45±0.05 (glass beads, R  >  2)
	∝1.3 R0.70±0.05 (glass beads, R  >  2)
Axisymmetric
Lajeunesse et al (2004)	∝R (R  <  3)	R (R  <  0.74)
	∝R1/2 (R  >  3)	0.74 (R  >  0.74)
Lube et al (2004)	1.24 R (R  <  1.7)	R (R  <  1)
	1.6 R1/2 (R  >  1.7)	0.88 R1/6 (1.7  <  R  <  10)
		∝R−y (R  >  10)d
Roche et al (2011)	d  =  330 µm	d  =  330 µm
	1.40 R (R  <  1.6)	R (R  <  0.7)
	1.82 R1/2 (R  >  1.6)	~0.7 (0.7  <  R  <  7)
	1.48 R (fluidized, R  <  1.6)	∝R−1/2 (R  >  7)
	1.99 R1/2 (fluidized, R  >  1.6)	0.26 R1/3 (fluidized, R  <  5)
		R−1/2 (fluidized, R  >  5)
	d  =  80 µm	d  =  80 µm
	1.61 R (R  <  1.6)	∝R (R  <  1)
	1.84 R2/3 (R  >  1.6)	~0.95 (1  <  R  <  6)
	2.49 R (fluidized, R  <  1.1)	∝R−1/2 (R  >  6)
	2.53 R1/2 (fluidized, R  >  1.1)	~0.1–0.2(fluidized)

Experiments involved coarse particles of typical grain size d  >~0.3 mm, unless specified. In unidirectional experiments, w is the width of the channel. ^aw  =  1.2d, ^bslow gate withdrawal, ^cthe exponent n is unknown, ^drecalculated from h₀/h_f.

Other things being equal, parameters causing run-out increase (i.e. higher λ) include the aspect ratio, the slope angle (Hogg 2007, Mangeney et al 2010, Farin et al 2014), the volume of material released and the presence of on erodible granular substrate (Mangeney et al 2010, Farin et al 2014), the material polydispersity with possible related friction reduction effects (Roche et al 2005, Phillips et al 2006, Meruane et al 2012), and high initial pore fluid pressure (Roche et al 2011). For flows within a viscous ambient fluid, low initial particle volume fraction dramatically increases the runout distance because pore pressure is generated due to fluid expulsion as the flowing material contracts. In contrast, initially dense flows with volume fraction close to the maximum packing travel about the same distance as dry flows (figure 11(b) of Rondon et al 2011). The related mechanisms and their influence on the scaling laws are yet to be investigated in detail; in contrast, run-out decrease (lower λ) is caused by both slow material release (i.e. low inertia, Mériaux (2006), increased material friction angles (Balmforth and Kerswell 2005, Roche et al 2011) and interparticle cohesion (Mériaux and Triantafillou 2008). The presence of a thin layer of erodible material lying on an inclined bed may increase the maximum runout distance of a granular avalanche flowing down the slope by up to 40% and change the flow regimes at slopes higher than the critical angle θ_c (Mangeney et al 2010 and Farin et al 2014). In that case, as for granular flows on steep slopes, a slow, almost steady propagation phase occurs. Note that there could be a variation of a few degrees of this critical angle as a function of the initial and boundary conditions, but this has not been demonstrated yet.

3.4.2. Insight from analytical and numerical models.

Analytical solution of the thin-layer depth-averaged equations (Mangeney et al 2010, Faccanoni and Mangeney 2013) and numerical simulations show that the coefficient λ varies inversely with the difference between the tangent of the slope angle and that of the material friction angle (Balmforth and Kerswell 2005, Kerswell 2005, Mangeney-Castelnau et al 2005, Mangeney et al 2010):

$$\begin{eqnarray}&&\frac{{{x}_{\text{f}}}}{{{x}_{0}}}=\frac{2k}{\tan \delta -\tan \theta}\end{eqnarray} \tag{ 3.5 }$$

with k an empirical parameter (k ~ 0.5 see Mangeney et al 2010). This contradicts the simple model proposed by Lajeunesse et al (2005) where the runout distance is proportional to the friction coefficient μ  = tanδ. This dependence (3.5) is particularly well illustrated by gas-particle flows with high pore pressure (and hence low internal friction) whose run-out distance is larger (i.e. higher λ values) than that of their dry granular counterparts (Roche et al 2011 and figure 10(B)).

While simulations using DEM and depth-averaged thin layer models reproduced the above-mentioned experimental scaling laws, very few quantitative comparisons with the deposit extent (i.e. exact value of λ) and the dynamics of the spreading have been performed. Though thin-layer models reproduced the deposit extent with friction angles only a few degrees higher than the typical friction angle of the material involved (Kerswell 2005, Mangeney-Castelnau et al 2005) they overestimate the spreading velocity, in particular during the first instants after release. This is notably due to the neglect in thin-layer models of the vertical acceleration, which is particularly important in the beginning of the spreading. On the other hand, 2D DEM (contact dynamics) strongly overestimates the runout distance (up to 40%) even with high grain/grain friction (Staron and Hinch 2005, 2007). To overcome this issue, 3D DEM (molecular dynamics) by Girolami et al (2012) included a dissipation related to the resistance to rolling to match the experimental results. Lacaze and Kerswell (2009) were able to reproduce quantitatively the runout extent using 3D DEM with a grain/grain friction coefficient of 0.5. Using a set of parameters typical of the 2D visco-plastic flow law μ(I) to compare with 2D DEM, Lagrée et al (2011) overestimated the runout distance of granular collapse observed experimentally. However, 2D continuum simulation using pressure and/or strain-rate dependent viscoplastic (Ionescu et al 2015) or elasto-viscoplastic (Crosta et al 2009) flow laws were able to reproduce quantitatively the deposit of granular columns on horizontal and inclined planes, and even the detailed dynamics by taking into account the lateral wall friction (Ionescu et al 2015, Martin et al 2016). These simulations show that a visco-plastic flow law derived from the μ(I) rheology makes it possible to reproduce quantitatively the transient regimes characterizing granular collapses over horizontal and inclined beds with rheological parameters derived from the experiments, without any fitting procedure. However, it has to be mentioned that for collapses of small aspect ratio typical of real landslides, a constant viscosity (not dependent on pressure and strain rate) defined in the sense of Ionescu et al (2015) gives results very similar to those obtained with the variable viscosity derived from the μ(I) rheology (Ionescu et al 2015). Note that the initial opening of the gate has a significant effect on the spreading while it does not change the final deposits (Ionescu et al 2015), as also observed for solid-fluid mixture dam-break flows (Iverson and George 2014). This has to be taken into account for quantitative comparison between simulation and experiments.

During granular collapse over rigid or erodible beds, part of the material is flowing over static material. This is qualitatively reproduced by DEM and continuum 2D models (e.g. Lagrée et al 2011, Ionescu et al 2015) as well as by a multilayer shallow models that was able to reproduce for the first time the increase of runout distance of granular flows with increasing thickness of the underlying erodible bed (Fernandez-Nieto et al 2016). Analytical solution and numerical simulation of the time change of the static-flowing interface shows that its dynamics is closely related to velocity profiles and to the rheological properties of the material (yield, viscosity) (Lusso et al 2016a, 2016b) as discussed in section 2.5. In particular, the penetration of the static-flowing interface within the erodible bed is qualitatively reproduced in the continuum 2D models (Mangeney et al 2007b, Crosta et al 2009, Lusso et al 2016a). Quantitative comparison is missing, however.

Waves that look similar to surges reported for debris flows (e.g. Iverson 1997, Zanuttigh and Lamberti 2007) have long been observed in dry granular flows (Savage 1979) (figure 9). These individual waves or series of waves have different characteristics, depending on the initial and boundary conditions and on the material involved.

Zoom In Zoom Out Reset image size

Soliton-like waves are observed for dry granular flows on rigid or erodible slopes with inclination angles close to the friction angles of the material involved (figure 9). Pouliquen and Forterre (2002) showed that, when a 3D granular cap is released on such an inclined plane covered by a thin erodible bed of thickness h_stop made of the same granular material, a triangular shaped traveling wave that propagates at a constant downslope velocity puts into motion static material from the erodible bed at the front (erosion) and leaves material at the rear (deposition) (figure 9(a)). Triangular shaped erosion/deposition waves have also been observed when a perturbation (Daerr 2001, Malloggi et al 2006 (in air or water)) or a constant flux of material (Börzsönyi et al 2005, 2008) are imposed on top of a thin granular layer lying over an inclined plane. The morphology and dynamical properties of these waves or avalanches are qualitatively and quantitatively different for smooth glass beads compared to irregular-shaped granular materials such as sand (Forterre and Pouliquen 2003, Börzsönyi et al 2005, 2008) (figures 9(d) and (e)). In particular, Börzsönyi et al (2008) showed that distinct, well-resolved avalanches made of rough non-spherical grains are bigger and travel faster than those made of spherical glass beads. For inflow conditions close to the minimum depth h_stop required for steady uniform flow to exist, a series of steadily travelling erosion/deposition waves develop (Edwards and Gray 2015). Forterre and Pouliquen (2003) showed that sand flows are much more unstable than glass beads flows in their setup configuration and that they generate roll-waves (figure 9(f)). These waves have smaller amplitudes than erosion/deposition waves, form in slightly deeper flows, and develop at a much larger distance from the perturbation. The most unstable modes (with large amplitude) of roll-waves occur close to the threshold conditions of the flow (low inclinations and thin flows). The generation of large-amplitude waves near the threshold of the flow may have strong implication for natural flows (see section 4).

These waves provide a unique test condition for rheological models. Indeed, the precise characteristics of the wave development (amplitude, shape, velocity, stability threshold, growth rate, cutoff frequency, etc) dramatically depend on the rheological properties of the material. Using the 'partial fluidization model' proposed by Aranson and Tsimring (2002) that involves a viscosity proportional to the square root of the pressure and that describes explicitly the static-flowing transition through a phenomenological order parameter, Mangeney et al (2007b) reproduced qualitatively the steady traveling erosion/deposition wave obtained experimentally by Pouliquen and Forterre (2002). Using the same model, Aranson et al (2006) identified a family of propagating soliton-like avalanches with shape and velocity controlled by the inclination angle and the thickness of the substrate. At high inclination angles, the solitons display a transverse instability, followed by coarsening and fingering similar to the experimental observation of Malloggi et al (2006). These simulations, however, only reproduce qualitatively the experimental results (see e.g. figures 9(b) and (c)). Analysis of these waves has also been done with depth-averaged thin-layer models, where the whole layer of grains at a given downslope position is implicitly assumed to be either completely mobile or completely static. Pouliquen and Forterre (2002) showed that depth-averaged thin layer models with the μ(I) rheology implemented at leading order does not reproduce the erosion wave they observed experimentally (figure 9(b)). Using the same model, Forterre and Pouliquen (2003) showed that roll waves observed for glass beads and for sand particles result from the same instability mechanism and that their different behaviour could be explained by the quantitative differences in the friction coefficient of the μ(I) rheology. This model, however, is unable to predict the observed cutoff frequency. The cutoff frequency results from the dissipation due to the longitudinal velocity gradients (the extensional viscosity term) that are neglected in leading order thin-layer models (Forterre 2006). By adding a depth-averaged viscous term (Gray and Edwards 2014), based on the μ(I)-rheology, into the model, Edwards and Gray (2015) showed that it was possible to predict the formation of erosion-deposition waves similar to those shown in figure 9, which have completely static zones between steadily traveling mobile regions. In particular, this theory was able to accurately capture the growth, coarsening dynamics and typical wave amplitudes and lengths, which is an important step forward in the goal to incorporate erosion and deposition into depth-averaged theories. Börzsönyi et al (2005, 2008) use also depth-averaged thin layer models with the leading order μ(I) rheology to explain the difference between the avalanches generated by glass beads or sand, no quantitative comparison with the avalanche shape and velocity was achieved.

Waves at the interface between a flowing granular mass and an initially static erodible bed were also observed experimentally (Mangeney et al 2010, Farin et al 2014). They could have tremendous effect on the amount of bed material entrained (Rowley 2011). Modelling these waves could be a unique way to investigate rheological laws near the static-flowing transition. Furthermore, a challenging issue is to quantify the role of interstitial fluid and grain polydispersity on erosion/deposition and roll waves development to ultimately compare them with natural surges (e.g. Iverson and Denlinger 2001, Zanuttigh and Lamberti 2007, Zanuttigh and Ghilardi 2010). Indeed, natural surges are strongly heterogeneous as they are generally composed of a first surge, mainly containing boulders, and of possible secondary waves that are generally muddier as shown in large scale experiments (see section 3.3 and figure 7).

Experiments of granular flows over complex topographies provide useful tests for numerical models (see Pudasaini et al 2007 for a review). Simulations of these experiments attempted to investigate the anisotropy of normal stresses in granular flows. Considering this anisotropy, however, did not permit reproducing quantitatively the experiments (Pirulli et al 2007), while adding different coefficient for along-slope and transverse normal stress anisotropy seemed to improve the results (McDougall and Hungr 2004). Quantitative agreement with the experiments instead required to add an ad-hoc variation of the friction angle towards the tail of the flow (Gray et al 1999). In similar experiments by Iverson et al (2004), avalanches of dry sand were confined within a sloping irregular 'valley' before they stopped on a horizontal plane. Using a depth-averaged thin layer model that took into account anisotropy of normal stresses and an approximation of the vertical acceleration that compensates the thin-layer approximation made in the horizontal/vertical reference frame (see section 2.6), they were able to simulate some characteristics of the flows, with friction coefficients measured independently with a tilt table apparatus. Details of the simulated dynamics and deposit were however quite different from those in the experiments (Iverson et al 2004, their figures 9 and 10). Using a similar model that took into account all the curvature terms of the topography derived from the thin-layer approximation applied in a reference frame tangent to the topography, Lucas et al (2014) (see their figure 6 in the supplementary material) also reproduced these experiments but without including anisotropy of normal stresses. As a result of the complex interplay between topography and rheological effects and because of the difficulty of properly taking into account topography in depth-averaged thin-layer models, more detailed comparison of the simulated and observed dynamics and deposit should be performed to be able to really discriminate relevant rheological laws.

Granular flows around obstacles (Tai et al 2001, Gray et al 2003, Hákonardóttir and Hogg 2005, Gray and Cui 2007, Vreman et al 2007, Faug et al 2008, Benito et al 2012) are also of paramount importance given the applications to defensive structures (e.g. deflecting and catching dams) aiming to protect people and infrastructure from hazardous natural granular flows. As an example, Benito et al (2012) showed that presence of a forest of cylinder obstacles significantly increase the stability of a granular layer. Moreover, such kind of flows is also interesting from a fundamental point of view. Indeed they generate shock waves at which there are rapid changes in the flow properties (thickness, velocity etc). Yet, they can be used as benchmarks to test theories, models or numerical simulations aiming to describe granular flows. It has been recently shown that a thin-layer depth-averaged model assuming isotropy of normal stresses and a constant friction coefficient (Cui and Gray 2013) captures experimental results in a satisfactory way. A precise description of such flows is out of the scope of the paper and the interested reader should refer to the aforementioned references for more details.

Field observations provide important constraints on the initial and final conditions (volume and nature of the material involved, type and duration of the mass release, triggering process, deposit morphology...), on the material properties, on the boundary conditions (nature of the ambient fluid, underlying topography, type of substrate...), and on the mechanisms that control the dynamics and deposition of geophysical flows.

Field data mostly consist of measurements of the morphology and extent of the deposits, provided surface erosion was not significant since the time of deposition (table 4). Very high resolution digital elevation models (DEM), i.e. topography data, with precision down to 1 meter are obtained through the increasing precision of satellite imagery on Earth and on the planet Mars (e.g. Huggel et al 2005, Lucas et al 2011, 2014). In fact, free satellite images of the surface of Mars that has a small erosion rate and no vegetation provide unique data on landslide deposits that are hardly obtained on Earth (Quantin et al 2004, Mangold et al 2010, Lucas et al 2011, 2014, Jouannic et al 2012). In subaerial environment on Earth, the shape of the front deposit as well as lateral levées may be measured using terrestrial or airborn scanner laser measurements (LiDAR) reaching a precision of about 2 cm at about 1 km of distance (Conway et al 2010, Sovilla et al 2010, Jaboyedoff et al 2012, Jessop et al 2012). In submarine environment, bathymetry data with a resolution down to 50 cm have been obtained locally with a multi-beam echo sounder mounted on a Remotely Operated Vehicle (ROV) and have revealed detailed morphological features (figure 2(d) and Cannat et al 2013). Sidescan sonar is well suited to identify subaqueous morphological features of deposits of turbidity currents (Piper et al 1985). On planetary surfaces, thermal infrared imagery is used to discriminate landslide deposits from the surrounding substrate, especially at night time where the temperature largely depends on the thermal inertia of the surface that is mainly related to grain size and degree of induration. The volume of material released is determined readily if the deposit is well preserved and if the underlying topography can be constrained accurately. In cases of severe erosion and/or the topography is unknown, the volume can be inferred from the geometry of the source zone from where the flow initiated. Accurate determination of this geometry, however, is commonly challenging because of multiple collapse events (e.g. Lucas et al 2011). Successive topography or bathymetry measurements before and after specific events makes it possible to better constrain the volume involved and to precisely measure the amount and spatial distribution of erosion and deposition (Conway et al 2010, Le Friant et al 2010).

Direct measurement of the flow properties during propagation, however, is very difficult because of the unpredictability and destructive power of the natural phenomena. Yet, parameters such as flow velocity and basal stresses have been quantified from direct measurements (e.g. Arattano and Marchi 2008, McCoy et al 2013) or inversion of seismic data (e.g. Brodsky et al 2003, Favreau et al 2010) as discussed in the following. Subaerial flows are obviously more investigated than their subaquous equivalent. Debris flows and snow avalanches are the most documented events (e.g. Marchi et al 2002, McArdell et al 2007, Vriend et al 2013). However, the number of studies on rockfalls, rock avalanches and pyroclastic density currents involving acoustic or seismic signals is increasing (e.g. Deparis et al 2007, Ripepe et al 2009, Hibert et al 2011, 2014a, 2014b, Levy et al 2015, Farin et al 2015). National observatories record series of complementary data (video, seismic, acoustic) on gravitational flows coupled with pluviometry, ground deformation, and soil properties measurements that are made freely available (e.g. the Observatoire Volcanologique du Piton de la Fournaise (OVPF) where more than 20 years of data are available, or the Observatoire Muldisciplinaire des Instabilités de Versants (OMIV)). Table 4 summarizes the parameters that can be measured in the field and the techniques used for such measurements. Examples of field measurement of parameters are presented below, and the principles and main limitations of the methods employed are highlighted.

The flow velocity is measured by various means. A direct approach to determine the flow surface and front velocity consists of analysis of videos (Turnbull and McElwaine 2007, Doyle et al 2011). However, dense basal pyroclastic flows and snow avalanches, for instance, are often obscured by an upper dilute ash or powder cloud, respectively, which prevents measurement of the underlying dense flow velocity. Doppler speedometers are used to measure the surface velocity whereas the mean front velocity can be calculated from pairs of ultrasonic sensors or geophones at a known distance along the pathway, as done for instance for debris flows (see Arattano and Marchi 2008 and references therein). The velocity of snow avalanches can be deduced from pressure data obtained with pitot-type devices (Nishimura and Ito 1997). In that case, however, a correction must be made according to the size of the piezometer tube. New high resolution radar measurements permit to penetrate through the powder cloud and measure the velocity of the dense snow avalanche front over its entire path with spatial resolution smaller than 1 m (Vriend et al 2013). Correlation of the signals from an array of infrasound sensors or of seismometers makes it possible to localize the source of radiated energy as a function of time and to deduce the flow velocity of pyroclastic density currents, snow avalanches, and rockfalls or individual blocks (Vilajosana et al 2007, Ripepe et al 2009, Lacroix et al 2012). As strong winds may generate noise when using infrasound sensors, these have to be buried and equipped with pipes connected to the atmosphere (Ripepe et al 2009). Flow duration can be deduced from seismic data. Indeed, combined photogrammetric and local seismic data (source-station distances of ~1 km) have shown that the flow duration of rockfalls of volume of 10¹–10⁴ m³ almost corresponds to the duration of the generated seismic signal (Hibert et al 2011, 2014a). Similar observations have been made at larger scale (source-station distances of 30–300 km) for bigger events (1–50  ×  10⁶ m³, Favreau et al 2010). Various techniques have been elaborated to measure the flow depth. These include particularly radar and wire sensors (Arattano and Marchi 2008) as well as laser devices for debris flows generated in large-scale experimental facility (Major and Iverson 1999, Iverson et al 2010).

Basal forces, related to normal and shear stresses as well as particle impacts, are measured with load cells or piezoelectric sensors in debris flows (Major and Iverson 1999, McArdell et al 2007, Iverson et al 2010, McCoy et al 2013) or in snow avalanches (Nishimura and Ito 1997, McElwaine and Turnbull 2005), respectively. The basal normal stress permits to infer the flow depth or the bulk density if the other parameter is known, whereas the ratio of the normal to the shear stress gives the Coulomb basal friction while the impact force indicates the flow speed if a flow density is assumed. Note, however, that inertial forces related to topography curvature, if present, can significantly affect the basal stresses making it more difficult to relate them to the flow depth or bulk density (Favreau et al 2010, Moretti et al 2015). The interstitial fluid has a fundamental influence on the dynamics of many geophysical flows. For instance pore pressure measurements reveal a stagnation pressure in front of snow avalanches while the flow head generates a relative underpressure, less than that of the ambient atmosphere (Nishimura et al 1995, McElwaine and Turnbull 2005). Basal excess pore fluid pressure measured with appropriate sensors indicates either the amount of the weight of the particles supported in flows close to maximum packing or the particle concentration in more dilute mixtures. Excess pore pressure generated by soil contraction (Iverson et al 2000) and/or heavy rain falls (Montgomery et al 2009) can be responsible of flow initiation as frictional forces are severely reduced. High pore pressure is also measured during emplacement of debris flows, hence providing an explanation for their high mobility (Major and Iverson 1999, McArdell et al 2007, Iverson et al 2010, Doyle et al 2011). The temporal evolution of the basal forces can also be recovered from seismic data. Indeed, as subaerial or submarine geophysical flows travel down slope they apply stresses to the ground, which generate seismic waves that can be recorded by regional seismic networks and even global networks for very big events (e.g. Kanamori and Given 1982, Kawakatsu 1989, Brodsky et al 2003, La Rocca et al 2004, Favreau et al 2010, Lin et al 2010). As an example, the 50x10⁶ m³ Mount Steller landslide, Alaska, generated a vertical ground velocity of about 10⁻⁴ m s⁻¹ at 37 km from the source and was recorded by stations located more than 1500 km away (Moretti et al 2012, Zhao et al 2015). The low frequency (0.01–0.05 Hz) basal force can be calculated by inversion of seismic data using the Green's function of the Earth between the source and the seismometer (e.g. Moretti et al 2012, 2015, Yamada et al 2012, 2016, Allstadt 2013, Ekström and Stark 2013). As an example, the maximum force associated with the Mount-Steller landslide is 2  ×  10¹¹ N. Note that broadband seismometers are required for source inversion because short-period instruments only measure frequencies higher than 1 Hz.

When a sufficient number of events can be recorded at a same site, the ratio between the loss of potential energy and the high frequency (>1 Hz) radiated seismic energy can be estimated (Hibert, et al 2011, Levy et al 2015). This energy ratio typically ranges in between 10⁻³ and 10⁻⁵, in agreement with experimental results of acoustic emission generated by block impacts (Farin 2015, Farin et al 2015, 2016). As a result, estimates of the rockfall volume can be calculated from the radiated seismic energy using simple formula derived from granular flow or impact modeling (Hibert et al (2011) for granular flows and Farin et al (2015) for rock impacts). For a given site, purely empirical relations can be also derived from rockfalls seismic signals and independent measurements of the volume involved (Deparis et al 2007, Dammeier et al 2011). The parameters of these relations depend on the material involved and on the local topography and ground structure. Attempts have also been made to relate the seismic signal frequency to interparticle and flow-substrate interactions. The frequency of the seismic signals associated to debris flows is assumed to reflect their solid concentration so that sliding frictional dense flows should produce lower frequency signals than more dilute flows with more collisions (see Doyle et al (2011) and references therein). For debris flows, McArdell et al (2007) have inferred the number of impacts of particles with a diameter typically larger than a few centimetres. This method, however, is very difficult to apply to dense dry granular flows as suggested by the experiments of Farin (2015).

Other parameters can be measured in the field. Substrate erosion by the flows is investigated using buried radars (Sovilla et al 2006, 2010), scour sensors consisting of an array of erodible elements (Berger et al 2011) or stress and pore fluid pressure sensors (McCoy et al 2012), whereas sampling in debris flows and powder snow avalanches is a direct way to know the local particle concentration and/or the grain size distribution (Doyle et al 2011, Rastello et al 2011). Scanner laser measurements may also provide the size and shape of individual grains as well as rugosity of the surface deposit (Haas et al 2012). Grain size, particle density and morphology, componentry and textural analysis of the deposits in aerial environment allow to investigate flow physical processes such as particle segregation, fragmentation and abrasion, as well as substrate erosion and melting due to friction (e.g. Calder et al 2000, Goren et al 2010, Imre et al 2010, Manga et al 2011). Petrographic analysis at both macroscopic and microscopic (e.g. x-ray diffractometry, Mössbauer spectroscopy on thin sections) scales makes it possible to quantify the conditions and characteristics of frictional melting of the material involved (Weidinger et al 2014).

Further field measurements are required to better understand, for instance, the rate of entrainment of the surrounding fluid, the rate of mass exchange between a dense underflow and its upper dilute parts, the amount of erosion of the substrate, and the onset of flow deposition. Substrate erosion is particularly critical in controlling the flow runout distance, as shown by recent experimental studies on either dry or water-saturated flows (Mangeney et al 2010, Iverson et al 2011, Farin et al 2014), and this issue should become the focus of more studies in the next future. Measurement of particle fragmentation in natural flows is also an open and key issue. Beyond the measurements of the properties of individual landslides, seismic data provide a unique tool to detect, localize and characterize gravitational flows. This provides a new way to quantify and understand the potential link between gravitational activity and external forcing. As an example, recent analysis of rockfall events together with volcanic seismicity at Piton de la Fournaise (La Réunion island) suggests that the volume of rockfalls increases just before magma intrusions reach the surface, hence providing an indicator for eruptions (Hibert et al 2014a). Satellite imagery of series of landslides provide complementary data that makes it possible, for example, to quantify links between gravitational instabilities and large earthquakes (Meunier et al 2007).

Before illustrating how laboratory experiments can be used to interpret field observation, let us recall the fundamental issues related to this approach. First, landslide dynamics and deposits are strongly affected by the natural topography, which is generally complex and variable from site to site. On the contrary, laboratory experiments are generally performed in simple configurations. As a result, quantitative interpretation of natural flows in terms of physical processes and rheological behaviour based on laboratory experiments is commonly very difficult. One way to make the link between laboratory experiments and field data is to perform numerical models. The aim is to design models that reproduce quantitatively laboratory experiments and then to apply these models to natural flows by taking into account the real local topography. This approach is illustrated in section 4.2.1 below through the challenging question of the high mobility of large landslides. However, both laboratory and numerical models always represent a huge simplification of natural processes. Indeed, natural flows result from a complex interplay between the effects of external forcing, heterogeneity and polydispersity of the material involved, presence of fluids, initial and boundary conditions, nature and entrainment of the substrate, etc. Even when the physical laws describing these processes may be estimated, real field conditions defining the parameters involved in these laws are rarely known and may also change with time during one landslide (e.g. fluid content or grain size distribution of the particles). As an example of the natural complexity, during landslide initiation stage, the mass generally consists of a cohesive rock that destabilizes due to different processes such as the increase of pore pressure due to rain falls or changes in the surrounding stresses due to seismic or volcanic events or groundwater table variations, erosion processes, aging, etc (e.g. Chigira et al 2003, Dunning et al 2007, Helmstetter and Garambois 2010, Kean et al 2013, Cappa et al 2014). Once destabilized the mass gets progressively fragmented in a way that is essentially not understood nor modelled experimentally or numerically (e.g. De Blasio and Crosta 2014). Finally, the respective contribution of the physical processes (e.g. pore-fluid viscosity, capillarity forces, etc) may be very different at small and large scale (e.g. Iverson 2015). The use of appropriate parameters for each scale makes it however possible for the same model to describe both laboratory and field flows (i.e. the terms in the equations will be the same but their relative magnitude will be different).

Despite this natural complexity, 'simple' laboratory experiments can be used to test if precise field observation (deposit extent, morphological features such as levees or fronts) and associated scaling laws may be qualitatively or quantitatively reproduced with or without specific physical processes (e.g. presence of fluid, particle segregation). Combining these experiments with analytical and numerical models makes it possible to understand the origin of these field data and scaling laws and to relate the parameters involved to the properties of the flowing material. Experimental and numerical models may also highlight specific behaviours, suggesting to look for specific features in the field or to revisit field data interpretation as illustrated in the following.

4.2.1. Effective friction and mobility of small to large landslides.

The initial aspect ratio of the granular mass released to create natural landslides is generally smaller or of the order of one (R  ⩽  1) (e.g. Legros 2002, Lucas et al 2014). Experimental dam-break experiments and their analytical and numerical modelling described in section 3.4 show that the scaling law ${{x}^{\ast}}=\left({{x}_{\text{f}}}-{{x}_{0}}\right)/{{x}_{0}}=\lambda R$ is observed on horizontal to moderate slopes (smaller than about half the repose angle of the material involved) and that the parameter λ could be approximated by k/(tan δ–tan θ), where δ is the effective friction angle, θ the slope angle and k a constant (equation (3.5)). For higher slopes, this parameter depends on the volume involved (Farin et al 2014). This scaling law is recovered by numerical thin-layer depth-averaged models for different shapes of the initial released mass (figure 10(C) and Lucas et al 2011). The effective friction experienced by the grains μ_eff  =  tan δ_eff can be deduced from (3.5) as a function of the slope and geometrical characteristics of the initial mass and deposit:

$$\begin{eqnarray}&&{{\mu}_{\text{eff}}}=\tan \theta +\frac{{{h}_{0}}}{{{x}_{\text{f}}}-{{x}_{0}}}\end{eqnarray} \tag{ 4.1 }$$

Zoom In Zoom Out Reset image size

As shown in figure 10(B), λ increases for flows with interstitial pore fluid pressure (obtained by fluidizing fine particles in laboratory experiments), reflecting the lower effective friction experienced by the grains (Roche et al 2011). Using (4.1) gives values of δ_eff  =  25°–28° for flows of coarse particles that cannot retain pore pressure, which is of the order of the repose angle of the particles involved, and δ_eff  =  14°–17° for flows of fine particles in which pore pressure decreases slowly (figure 10(B) and Roche et al 2011).

When trying to recover this scaling law for a series of about 40 well constrained landslides of different sizes in various environments and flowing over different slopes on Earth and on other planets, Lucas et al (2014) show that there is a very large dispersion of the data (figure 10(A)). Using (4.1) with the real values of h₀, x_f, x₀ and θ for each landslide leads to a wide range of values for μ_eff. This may explain the dispersion of the data shown in figure 10(A) as detailed below. On the contrary, when cxonsidering only very large landslides that occurred in the same environment (Valles Marineris on Mars) over gentle slopes of about 2°, the data obey to the scaling law given by (3.4) (figure 10(B)). The difference between the two data sets of Valles Marineris landslides represented in figure 10(B) highlights the difficulty of obtaining field measurements. Indeed, the data are different because Lucas et al (2014) only considered landslides that were not stopped by the opposite valley wall and used more precise satellite data than those available up to 2006. Furthermore, each point in Lucas et al (2014) represents a landslide while the points of Lajeunesse et al (2006) result from the averaging of several landslides. However, whatever the data set, the coefficient λ of the scaling law is much smaller than that of lab-scale dry granular flows. This suggests a much higher mobility of large Martian landslides than experimental granular flows in contrary to the conclusion of Lajeunesse et al (2006). Simulating these landslides based on precise Digital Elevation Model deduced from satellite data with a thin-layer depth-averaged model with only one free parameter, the effective friction coefficient, leads to very small friction coefficients between 0.1 and 0.16 (friction angles of about 6°  ⩽  δ  ⩽  9°) (Lucas et al 2014, table 1 of the Supplementary), which are similar to those in experiments on granular flows with reduced friction caused by pore fluid pressure (Roche et al 2011). Interestingly, the friction angles making it possible to reproduce numerically the extension of these Martian landslides occurring in the same environment and with a similar range of volumes (10¹⁰ m³  ⩽  V  ⩽  2.5  ×  10¹² m³) are almost the same, i.e. varying by less than 3° (Lucas and Mangeney 2007, Lucas et al 2011).

Calculating the effective friction from field data using (4.1) or using numerical simulation based on thin-layer depth-averaged model shows a clear decrease of the effective friction with increasing volume of the flowing mass (figure 11 and Lucas et al 2014). A distinct trend is observed despite the very different nature and location of the investigated landslides that occur in wet or dry environments on planets with or without atmosphere. Most of the numerical analyses of individual landslides in the literature lead to calibrated friction coefficient in agreement with figure 11 (e.g. Hungr and Evans 1996, Kelfoun and Druitt 2005, Pirulli and Mangeney 2008). The numerical models can be viewed here as an empirical tool making it possible to recover the effective friction while taking into account the geometry of the released mass, the dynamics of the flow (spreading, sliding, etc) and the topography effects. It represents a much better estimate of the effective friction than the H/L Heim's ratio widely discussed in the literature (e.g. Legros 2002, Lucas et al 2014 and references therein). Even though the origins of this friction weakening with volume (or with other related properties such as velocity, shear rate, granular temperature, fluids, flash-heating, etc) are still open issues, the proposed physical processes (e.g. Campbell et al 1995, Brodu et al 2015, Johnson et al 2016) should be quantitatively compatible with these field observations.

Zoom In Zoom Out Reset image size

4.2.2. Physical processes involved.

Going further than just comparing experimental or numerical deposit extent with field observation is very challenging. In particular, all the results described in the previous section are based on comparison between characteristics of experimental, simulated and observed deposits that represent only the final flow stage. However, deposits can be reproduced with models that do not accurately describe the flow dynamics. The effect of friction weakening laws or of other physical processes such as pore fluid pressure or substrate erosion that may significantly increase the flow runout or generate surges should be constrained with 'dynamic' flow data such as stress and pressure gages or seismic data as discussed in section 4.1. Figure 12 illustrates how seismic records in different frequency ranges can be used to extract information on the interaction between the flow and the topography. This interaction is strongly related to the flow rheology and to the physical processes involved (presence of fluid, entrainment, etc). As an example, combined seismic data analysis and numerical modelling of 200 rockfalls and pyroclastic flows in the same volcanic environment showed that taking into account friction weakening with increasing volume is necessary to reproduce seismic observations (Levy et al 2015). Recent laboratory experiments on acoustic waves generated by granular column collapses also show how the seismic energy reflects the different flow regimes. Such experiments permit to quantify the balance between loss of potential energy, frictional work and seismic energy (Farin 2015, Farin et al 2015, 2016). Another interesting example is the dynamic measurement of the position of the static-flowing interface within an erodible bed for snow avalanches (Sovilla et al 2006, 2010). The temporal variation of this interface is closely related to the rheology of the flow and substrate material involved (e.g. Lusso et al 2016a). Quantitative comparison with such measurements would make it possible to identify the difference between erosion rate with and without the presence of interstitial fluids. Indeed, laboratory experiments on dry granular flows showed that even without fluids, the presence of an erodible bed can significantly affect the flow dynamics and deposit (e.g. Mangeney et al 2010, Farin et al 2014).

Zoom In Zoom Out Reset image size

As described in section 4.1, a lot of detailed data on natural deposits (morphometric characteristics, size, shape and nature of the grains, etc) are commonly available. Let us discuss below a few examples of the use of laboratory experiments to interpret these field data and highlight the related issues. As described in section 3, laboratory experiments and their simulation have shown that neither grains polydispersity nor fluid effects are necessary to create self-channeling flows and levees-channel deposit, while these processes change the detailed shape of the flow and deposit. This has strong implication for planetary science and in particular on Mars (e.g. Mangold et al 2003, 2010) where the levee-channel landslide morphologies were first interpreted as indicating the presence of water during emplacement (Malin and Edgett 2000). Application to natural pyroclastic or submarine flows of the scaling laws described in section 3 relating the morphometric characteristics of the levee-channel deposit to the flow velocity and rheological properties may help to constrain the flow dynamics during emplacement (e.g. Jessop et al 2012, Cannat et al 2013)). Based on this approach, Jessop et al (2012) suggest that the quasi-static lateral borders of the Lascar pyroclastic flows behaved like dry granular material whereas within the central channel the flow had a high pore fluid pressure. The physical processes (particle segregation, pore fluid pressure, etc) that may lead to these differential behaviours are well explained by laboratory experiments and by their simulations as described in section 3. There is however a strong scattering in the scaling laws derived from natural data because the slope of the topography changes from more than 20° to less than 10° along the path of each pyroclastic flow, the flux that generated these flows is expected to change in time, and erosion of the substrate made of former deposits probably occurred. Even in the same geological setting, it is difficult to compare the detailed morphological features or composition of a flow with another to infer general flow behaviour because flows generally travel on different paths and were probably generated by different upslope conditions (e.g. Brand et al 2014). Other morphological features could be interpreted using laboratory experiments and/or numerical modelling. As an example, the large scale abrupt change in deposit thickness observed in the Socompa debris avalanche deposit in Chile (Kelfoun and Druitt 2005) or on the Skollahvilft snow avalanche in Iceland (Cui et al 2007) were shown to be related to the deflection of the flow by the local topography that generated a wave that froze when the avalanche stopped. These studies improve the understanding and prediction of maximum runup height on topographies or dams, with strong implication for hazard assessment and design of deflecting dams.

An open issue is the origin of the hummocky surface of large debris avalanches (figure 1(c)). The question as to whether they may result from flow instability or waves such as those described in section 3.5 has not been investigated in detail yet. Laboratory experiments, in particular, can provide fundamental insights into flow-substrate interaction processes. Experiments by Roche et al (2013b) have reproduced the counterintuitive modes of interaction between pyroclastic flows and horizontal granular substrates. Field observations show that pyroclastic flows whose matrix consists of ash of size less than 2 mm are commonly unable to rework unconsolidated substrates of fine particles of similar grain size whereas they can entrain large blocks commonly larger than 10 cm in diameter (Roche et al 2016). The experiments reveal that a flow actually slides on a finely grained substrate. In contrast in case of a substrate of coarse grains the flow particles first penetrate into the substrate interstices and basal shear promotes extraction of the grains. Substrate-derived particles are then first dragged slowly at flow base before being uplifted because of a dynamic pressure gradient at the flow-substrate interface that scales with the square of the flow velocity. In case a granular substrate is inclined, experiments by Rowley (2011) and Farin et al (2014) show that shear imposed by the flow on the underlying material generates features that resemble Kelvin–Helmholtz instabilities typical of pure fluids. The material entrained from the substrate can then be mixed thoroughly with the flow particles and may represent a significant proportion of the basal part of the deposits.