Unifying Suspension and Granular flows near Jamming

Rheological properties of dense flows of hard particles are singular as one approaches the jamming threshold where flow ceases, both for granular flows dominated by inertia, and for over-damped suspensions. Concomitantly, the lengthscale characterizing velocity correlations appears to diverge at jamming. Here we review a theoretical framework that gives a scaling description of stationary flows of frictionless particles. Our analysis applies both to suspensions and inertial flows of hard particles. We report numerical results in support of the theory, and show the phase diagram that results when friction is added, delineating the regime of validity of the frictionless theory.


Introduction
Microscopic description of particulate materials such as grains, emulsions or suspensions is complicated by the presence of disorder, and by the fact that these systems are often out-of-equilibrium. One of the most vexing problems is how these materials transition between a flowing and a solid phase. When this transition is driven by temperature, it corresponds to the glass transition where a liquid becomes a glass, an amorphous structure that cannot flow on experimental time scales. Here we focus instead on athermal systems driven by an imposed stress, such as granular flows, and consider both the case where inertia is important (such as in aerial granular flows) or not (such as over-damped suspensions). We focus primarily on the case of hard particles.
Empirical constitutive relations have been proposed to describe such dense flows in the limit of hard particles [1][2][3]. Two important dimensionless quantities are the packing fraction φ and the stress ratio µ ≡ σ/p (also called the effective friction), where σ is the applied shear stress and p the pressure carried by the particles. For inertial flow, dimensional analysis implies that both quantities can only depend on the strain rate˙ , p, the particle diameter D and the mass density of the hard particles ρ via the inertial number I ≡˙ D ρ/p. One finds empirically that the constitutive relations µ(I) and φ(I) converge to a constant as I → 0, corresponding to the jamming transition where flow stops. We define µ(0) ≡ µ c and φ(0) ≡ φ c , which are system-specific and depend on particle shape, poly-dispersity, friction coefficient, etc. Near jamming, the constitutive relations are observed to be singular with: e-mail: eric.degiuli@epfl.ch e-mail: matthieu.wyart@epfl.ch As jamming is approached the dynamics becomes increasingly correlated in space [4,5]. By considering the dominant decay [6] of the velocity correlation function, one can define a length scale c : Similar dimensional arguments have been made for dense suspensions of non-Brownian particles [7,8]. In that case the relevant dimensionless number is the viscous number J = η 0˙ /p, where η 0 is the viscosity of the solvent. Empirically one finds similar relations: These relations imply that the viscosity η = σ/˙ of the suspension diverges as jamming is approached. Indeed Eq.(4) implies that σ ∼ p near jamming (in our scaling arguments below we may thus exchange freely σ and p), so that J ∝ η 0 /η. Eq.(5) then implies that: When both viscosity and inertia are present, a transition from viscous to inertial flow occurs as strain rate˙ is increased at fixed volume fraction [9][10][11]. This defines a cross-over strain ratė where the prefactors follow from a dimensional analysis. Empirical values found for the exponents in Eqs.(1-7) are reported in Table 1. They seem not to depend on dimension, which we thus did not report in our table. In the case of inertial flow they appear to depend on the presence of friction, whereas for suspended particles exponents appear to be similar with and without friction. In this work we focus initially on theory for frictionless particles, and in a second section, discuss the effect of friction.
Currently there is no accepted microscopic theory describing quantitatively these singular behaviors, in particular Eqs. (1)(2)(3)(4)(5)(6)(7). Observations support that as jamming is approached, particles form an extended network of contacts, and that the stress is dominated by contact forces [2,8,13]. In this work we review a framework to describe flow in such situations. The detailed arguments are presented in [23]; here we discuss the essential ingredients of the theory, and focus on the effect of friction, discussed in more detail in [18,24].
We attack the problem in two steps. First, we isolate the microscopic quantities that control flow. Then, we compute the scaling properties of these quantities by performing a perturbation around the solid phase. The idea is to consider the solid in the critical state, i.e. carrying the maximal stress ratio possible µ = µ c , corresponding to a packing fraction φ = φ c . Next, one adds an additional kick to the system, corresponding to a small additional stress ratio δµ. As a result, some contacts between particles will open, forces will be unbalanced, and the system will start to flow (see Fig. 1). Our key assumption is that flowing configurations are similar to a solid that is thus destabilized. This approach enables us to propose a full scaling description of the problem, and to predict the exponents entering Eqs. (1)(2)(3)(4)(5)(6)(7) in good agreement with observations in the absence of static friction. Moreover, our approach predicts several other properties singular near jamming: the speed of the particles, the strain scale of velocity decorrelation, and the coordination of the contact network. The first two quantities are accessible experimentally, and provide an additional experimental test of our views.
When an additional stress ratio is imposed to a marginal solid, naively one expects the contacts carrying the smallest forces to open. However, only a vanishingly small fraction of weak contacts are significantly coupled to external stresses [28]. We call them extended contacts, because perturbing such contacts mechanically lead to a  [28,29]. In a packing only those contacts lead to plasticity when stress is increased, or when a shock (say a collision) occurs in the bulk of the material [28,30]. The density of extended contacts as a function of the force f in the contact follows Numerically it is found that θ e ≈ 0.44 both in two and three dimensions [28,29], suggesting that this quantity may be independent of dimension. Figure 2. Extended vs localized contacts. When a contact is opened from an isostatic packing, the resulting deformation (arrows) can either be extended, as shown at left, or localized, as shown at right. Localized contacts are more numerous, but only extended contacts couple strongly to an imposed shear stress.
Reproduced from [28] by permission of The Royal Society of Chemistry (RSC).
Moreover, its value does not depend on the preparation protocol of the isostatic state: up to error bars, equal values are found from compression of hard spheres [28], shear-jammed hard disks [31], and decompression of soft spheres [29,32] 1 . The exponent θ e can be shown to control the stability of the solid phase [28,33,34]. Recently replica calculations in infinite dimension on the force distribution [32,35,36] led to the prediction [29]: within the error bar of our measurements. In our proposed scaling description all exponents can be expressed in terms of θ e , in particular: Empirically it was noticed that γ µ = γ φ and that α µ = α φ , which our arguments rationalize.

