Friction law and hysteresis in granular materials

TL;DR

This paper investigates the velocity-weakening friction law and hysteresis in granular materials, revealing the role of acoustic noise and soft spots, and provides a microscopic theory with experimentally testable predictions.

Contribution

It introduces a microscopic theory explaining the non-monotonic friction law and hysteresis in granular materials, highlighting the role of acoustic noise and soft spots.

Findings

Velocity-weakening friction law exists even with identical microscopic friction coefficients.

Soft particles eliminate the velocity-weakening behavior.

Predicted scaling laws for acoustic noise and sliding contacts are confirmed numerically.

Abstract

The macroscopic friction of particulate materials often weakens as the flow rate is increased, leading to potentially disastrous intermittent phenomena including earthquakes and landslides. We theoretically and numerically study this phenomenon in simple granular materials. We show that velocity-weakening, corresponding to a non-monotonic behavior in the friction law $\mu(I)$, is present even if the dynamic and static microscopic friction coefficients are identical, but disappears for softer particles. We argue that this instability is induced by endogenous acoustic noise, which tends to make contacts slide, leading to faster flow and increased noise. We show that soft spots, or excitable regions in the materials, correspond to rolling contacts that are about to slide, whose density is described by a nontrivial exponent $\theta_s$. We build a microscopic theory for the non-monotonicity of $\mu(I)$, which also predicts the scaling behavior of acoustic noise, the fraction of sliding contacts $\chi$ and the sliding velocity, in terms of $\theta_s$. Surprisingly, these quantities have no limit when particles become infinitely hard, as confirmed numerically. Our analysis rationalizes previously unexplained observations and makes new experimentally testable predictions.

Figures