General Approach
We argue that several dimensionless quantities that characterize the microscopic dynamics under flow critically affect rheological properties. As jamming is approached, the assembly of particles acts as a lever: due to steric hindrance, the typical relative velocity between adjacent particles V r becomes much larger than the characteristic ve-locity˙ D where˙ is the strain rate and D the mean radius of the particles [37,38]. We thus define the amplitude of this lever effect L as: Another fundamental quantity, particularly relevant for inertial flow, is the strain scale v beyond which a particle loses memory of its direction relative to its neighbors.
v can be extracted from the decay of the autocorrelation where the average is made over all pairs of adjacent particles α. As the packing fraction φ increases toward jamming, collisions are more frequent per unit strain (due to the increase of relative particle motion L), and each collision affects the motion of the particles on a growing length scale. These two effects imply that v vanishes rapidly near jamming.
We now argue that dissipation is entirely governed by L in overdamped suspensions, and by both L and v in inertial flows. In both cases the power injected into the system at the boundaries, which is simply P = Ωσ˙ at constant volume, must be dissipated in the bulk.
In a dense suspension we expect dissipation to be governed by local mechanisms such as lubrication. Lubrication forces are singular for the ideal case of perfectly smooth spheres, but not for rough particles where they must be cut off. Thus the viscous force exchanged by two neighboring particles must dimensionally follow F ∼ η 0 V r D d−2 , leading to a power dissipated P/N = Cη 0 V 2 r D d−2 where d is the spatial dimension, C is a dimensionless constant that depends on the particle shape and roughness, and N is the number of particles. Equating the power dissipated to the power injected, one gets that for a given choice of particles: implying that the divergence of viscosity is governed by L. This result holds by construction in simple models of dissipation in suspension flows [37][38][39].
Concerning inertial flows, we suppose that the restitution coefficient characterizing a collision between two particles is smaller than one, and that collisions dominate dissipation. Then each time two neighboring particles change relative direction, a finite fraction of their relative kinetic energy E c ∼ MV 2 r must be dissipated, where M is the particle mass. Then the total power dissipated must follow P ∝ N˙ E c / v . Using Eq.(14) and balancing power injected and dissipated, one gets σ/(˙ 2 D 2 ρ) ∼ L 2 / v where ρ is the mass density of the particles, so that the inertial number I follows: To our knowledge Eq. (16) has not been proposed before, and could be tested empirically.

Perturbation around the solid
To obtain a complete description of flow, one must therefore express L and v in terms of control parameters such as δµ or δφ. To achieve this goal, we make the assumption that the contact network of configurations in flow is similar to that of jammed configurations at µ c immediately after increasing the stress ratio by δµ > 0. The detailed arguments in this section are presented in [23]; here we discuss the physical ingredients of each argument and summarize the results. The coordination z of the network of contacts is a key microscopic quantity that distinguishes flowing from jammed configurations. At jamming the coordination is just sufficient to forbid motion, corresponding to z c = 2d for frictionless spheres [40][41][42]. As illustrated in Fig.(1), the kick of amplitude δµ opens a fraction δz ≡ z c − z of the contacts, allowing collective motions of the particles for which particles do not overlap, but simply stay in contact, the so-called called floppy modes. Using a virtual work theorem to compute the work done after a contact is opened, we find [23] that the lever amplitude is directly related to the density of floppy modes δz, in particular To see how many contacts are opened when an increment of stress ratio δµ is added to an isostatic packing, we can use the known behavior of the shear modulus near jamming. This results in Jointly Eqs. (17,18) predict a relationship between lever amplitude and stress ratio: Eqs. (15,19) lead to a prediction for the exponent γ µ entering in the constitutive relation µ(J). Together with previous results showing that c ∼ 1/ √ δz [6,43], we obtain expressions for γ and α , corresponding to: both for inertial and viscous flows. We now consider the characteristic strain scale v at which velocities decorrelate. This is fixed by a geometrical argument: we use the fact that in dense flows, when a grain has an unbalanced net contact force, F, the ensuing motion will tend to make the remaining contacts of the grain align along F (see Fig.3). Since forces are repulsive, this further increases the unbalanced contact force. The increase in force is proportional to the typical contact force, pD d−1 , as well as to the rotation of the contacts, of magnitude ∼ Ld , thus where L is the dimensionless magnitude of the velocity fluctuation. This equation can also be derived formally, see [23]. In inertial flow, unbalanced forces are proportional to accelerations, F = pI 2 dL/d , which leads to Eq.22 indicates that there is a characteristic strain scale v ≈ I in which a velocity fluctuation grows by an amount proportional to its initial magnitude. In steady flow, such growth must be destroyed by collisions on the same strain scale, since the latter reorganize the direction of particle motion. Hence v is indeed the scale of decorrelation of particle velocities. Together with Eqs. (16,19) this result leads to a prediction for the exponent α µ characterizing the constitutive relation µ(I).
For viscous flows, we still have Eq. (21), but now unbalanced contact forces are equal to the viscous drag forces, since the total net force vanishes. In this case By a similar argument as above, this implies a characteristic strain scale v ≈ J in which velocity fluctuations are created and destroyed. One missing link to obtain a full scaling description of the problem is how the packing fraction depends on other control parameters. If we make the additional assumption that isotropic packings of frictionless particles in the thermodynamic limit have a finite (although presumably small) dilatancy, then this implies the scaling relation δφ ∼ δµ, known to agree well with observations [23]. This result enables us to predict the exponents α φ and γ φ entering the constitutive relation for φ(I) and φ(J), leading to a complete scaling description of rheological properties near φ Figure 3. Illustration of geometrical nonlinearity. If the central grain has an unbalanced contact force as indicated by the arrow, then the ensuing flow will tend to align the contact normals of the dominant contact forces (thick lines), i.e. the angle φ will increase. Geometrically, dφ/dt ∝ V, the velocity of the particle.  Figure 4. Phase diagram of dense homogeneous inertial frictional flow. In the frictionless and rolling regimes, most energy is dissipated by inelastic collisions, while in the frictional sliding regime energy dissipation is dominated by sliding. Along the phase boundary, grains dissipate equal amounts of energy in collisions and in sliding. For I 0.1, one enters the dilute regime [16]. The dashed line has slope 2.
jamming for frictionless particles. In particular the divergence of viscosity with packing fraction in suspensions is expected to follow: Finally, by a straightforward analysis one can extract the strain rate when one expects a transition from viscous to inertial flow, leading to Eq.(8) [23].
Our results are compared to previous empirical and numerical observations in Table 1. Overall we find a very good agreement between observations and predictions for frictionless particles. However, friction appears to change exponents for inertial flows, and may also affect viscous flows.