c in the main text. Defining $\Delta\mu_{hyst}\equiv\mu(I=2\times 10^{-5})-\mu(I_{*})$, we show directly the depth of the minimum for all data in Fig. [A.1](#S6.F1)a. Error bars are determined as follows: splitting the total strain into halves, for each run we compute $\Delta\mu_{hyst}$ in each half of strain. The error bar is the standard deviation of these values.

Equations

$$\mu (I)\approx {\mu }_c+c_1{I}^{{\alpha }_{\mu }}.$$

$$\mu (I)\approx {\tilde{\mu }}_c+c_{1}I={\mu }_c+{\mu }_c(g(\chi )-1)+c_{1}I.$$

$$\mathcal{R}\equiv \frac{E_{noise}}{E_p}\propto \frac{{\mathcal{D}}_{coll}}{{\Delta }^2}\frac{\tau }{k{D}^2}.$$

$$\mathcal{R}\sim Q\frac{{\mathcal{L}}^2{I}^3}{{\Delta }^{1/2}{ϵ}_v}\sim Q{\mathcal{L}}^2{\left(\frac{I}{{\Delta }^{1/4}}\right)}^2,$$

$${\mathcal{L}}_T\propto 1/\chi .$$

$$\delta \tilde{f}=p{D}^{d-1}{\mathcal{R}}^{1/2}.$$

$$\chi -{\chi }_c\approx \frac{1}{2}{\int }_0^{\delta f}dxP(x),$$

$$\displaystyle\chi\propto\mbox{$\frac{1}{2}$}{\cal R}^{\frac{1+\theta_{s}}{2}}.$$

$$\chi$$

$${\mathcal{L}}_T$$

$$\mathcal{R}$$

$$\delta f$$

$${\mathcal{R}}_{tot}\sim Q{\mathcal{R}}_{tot}^{-1-{\theta }_s}{\left(\frac{I}{{\Delta }^{1/4}}\right)}^2+\frac{{\mathcal{V}}^2}{{\Delta }^2}.$$

$$\chi =\frac{\sqrt{a_1+a_2(I/{\Delta }^{1/4}{)}^{2a_3}}}{1+a_4(I/{\Delta }^{1/4}{)}^{a_3}}$$

$${\chi }^2(\chi -{\chi }_c{)}^{\frac{2}{1+{\theta }_s}}\propto Q{\left(\frac{I}{{\Delta }^{1/4}}\right)}^2$$

Full text

\verticaladjustment

-2pt

\dropcap

A central property of particulate materials is their macroscopic friction (1, 2). In velocity-strengthening materials, friction grows with the deformation rate and flow is stable. By contrast, velocity-weakening materials are susceptible to instabilities under loading, including stick-slip. This classification is believed to distinguish faults that are prone or not to earthquakes (3), which ultimately corresponds to the shear of a granular material made of rocks and debris – the gouge – contained within the narrow fault between two tectonic plates. The Rice-Ruina rate-and-state model can describe both kinds of frictional behavior (4). It is, however, a heuristic model, and building a microscopic theory to justify a priori which materials weaken or strengthen under flow remains a challenge. To make progress, it is natural to consider well-controlled granular materials such as glass beads or sand, which have received considerable attention in the last decades. In these systems an important result comes from dimensional analysis: assuming that grains are strictly hard and that this limit is not singular, the macroscopic friction or stress ratio $\mu\equiv\sigma/p$, where $\sigma$ is shear stress and $p$ is pressure, can be shown to depend on strain rate $\dot{\epsilon}$ and pressure only via the dimensionless inertial number $I\equiv\dot{\epsilon}D\sqrt{\rho/p}$. Here $D$ is mean grain diameter and $\rho$ is grain density (5, 6, 7). The constitutive relation $\mu(I)$ has been extensively studied in the range $I\geq 10^{-3}$ and is found to be a growing function of $I$, corresponding to a velocity-strengthening material. However, several recent experiments have reported non-monotonic behavior of $\mu(I)$ for very small $I$(8, 9, 10), corresponding to velocity weakening. There is currently no microscopic theory rationalizing these observations.

As depicted in Fig.1a, the non-monotonic behavior of $\mu(I)$ must lead to hysteresis effects when $\mu$ is repeatedly cycled around its quasi-static value. Such hysteresis is well-known to characterize the jamming transition of granular materials. Indeed, the maximum angle of repose of a granular layer $\theta_{start}$ (corresponding to a macroscopic friction $\mu_{start}=\tan(\theta_{start})$) is larger than the angle $\theta_{stop}$ where avalanches stop (corresponding to $\mu_{stop}=\tan(\theta_{stop})$) (11). (Note that $\theta_{start}$ is measured by increasing the angle from a configuration that had stopped flowing at $\theta_{stop}$, otherwise it is not well-defined and depends on system preparation.) Several properties of this hysteresis should be explained by a microscopic theory: (i) externally applied vibrations (8, 12, 13, 10) can eliminate hysteresis, if the vibration velocity passes a threshold amplitude (13). (ii) Hysteresis appears to become very small if inertial effects are negligible (14). (iii) There is no evidence of velocity weakening for frictionless particles (15), as confirmed below 111This statement does not contradict the observation that the macroscopic friction at which yielding occurs depends on preparation (16). However, true velocity weakening in steady state without friction has only been observed for discontinuous, unphysical dissipation mechanisms (17)..

In this letter we build a microscopic theory that explains these observations, and justify why for hard enough frictional particles $\mu(I)$ is non-monotonic at small inertial numbers. Our results, which we test using the discrete element method, hold even if the friction between two particles is assumed to be a simple Coulomb law with identical static and dynamic friction coefficients. This demonstrates that hysteresis emerges as a collective effect, even when absent at the contact level. Our approach is based on two previous fundamental insights. First, acoustic emissions generated by particle collisions during flow can fluidize the granular material, reducing dissipation. This idea was proposed by Melosh (18, 19, 20) to explain why the macroscopic friction $\mu$ can be reduced by a factor of order $10$ from its static value during earthquakes. However, the quantitative treatment of this effect has been criticized (21), and the role of self-fluidization in regular granular materials is unclear, although self-fluidization was experimentally shown to affect the density of the flowing material (22).

To quantify this role, we will use a second, recent insight: elementary excitations are a key feature of amorphous solids, in particular their density affects a material’s stability and its susceptibility to plastic flow (23, 24, 25). These excitations correspond to soft spots – equivalent to dislocations in metals – where local rearrangements can be triggered if the stress is increased, or if mechanical noise is present. For hard frictionless particles, excitations correspond to contacts carrying a small force, which are thus susceptible to open (26, 27, 23). Here we will argue that for frictional granular materials, another kind of excitation appears, associated with contacts that are very close to sliding. Denoting by $\delta f$ the distance to sliding of a contact, we find that the distribution of excitations follows $P(\delta f)\sim\delta f^{\theta_{s}}$ with $\theta_{s}\approx-1/3$. From this knowledge, we can compute the effect of mechanical noise on the structure, which in turn affects the noise itself. Ultimately this leads to an experimentally testable microscopic theory of velocity weakening and mechanical noise. In addition, several non-trivial exponents are predicted that characterize noise and other microscopic quantities, which we find to be in good agreement with our numerical results.

1 Numerical evidence for non-monotonic flow curves

The macroscopic friction $\mu(I)$ has been studied extensively for inertial numbers $I\geq 10^{-3}$ (5, 6, 7, 28). In the range $0.1\geq I\geq 10^{-3}$ it follows

[TABLE]

This velocity-strengthening behavior is caused by an increase in the dissipation induced by collisions as $I$ increases (29). There is, however, no microscopic theory for the value of $\alpha_{\mu}$, except for frictionless particles (23). For frictional particles one observes $\alpha_{\mu}\approx 0.85$ (30), often approximated as $\alpha_{\mu}=1$.

Here we seek to numerically study if $\mu(I)$ becomes non-monotonic in the quasi-static regime $I\leq 10^{-3}$. However, if velocity-weakening is indeed present, then homogeneous flow is unstable, and shear bands or chaotic behavior are expected. Experimental reports of velocity-weakening use a set-up with stress gradients (8, 9, 10), which has been argued to stabilize flow, but makes quantitative analysis challenging. Here instead we use the fact that shear bands have a minimal thickness (typically of order ten grain diameters wide), and that flow can be stabilized by considering small systems. Then, if $\mu$ has an instability, the non-monotonicity will be most pronounced for small $N$, and decay as larger systems are considered.

In practice, our numerical results are created with the standard Discrete Element Method. We work in spatial dimension $d=2$ (our theoretical conclusions do not depend on $d$). The frictional disks follow static Coulomb friction: at each contact, the transverse component $f_{T}$ and normal component $f_{N}$ must satisfy $|f_{T}|\leq\mu_{p}f_{N}$, where we take the Coulomb friction coefficient to be $\mu_{p}=0.3$. This coefficient is the same for static and dynamic motion. Systems are simply sheared by horizontal motion of walls at the top and bottom domain edges, on a periodic domain. The confining pressure is controlled by fixed vertical forces exerted upon the walls. Particles have linear elastic-dashpot interactions, modeling a finite restitution coefficient $e<1$ and a particle stiffness $k$. $\Delta=p/k$, which characterizes the typical contact deflection relative to the grain diameter, is varied from $\Delta\approx 10^{-2}$ to $\Delta\approx 10^{-5}$.

Our results for $\mu(I)$ are shown in Figs. 1b,1c. A crucial finding is that non-monotonicity is indeed observed, but disappears if particles are too deformed (large $\Delta$). For $\Delta\approx 10^{-4},10^{-3}$ we find that the minimum occurs at $I_{*}\approx 10^{-3}$, the same inertial number below which intermittency and large fluctuations were previously observed in simulations (6). We measure the amplitude of non-monotonicity as $\Delta\mu_{hyst}\equiv\mu(I_{0})-\mu(I_{*})$ where $I_{0}=2\times 10^{-5}$ is the smallest inertial number we probe. As expected, $\Delta\mu_{hyst}$ increases as smaller systems are considered, but $I_{*}$ does not vary significantly, as shown in Fig.1c and SI.

By contrast, for frictionless particles velocity weakening is essentially absent, or at least much weaker than for frictional particles, as shown in Fig.A.1,A.2 of SI.

The existence of a minimum at any finite $I_{*}$ implies a dramatic behaviour as the shear stress is increased from the solid phase, as illustrated in Fig. 1a, in particular flow is predicted to start at a finite inertial number $I\approx 2\times 10^{-3}$. Essential questions include what microscopic mechanism leads to a minimum of $\mu(I)$ at $I_{*}$, and why does this instability disappear when $\Delta$ increases?

2 Mechanism for instability in granular media

We argue that these questions are naturally explained if one considers the role of the acoustic noise endogenously generated in flow, as measured for example in (22). In dense flows of hard particles, a network of contacts is formed that constrains motion: particles cannot overlap and cannot slide if the considered contact satisfies the Coulomb criterion. For infinitely hard particles, the dynamics only occurs along floppy modes for which these constraints are satisfied (31). For any finite $\Delta$, there will exist small-amplitude motion orthogonal to the floppy modes, which corresponds to vibrations of the contact network. The formation of new contacts through collisions pumps energy into these vibrations, which eventually decays due to grain visco-elasticity. Henceforth we call this vibrational energy per particle “mechanical noise” and denote it $E_{noise}$. If this noise is not negligible with respect to the characteristic potential energy in the contact, $E_{p}\propto\Delta^{2}$, it will affect the contact network and “lubricate” it (22, 32). In turn, this lubrication will facilitate particle motion, leading to faster flow, stronger collisions and larger noise. Our contention is that this positive feedback is responsible for instability in the flow curve and its associated hysteresis. This view naturally explains (i) why increasing particle deformability diminishes hysteresis, since it increases the potential energy which makes contacts less sensitive to noise, (ii) why immersing the grains in a viscous fluid reduces hysteresis, since it damps vibrations faster and (iii) why externally applied vibrations can eliminate hysteresis, since if the external noise is much larger than the endogenous noise, the positive feedback becomes negligible.

We now argue more specifically that the dominant “lubricating” effect of mechanical noise is to intermittently remove rolling constraints, thus increasing the fraction of sliding contacts $\chi$. In previous work, such an increase was indeed observed in response to a transient forcing (33). It is important to note that, with respect to noise, rolling constraints and no-overlap constraints are qualitatively different. Noise can lead to fluctuations in the distance between two particles in contact, but will never allow them to eventually penetrate each other. By contrast, vibration leading to intermittent stick-slip motion at a contact yields net relative motion between particles. This is illustrated by an object sitting on an inclined plane, which will eventually slide down if the plane is vibrated sufficiently. In our approach, this distinction explains why frictionless particles do not appear to present noticeable hysteresis, as shown in SI (15). (For very soft and elastic particles, other mechanisms than the removal of constraints described here could lead to a finite amount of hysteresis 777For highly deformable particles flow occurs via a succession of saddle-node bifurcations in the energy landscape, called shear transformations. Mechanical noise can then speed up flow by allowing jumps over barriers of potential energy, before these are destroyed by shear. This mechanism can cause hysteresis for very weak damping $e\approx 1$, see (34) for a recent discussion. In our framework such effects can be accounted for by replacing $\chi$ by the relative mechanical energy ${\cal R}$ in Eq.2. However for realistic restitution coefficient and particle hardness this effect appears negligible, since it should exist in frictionless systems in which velocity weakening is essentially absent.).

If $\chi$ is indeed controlled by mechanical noise, then it should be a function of both $I$ and $\Delta$. This prediction is confirmed in Fig. 2, which shows that $\chi$ is a function of $I/\Delta^{1/4}$. Intriguingly, this result implies that the hard sphere limit $\Delta\rightarrow 0$ is singular, at least for some microscopic properties. In the next section we will build a quantitative theory of mechanical noise and explain this scaling property. First, we argue why this behavior of $\chi$ can indeed generate hysteresis.

Consider the following Gedankenexperiment: the stress is increased in a granular solid. According to Eq.1, flow will start when $\mu$ reaches $\mu_{c}$. Now, consider the same protocol, but with vibrations imposed on the sample. Because of this noise, some contacts that were close to the Coulomb cone, as pictured in Fig.3a, will intermittently slide, leading to an overall increase of the fraction of sliding contacts, $\chi$. Being less constrained, the material will be less stable and start flowing earlier, for some $\tilde{\mu_{c}}<\mu_{c}$. We write $\tilde{\mu}_{c}=\mu_{c}g(\chi)$, where $g$ is a decreasing function and $g(0)=1$. The same reduction of macroscopic static friction must apply when noise is endogenously generated in flow. This effect can be included in Eq.1, which now becomes (taking $\alpha_{\mu}=1$ for simplicity):

[TABLE]

In this equation the last term is the usual velocity-strengthening effect due to the increased dissipation induced by collisions. The second-to-last term instead is velocity-weakening. To illustrate this point we consider a linear model $g(\chi)=1-b\chi$. Qualitative results do not depend on the parameter $b$; here we choose $b=0.3$. Other parameters entering Eq.2 are fixed by previous observations to $c_{1}\approx 1.4,\mu_{c}\approx 0.2$. Using fits to $\chi$ from Fig.2 (described in SI), Eq.2 predicts the curves shown in Fig. 1d. We see that the qualitative features of $\mu(I)$ are captured: $\mu$ has a minimum, whose depth, but not its location, depends strongly on $\Delta$.

3 Microscopic description of mechanical noise

Our approach to compute the mechanical noise and its effect on $\chi$ is illustrated in Fig. 3b. We proceed in three steps. First, we express the mechanical noise in terms of two kinetic properties: the dimensionless velocity fluctuations ${\cal L}=V_{r}/(\dot{\epsilon}D)$ (where $V_{r}$ is the typical relative velocity between neighboring particles, $D$ the particle size and $\dot{\epsilon}$ the strain rate) and $\epsilon_{v}$, the characteristic strain at which particles collide and change direction. Second, we use a previous argument based on energy balance relating $\chi$ to the normalized typical sliding velocity ${\cal L}_{T}$. Finally, we compute the density of sliding contacts $\chi$ induced by mechanical noise. This requires knowing the density of rolling contacts at a small distance $\delta f$ to sliding. From these three coupled equations quantities of interest are derived.

Along the way we make two approximations. First, we assume that the characteristic velocity ${\cal L}$ with which particles collide is identical to the characteristic sliding velocity ${\cal L}_{T}$. Assuming that there is only a single velocity scale in flow is essentially a mean-field approximation, known to be correct for frictionless particles (23) but not exact for frictional ones (29). Second, we assume that the vibrational noise is damped at a rate $1/\tau$ independently from the distance to jamming, which allows us to estimate the scaling of $E_{noise}$ directly in terms of the energy dissipated in collisions. The utility of these approximations is supported by the good comparison between predictions and numerics.

3.1 Computing mechanical noise from kinetics

$E_{noise}$ characterizes the density of vibrational energy of the contact network. The power injected in vibrational modes is supplied by collisions at a rate ${\cal D}_{coll}$. It can be estimated as the product of the kinetic energy of a particle $\sim{\cal L}^{2}\dot{\epsilon}^{2}mD^{2}$ times the collisional rate $\sim\dot{\epsilon}/\epsilon_{v}$, leading to ${\cal D}_{coll}\sim{\cal L}^{2}\dot{\epsilon}^{3}mD^{2}/\epsilon_{v}$ (23, 29). Introducing a time scale $\tau$ at which this vibrational energy is dissipated into heat, we get the estimate $E_{noise}=\tau{\cal D}_{coll}$. We make the simplifying assumption that $\tau$ does not depend on $I$, characterizing the distance to jamming. Then dimensional analysis implies that $\tau=C_{0}m/\eta_{N}$ where $\eta_{N}$ the damping coefficient of the particle interaction. In our approximation $C_{0}$ is a constant, which we expect to be rather large 777In a network of interacting particles with linear dashpot interactions, the rate $1/\tau_{\omega}$ at which the vibrational energy in a mode of frequency $\omega$ decays follows $\tau_{\omega}\sim 1/\omega^{2}$. Thus the energy in low-frequency modes takes a much longer time to decay. It is well known that both amorphous solids (23) and contact networks in flow (31) display an abundance of low-frequency modes. Their frequency scale can be $\sim 10$ times smaller than the characteristic frequency of a single contact, suggesting the order of magnitude $C_{0}\sim 10^{2}$.. Normalizing $E_{noise}$ by the typical potential energy in a contact $E_{p}=\frac{1}{2}k\Delta^{2}D^{2}$, we define their ratio

[TABLE]

Neglecting pre-factors, we obtain the scaling behavior:

[TABLE]

where $Q\equiv\sqrt{mk/\eta_{N}^{2}}$ is the bare quality factor of the grains. We used the definition $I\equiv\dot{\epsilon}D\sqrt{\rho/p}$ as well as the previous result $\epsilon_{v}\propto I$ observed experimentally (35), numerically (23, 29) and justified theoretically in (23, 29).

3.2 Constraint from energy balance

For realistic microscopic friction coefficient, it is found in dense flows that most of the dissipation is induced by sliding at frictional contacts (29) (the collisional dissipation ${\cal D}_{coll}$ is sub-dominant, but is the only one contributing to the mechanical noise). On average, the power dissipated in steady flow must balance the power injected by the shear stress; by a straightforward estimation of the sliding dissipation rate in terms of the typical sliding velocity (29), this implies that 111The power injected is $\Omega\sigma\dot{\epsilon}$, where $\Omega$ is the system volume. At each sliding contact, the power dissipated is $f_{T}u_{T}$, where $f_{T}$ is the transverse force and $u_{T}\sim{\cal L}_{T}D\dot{\epsilon}$ the transverse velocity. Equating the total power dissipated $\sim N\chi f_{T}u_{T}$ to $\Omega\sigma\dot{\epsilon}$, we find $ND^{3}\mu p\dot{\epsilon}\sim N\chi\mu_{p}f_{N}{\cal L}_{T}\dot{\epsilon}D$, or ${\cal L}_{T}\chi\sim\mu/\mu_{p}$. For small $I$, $\mu$ tends to a constant, so that ${\cal L}_{T}\chi\sim 1$.:

[TABLE]

where ${\cal L}_{T}$ is the characteristic dimensionless sliding velocity of sliding contacts. In what follows we assume ${\cal L}_{T}\sim{\cal L}$, which is an approximation (29, 36).

3.3 Noise-induced sliding

For harmonic grains, the mechanical noise corresponds to a characteristic force scale

[TABLE]

$\delta{\tilde{f}}$ characterizes the fluctuations of forces at contacts. Such fluctuations will induce contacts near the edge of the Coulomb cone, shown in red in Fig.3, to slide intermittently. Thus, it is of crucial importance to determine the effect of noise on the density of contacts at small distances ${\tilde{x}}=f_{N}-|f_{T}|/\mu_{p}$ from the Coulomb cone. Fig.4 shows the distribution $P(x)$ of the dimensionless quantity $x={\tilde{x}}/pD^{d-1}$. As $I$ decreases, $P(x)$ develops a divergence $P(x)\sim x^{\theta_{s}}$ where $\theta_{s}\approx-1/3$. At any time, we expect that about half of the contacts within a dimensionless distance $\delta f=\delta{\tilde{f}}/pD^{d-1}$ of the Coulomb cone will be induced to slide, implying

[TABLE]

where $\chi_{c}$ corresponds to the density of sliding contacts in the system when flow stops and noise is absent. As shown in SI, our numerics support that $\chi_{c}\rightarrow 0$ in the thermodynamic limit, but is finite for finite $N$. In what follows we therefore consider $\chi_{c}=0$, but our theory can readily be extended to the case $\chi_{c}>0$ to capture some finite size effects, see S.I. Inserting the scaling behavior $P(x)\sim x^{\theta_{s}}$ in Eq.7 leads to

[TABLE]

3.4 Results

Combining Eqs.4,5,8 and keeping only on the dependence on $I$ and $\Delta$ leads to222The $Q$ dependence is readily obtained in the Eqs.9,10,11,12 by replacing $\Delta$ by $\Delta/Q^{2}$:

[TABLE]

where $\alpha=\frac{1+\theta_{s}}{2+\theta_{s}}\approx 2/5$, $\gamma=\frac{2}{2+\theta_{s}}\approx 1.2$ and $\beta=\frac{1}{2+\theta_{s}}\approx 0.6$. In Eq.12 the characteristic dimensionless fluctuation of forces is obtained from Eq.6, and follows $\delta f\sim{\cal R}^{1/2}$. These fluctuations are expected to smoothen the distribution of the distance to the Coulomb cone $P(x)$ for $x\ll\delta f$. $\delta f$ can thus be extracted numerically by locating where the power law $P(x)\sim x^{\theta_{s}}$ breaks down, as shown in Fig.4b.

Eq.9 is tested in Fig.2, Eq.11 in Fig.5, Eq.12 in Fig.4 and Eq.10 in Fig.6. Note that testing Eq.11 requires measurement of $E_{noise}$, which is difficult to do directly in the numerics, as it would require identification of all the vibrational modes of the contact network, at all instants of time. To test Eq.11 we measure a proxy for $E_{noise}$: the characteristic kinetic energy in the relative motion between particles $m\langle u_{N}^{2}\rangle/2$, where $u_{N}$ is the normal velocity at contact, as shown in Fig.5. It is a lower-bound on $E_{noise}$ 777 Indeed for low-frequency vibrational modes, the kinetic energy is larger than the relative kinetic energy between neighboring particles..

Our prediction that the correct scaling variable for these macroscopic quantities is $I/\Delta^{1/4}$ is in perfect agreement with data. The predicted exponent $\alpha$ of Eq.10 is in excellent agreement with observations. This is also true for Eq.9 if finite size effects are included in the theory (see SI). Our measurements for the other exponents give $\gamma=1.6$ (Fig.5b) and $\beta=0.8$ (Fig.6b), overall in good agreement with predictions, considering the approximations made.

These results correspond to a detailed microscopic description of mechanical noise in granular materials. It is mean-field in character, but already appears to capture the essential aspects of the phenomenon. Obviously, the contact mechanics of real granular media can be much more complicated than the harmonic contact description discussed here (4, 37). For example, real grains in three dimensions are Hertzian, with an elastic potential $V(h)\propto E(h/D)^{5/2}$, where $E$ is the grain Young’s modulus. For such grains, the typical dimensionless contact deflection $\Delta\propto(p/E)^{2/3}$ and the characteristic stiffness scales as $k\propto ED(p/E)^{1/6}$. This entails that $\delta f\propto E_{noise}^{3/5}$. Likewise, the damping coefficient $\eta_{N}$ may also depend on $\Delta$. Our framework can readily be extended to include these different contact properties, which can slightly alter scaling exponents, but will not change the physical picture.

Another important possibility to consider is that the dynamic friction coefficient is smaller than the static one. This will not affect the description of mechanical noise leading to Eqs.9,10,11,12. However, it will cause the macroscopic friction to decrease even more with $\chi$. In our phenomenological model, it corresponds to an increase in the coefficient $b$ appearing in Eq.2, which will enhance the amplitude of hysteresis.

4 External noise in granular flows

Applying external vibrations is known experimentally to affect rheological properties and hysteresis (8, 37, 38, 32, 13, 10), as we now quantify. Even small typical noise can generate very rare events where noise is locally large, leading to rearrangements and creep flow. We do not seek to describe this effect here, which occurs at very small inertial number (10). Instead we focus on the dynamical range relevant to hysteresis and to our theory. When a sheared material is subject to external noise, the dimensionless noise ${\cal R}$ will get an additional contribution. Let $E_{e}$ be the characteristic energy per grain from the external noise, and define ${\cal R}_{ext}\equiv E_{e}/E_{p}\sim({\cal V}/\Delta)^{2}$, where ${\cal V}=\sqrt{E_{e}/(\frac{1}{2}kD^{2})}$ and ${\cal R}_{tot}={\cal R}+{\cal R}_{ext}$. Using Eqs.4,5 for ${\cal R}$, and Eq.8 with ${\cal R}_{tot}$ in place of ${\cal R}$, we find

[TABLE]

Thus external noise is relevant for $I\lesssim I^{\dagger}$, where $I^{\dagger}\leavevmode\nobreak\ \sim\leavevmode\nobreak\ {\cal V}^{2+\theta_{s}}\leavevmode\nobreak\ \Delta^{-7/4-\theta_{s}}\leavevmode\nobreak\ Q^{-1/2}$. The results of Eq.13 are shown in Fig. 6b for a range of ${\cal V}$ (with $\theta_{s}=-1/3$). There are two regimes: for weak noise, $\mu$ remains non-monotonic, but the value of $\mu_{c}$ is lowered. However, as noise is increased and $I^{\dagger}>I_{*}$, the non-monotonicity of $\mu$ disappears entirely. We thus predict that for sufficiently strong applied noise, flow becomes stable and hysteresis disappears.

In the slider-block experiments of (13), a block is pulled over a granular bed in the presence of external vibrations of amplitude $A$ and frequency $\omega$. A transition from stick-slip to steady motion is observed, which does not depend on $A$ and $\omega$ independently, but only on their product $A\omega$. This is expected from our analysis, since the external vibrational noise controlling the transition should be of order $E_{e}\sim\frac{1}{2}m(A\omega)^{2}$, indeed a function of $A\omega$ only.

Our framework can also predict the noise needed to trigger flow in a static configuration with $\mu=\mu_{c}-\delta\mu$. Yielding will occur when the noise-induced reduction of static macroscopic friction is larger than $\delta\mu$, which according to Eqs.2,8 occurs when \delta\mu<b\mu_{c}\chi\sim b\mu_{c}\big{(}{\cal V}/\Delta\big{)}^{1+\theta_{s}} or ${\cal V}>C\Delta(\delta\mu)^{1/(1+\theta_{s})}$, where $C$ is a constant. This prediction is consistent with the earlier observation that the noise threshold increases with the confining pressure (39, 40). It would be interesting to further test the dependence of this threshold on $\delta\mu$.

5 Conclusion & Outlook

We have shown that velocity-weakening in dry granular flows is a collective phenomenon, which emerges even if absent at the contact level. We have explained this observation based on a microscopic theory characterizing the endogenous mechanical noise induced by collisions, as well as its effect on the structure. To understand the latter a characterization of elementary excitations on granular packings was performed, argued here to correspond to contacts that are rolling, but close to sliding. This framework rationalizes several experimental observations, including the factors governing the strength of hysteresis and the effects of exogenous noise. It also makes detailed scaling predictions on the noise level, the fraction of sliding contacts and the typical sliding velocity, in good agreement with numerical observations. Several new predictions are made that can be experimentally tested, most obviously the fact that velocity weakening should become very small or even vanish if sufficiently deformable particles are considered. This could be checked in an inclined plane experiment. From a theoretical perspective, an important question for the future is what determines the exponent $\theta_{s}$ characterizing excitations in frictional packings – a question now well-understood for hard frictionless particles (26, 41).

It is interesting to compare our results to Melosh’s work on earthquakes (18, 19), where he proposed that acoustic emission fluidizes the fault gouge and causes the macroscopic friction to drop by a factor of order ten, as observed (42). Melosh assumes that most of the energy dissipated goes into vibrations, whereas we find numerically in granular flows that for inertial numbers of interest for velocity weakening (say $I=10^{-3}$), this is only true for a few percent of the injected power (29). Likewise, in granular materials the hysteresis is a small effect, 10% at most (1). Why such effects would be much larger in the context of a fault gouge is currently unclear, although it has been reported based on landslides that velocity weakening effects appear to grow with the length scale of the phenomenon considered (20).

Here we have focussed on $\mu$, but the particle volume fraction $\phi$ is also important in constitutive relations. Several works have observed that $\phi(I)$ can become non-monotonic, corresponding to anomalous compaction (22, 9, 43) 222In our numerical results we did not detect a noticeable non-monotonicity in $\phi(I)$. As indicated by (22), particle angularity may significantly enhance such anomalous compaction.. In (22) this compaction was argued to arise from mechanical noise, and shown to diminish under the action of external vibrations, in a scenario qualitatively similar to our theory. Future work should carefully probe both $\mu$ and $\phi$ together to see if this link can be strengthened.

Finally, it would be desirable to extend the present description to over-damped non-Brownian dense suspensions. Earlier experiments show that hysteresis increases with inertia, but do not rule out a finite amount of hysteresis in the viscous limit (14, 44). Over-damped simulations support that $\chi$ increases with the viscous number (36) (playing a role similar to the inertial number), but very weakly. Following the discussion of Eq.2, this behavior may be sufficient to lead to hysteresis. In that case, however, there is no microscopic theory for $\chi$ available.

\acknow

We thank M. Bouzid, E. Clément, Y. Forterre, X. Jia, B. Metzger, and V. Vidal for discussions. M.W. thanks the Swiss National Science Foundation for support under Grant No. 200021-165509 and the Simons Collaborative Grant “Cracking the Glass Problem" #454953.

\showacknow

\pnasbreak

6 Supplementary Material

6.1 Friction and system-size dependence of hysteresis magnitude

As described in the main text, when smaller systems are considered, the magnitude of hysteresis is expected to increase. This can be inferred from Fig.1c in the main text. Defining $\Delta\mu_{hyst}\equiv\mu(I=2\times 10^{-5})-\mu(I_{*})$, we show directly the depth of the minimum for all data in Fig. A.1a. Error bars are determined as follows: splitting the total strain into halves, for each run we compute $\Delta\mu_{hyst}$ in each half of strain. The error bar is the standard deviation of these values.

We verified that for frictionless systems, hysteresis is absent or significantly smaller than for frictional particles. Following the same protocol as for the frictional systems, we computed $\mu(I,\Delta,N)$ for $N=10^{2.2},10^{2.6},10^{3}$ and $\Delta\approx 10^{-4},10^{-3},10^{-2}$. In Fig. A.2 we show example curves of $\mu$ at fixed $\Delta=10^{-2.8}$ (Fig. A.2a) and at fixed $N=10^{2.6}$ (Fig. A.2b). The horizontal and vertical range is equal to curves Fig. 1b, Fig. 1c in the main text. A dip is hardly distinguishable, but can be determined algorithmically in some cases, as for the frictional data. The resulting hysteresis magnitude is shown in Fig. A.1b and is seen to be negligible compared to the frictional case considered in the main text.

6.2 Phenomenological Model Parameters

The parameters $\mu_{c}$ and $c_{1}$ in the phenomenological model can be estimated from existing data (45, 29). Fig.A.3a shows $\mu$ for a range of $\mu_{p}$. For realistic $\mu_{p}$ in the range $0.2<\mu_{p}<0.8$, we find $0.23<\mu<0.27$. In the main text we take $\mu_{c}=0.2$.

We also use the linear form of $\mu(I)$ at intermediate $I$. Experiments and simulations are often fit to a form $\delta\mu(I)=c_{1}I_{1}/(1+I_{1}/I)$, which behaves as $\delta\mu(I)\approx c_{1}I$ for $I$ in the dense regime, $I\ll I_{1}\approx 0.3$ (45). Here $c_{1}\approx 1.4$ independently of $\mu_{p}$ (45).

Finally, to show the predictions of Eq.2, we use a functional form of $\chi(I,\Delta)$ determined from data. As shown in Fig.A.3 for $\Delta=10^{-3.8}$ and various $N$, our data is very well fit by

[TABLE]

The constant $a_{1}$ vanishes as $N\to\infty$; we use this large $N$ version in Fig. 1d in the main text.

6.3 Finite-size effects in $\chi$

$\chi$ has a plateau value as $I\to 0$ that depends on $N$. Using the fit form described above, we find that all data can be collapsed onto a master curve with finite-size scaling, as shown in Fig. A.4a for $\Delta=10^{-3.8}$. The decay of $\chi_{c}=\chi(I\to 0)$ with $N$ is shown in Fig.A.4b.

If the existence of $\chi_{c}$ is taken into account in the theory, this changes Eq.(8) to $\chi-\chi_{c}\propto\frac{1}{2}{\cal R}^{\frac{1+\theta_{s}}{2}}$, giving

[TABLE]

Eq. A.15 is a closed equation for $\chi(I/\Delta^{1/4})$, and can be solved numerically. As in the main text, we predict that $\chi$ depends on $I/\Delta^{1/4}$, as verified in Fig.2. The effect of $\chi_{c}$ in Eq. A.15 is to saturate $\chi$ at small $I/\Delta^{1/4}$, as shown by the dash-dotted line in Fig.2.