Phase Diagram
In the theory for frictionless particles, the fact that the power injected and the power dissipated exactly compensate on average plays a crucial role in relating the strain rate to geometrical quantities, as expressed in Eqs. 15,16. Central to the arguments are the dominant dissipation mechanisms, assumed above to be viscous dissipation in over-damped dynamics, and grain inelasticity in inertial dynamics. However, when particles are frictional, energy will also be dissipated by sliding at frictional contacts, in addition to the other dissipation mechanisms above. As Coulomb friction coefficient µ p and shear rate are varied, it is not obvious a priori in which regimes the sliding dissipation rate D slid will dominate over inelastic collisions in inertial flows, and viscous dissipation in viscous flows; we call these latter sources of dissipation D 0 .
This question can be straightforwardly investigated with numerical simulations. To map out the phase diagram as µ p and shear rate are varied, we considered, separately, the parameter space (µ p , I) for inertial dynamics [18], and (µ p , J) for viscous dynamics [24]. and inertial dynamics, we find 3 phases: (i) a 'Frictionless' regime in which D 0 dominates, most contacts are sliding, and investigated quantities are consistent with the theory for frictionless particles, in particular L, v , δµ, and δφ; (ii) a 'Rolling' regime in which again D 0 dominates, but most contacts are rolling, not sliding. Here L, δµ, and δφ appear to be consistent with the frictionless theory, particularly in the limit µ p 1. Finally, in between these regimes we find (iii) a 'Frictional Sliding' regime in which D slid > D 0 .
Here v ∼ I(J) in inertial (viscous) dynamics, as predicted by Eqs.21,22,23, which hold in the presence of friction, but other quantities have new scalings: for example, L has a weaker divergence, and δµ has a much less dramatic behavior at small strain rate.
The crossover from 'Frictionless' to 'Frictional Sliding' regimes can be predicted from the theory above. Indeed, D slid is simply the sum over sliding contacts of f α T · U T α , where T indicates the tangential component. Denoting by χ the fraction of contacts that are sliding, we have D slid ≈ N C χµ p f N L˙ , where N C is the total number of contacts and f N ∼ p is the typical normal force. Inside the frictionless regime, nearly all contacts are sliding, χ ≈ 1. Therefore we find that sliding friction constitutes a fraction of the total energy dissipation rate (which must equal P in steady flow). Since L is diverging as jamming is approached (Eq.19) we thus expect a transition to a regime dominated by frictional sliding, both for inertial and viscous flows, as observed. From Eqs. (15,16) and v ∼ I, it follows that L ∼ I −1/2 and L ∼ J −1/2 in frictionless inertial and viscous dynamics, respectively, thus leading to a crossover 'Frictionless' to 'Frictional Sliding' at I ∼ µ 2 p and J ∼ µ 2 p . These scalings are plotted as dashed lines in Figs.4,5, and are consistent with the data.

Frictional Theory
Still missing is a complete theory of the Frictional Sliding regime. We expect Eqs.21,22,23 to be always valid, and indeed the associated prediction v ∼ I appears to hold in the Frictional Sliding regime as well. However, the velocity fluctuations and rheological properties present new scalings. In particular, in all cases the velocity fluctuations are large and their spatial correlations are long-range as jamming is approached. For frictionless particles, a single characteristic scale characterizes these velocity fluctuations. By contrast, in [24] we provided evidence that several scales are required to describe the kinetics of frictional particles. For example, the relative velocity between contacting particles differs from the amplitude of non-affine velocities. This additional complexity is plausibly caused by the intermittent and non-extensive localization of strain observed in several studies [44][45][46][47]. A similar localization is observed in the plasticity of soft particles [48], whose connection to hard frictional particles is, however, unclear. A possible factor causing the difference between frictional and frictionless particles relates to the criticality of the marginal solid: for hard frictionless spheres, Maxwell counting leads to a precise identification of the marginal state, which is exactly isostatic. By contrast, frictional solids at jamming are generally hyperstatic [49][50][51]. How these issues are related are outstanding problems.

Conclusion
We reviewed a theory for flow of hard, frictionless particles, both in viscous and inertial dynamics. In a first step, we related the power dissipated in flow of frictionless particles to certain microscopic kinetic quantities. The latter control singularities in the rheological properties near jamming. In a second step, we have computed these quantities, using a perturbation around the solid phase. Our main hypothesis is that configurations in flow are similar to jammed configurations at maximum stress ratio µ c , destabilized by an additional stress increment δµ. In this approach, the properties of the solid phase are central, in particular the fact that the density of contacts which can couple to external forces is singular at small forces, and characterized by a non-trivial exponent θ e . Our description of flow can thus be thought as that of a jammed solid, populated by elementary excitations corresponding to the opening of weak contacts, of density δz.
When friction is added, the theory holds in a finite region of parameter space, whose boundary is well predicted by the theory. Theory for the regime in which sliding friction dominates is, however, lacking. In our view, a central question for the future is what controls the stability of isostatic frictional systems, how these respond to an additional stress ratio δµ, and how the combination of finite softness, inertia and friction qualitatively affects the flow curves [52].