Flow, Ordering and Jamming of Sheared Granular Suspensions
Denis S. Grebenkov, Massimo Pica Ciamarra, Mario Nicodemi, Antonio Coniglio
aa r X i v : . [ c ond - m a t . s o f t ] F e b Flow, Ordering and Jamming of Sheared Granular Suspensions
Denis S. Grebenkov,
1, 2, ∗ Massimo Pica Ciamarra,
3, 1
Mario Nicodemi,
1, 4 and Antonio Coniglio Dip.to di Scienze Fisiche, Universit´a di Napoli “Federico II” and INFN, Naples, ITALY LPMC, C.N.R.S. – Ecole Polytechnique, F-91128 Palaiseau, FRANCE CNISM & Dip.to of Information Engineering, Seconda Universit´a di Napoli, Aversa, ITALY Complexity Science & Department of Physics, University of Warwick, UK (Dated: Received: December 1, 2018/ Revised version:)We study the rheological properties of a granular suspension subject to constant shear stressby constant volume molecular dynamics simulations. We derive the system ‘flow diagram’ in thevolume fraction/stress plane ( φ, F ): at low φ the flow is disordered, with the viscosity obeying aBagnold-like scaling only at small F and diverging as the jamming point is approached; if the shearstress is strong enough, at higher φ an ordered flow regime is found, the order/disorder transitionbeing marked by a sharp drop of the viscosity. A broad jamming region is also observed where, inanalogy with the glassy region of thermal systems, slow dynamics followed by kinetic arrest occurswhen the ordering transition is prevented. PACS numbers: 45.70Vn,83.50.Ax,83.10.TvKeywords: granular suspension, kinetic arrest, jamming, order/disorder transition
Under the effects of external drives granular media ex-hibit a variety of complex dynamical behaviors. Com-paction under shaking is a well studied phenomenon,characterized by slow relaxation followed by a jammingtransition at high volume fraction [1, 2, 3] with deepanalogies to thermal glasses [2, 3, 4, 5]. In addition shear-induced transitions from flowing to jammed states [6] orto ordered flowing states are observed [7, 8]. In such acomplex panorama, the rheology of the different phasesof a system under shear is far from being clarified, andeven a clear ‘flow diagram’ locating the regions with or-dered flow, disordered flow, and jamming, is missing.Within such a perspective, we investigate by moleculardynamics (MD) simulations the rheology of a dense gran-ular suspension subject to a constant shear stress in a boxof constant volume. Alike previous experiments on non-Brownian particles (see [9, 10] and references therein),we opt for a constant shear stress rather that a constantshear rate setup (see [8] and references therein) becausejamming may be precluded in the latter case. We derivethe system ‘flow diagram’ as obtained by shearing for along time a disordered assembly of grains. It includes anorder/disorder transition line, signaled by a sharp drop ofthe effective viscosity η , as well as a jamming region. Sim-ilarly to supercooled liquids, a crossover from a transientdisordered flow to a stationary ordered one is also found,jamming occurring when ordering is kinetically impeded.In the regime of disordered flow, the system viscosity, η ,is found to grow with the applied shear stress accordingto Bagnold scaling, but it departs from it when the shearstress increases. As the packing fraction grows towardsa critical value, η is found to diverge with a power law.We also show that, similarly to glassy thermal systems,the long time state of the system in the jamming regionis dependent on the dynamical preparation protocol.We consider monodisperse spherical beads of diameter d and mass m , initially randomly located between the up-per and lower plates of a box of a given size l x × l y × l z .The lower plate is immobile and the upper one may movealong x in response to a given applied shear stress σ ,i.e., under a constant force F = σl x l y (periodic bound-ary conditions are used along the x and y axes). Grainsare studied by MD simulations of a well known linearspring-dashpot model (L3 of Ref. [11]), including particlerotations and static friction: normal interaction betweentwo beads is characterized by the elastic and viscoelas-tic constants k n = 2 · and γ n = 50, given in stan-dard units [12]; static friction is implemented by keepingtrack of the elastic shear displacement throughout thelifetime of a contact. The tangential elastic constant is k t = 5 . · ; the static friction coefficient µ = 0 .
1; therestitution coefficient is 0 .
88. The box top and bottomplates are made of spherical beads, randomly ‘glued’ atheight z = l z − / ε and at z = 1 / ε , respec-tively, where ε is randomly chosen for each bead in theset [ − / , / m . Weexplored a range of box sizes, and below we illustrate theresults by considering the case with l x = l y = 16, and l z = 8 where finite size effects are absent. The smallvalue of l x and l y , taken here to allow simulations ontime scales long enough, is partially compensated by ourperiodic boundary conditions ( l z has a value typical forexperiments [7, 8, 10]). To avoid gravity-induced com-paction, gravitational acceleration of grains is neglectedas in density-matching liquid experiments [10].For each value of the grain volume fraction φ and ofthe applied force F we have performed 10 simulations us-ing different initial disordered states, always finding thesame typical velocities as well as the same long time state(we have performed simulations up to t = 8000): the sys- V s h ea r ( gd ) / t (d/g) φ = 0.614φ = 0.608φ = 0.593 (a)(b)(c) JammedDisordered OrderedDisordered Disordered
FIG. 1: Velocity of the system upper plate, v shear , as a func-tion of time, t , for the shown values of the volume fraction, φ ,sheared with a constant force F = 4 · . The system startsfrom a disordered initial pack. Panel (a)
When φ is smallenough the system flows in a disordered configuration (seetext and Fig. 2 left panel). Panel (b)
When φ is increasedthe disordered flow is only transient, as the system has anabrupt transition to an ordered flow with a reduced viscosity(see text and Fig. 2 right panel). Panel (c)
At even higher φ the transient disordered flow has a transition to a jammedconfiguration where v shear = 0. tem rheology is therefore characterized by the value of F and on φ (or equivalently on grain number N [19]). Inthe φ and F range of values we explored, which is closeto the hard-particle limit as the maximum deformationof a particle is δd/d < − , three qualitatively differ-ent regimes are usually observed in the system dynam-ics, corresponding to different behaviors of the measuredupper plate velocity, v shear ( φ, F ). They are summarizedin Fig. 1 which shows v shear as a function of time, t , atincreasing values of the volume fraction, φ = 0 . , . . F = 4 · .When φ is small enough (Fig. 1a), the system flows in astationary disordered state (as the one depicted in the leftpanel of Fig.2) from the initial random configuration and v shear fluctuates around a constant value depending on φ and F . The degree of ordering of the system is usuallyquantified by the amplitude of the peak of the structurefactor S ( k ) for, say, k = (0 , , π/d ); in the disorderedregion, S ( k ) has no peaks at all and takes values typicalto disordered arrangements, as shown in Fig.5b.At higher φ values, when F is strong enough (Fig. 1b),the system is trapped for a certain time in a transientstate where v shear keeps a constant value up to a momentwhen it suddenly jumps to a higher stationary plateau.Correspondingly, the system exits the disordered flowingstate (left panel of Fig.2) and develops a layered andpartially ordered structure (shown in the right panel ofFig.2). z zx xx xy y FIG. 2: Vertical ( ZX , top) and horizontal ( Y X , bottom)snapshots of sections of the system for the run of Fig. 1b.The left and the right panels show sections taken at t = 400,when the system flow is disordered, and at t = 800, whenthe flow is ordered. The ordered state is characterized by theformation of crystal-like layers in the XY plane, shifted along z . Dark beads form the boundaries (the shown circles havevarious sizes as they are different cuts of our monodispersespheres). The high φ scenario is drastically changed if the driv-ing force F is not strong enough (see Fig. 1c): after atransient flow the system jams in a state as disordered asthe initial one, and the upper plate velocity becomes zero.The above observations highlight that disordered statescan be either stationary, or transient. The latter statescan undergo transitions towards ordered flowing states,signalled by a marked increase of the shear velocity, oralternatively freeze towards jammed packs.As the velocity profiles appear to be approximatelylinear in the considered cases, the system resistance toflow can be quantified via an effective viscosity η = F/v shear ( φ, F ). In Fig. 5c we depict the ‘flow diagram’ ofthe suspension in the ( φ, F )-plane. We first consider theflow in the disordered region. Figure 3 shows the depen-dence of the effective viscosity on the applied shear force F , for several values of the volume fraction of the system, φ . At low φ and small applied shear stress, the viscos-ity increases to a good approximation with a power law, η ∝ F α , where the exponent is α ≃ .
5. Accordingly,the relation between shear stress ( σ ∝ F ) and shear rate( ˙ γ ∝ v shear ) is σ ∝ ˙ γ β , with β = 1 / (1 − α ) ≃
2, ingood agreement with the prediction of Bagnold scaling( β = 2).Within the disordered region, for a given force F , thedynamics of the system strongly slows down as φ in-creases. This is apparent from Fig. 4 where we showthat, at constant F , the effective viscosity diverges witha power law η ( φ ) ∝ ( φ c ( F ) − φ ) − b as the packing frac-tion increases. Such a relation recorded in the disorderedregime anticipates the presence of a jamming critical vol- F (mg) η ( m ( g / d ) / ) F (mg) η ( m ( g / d ) . ) φ = 0.605φ = 0.608φ = 0.611φ = 0.614 F F φ FIG. 3: The effective viscosity, η = F/v shear ( F ), as a func-tion of the applied external force F for different values of thesample volume fraction ( φ = 0 . − . F the viscosity increases, ata good approximation, as a power law in F α , with α ≃ . σ ∝ ˙ γ β , with β ≃ Inset:
Asimilar power law scaling is found for η ( F ) is the ordered flowregime, but almost no dependence with φ is observed. φ c - φ η ( φ ) / η ( φ c - . ) F = 4 10 F = 4 10 F = 4 10 F = 4 10 F = 1 10 F = 2 10 F φ c FIG. 4: At a given force F , the effective viscosity measuredin the disordered flow regime diverges as a power law as thevolume fraction of the sample approaches a jamming criticalthreshold φ c ( F ). Interestingly, data collected at several F values collapse on a master power law, suggesting that itsexponent is F independent. Here, we have scaled the datain such a way that η ( φ c − .
01) = 1.
Inset: the fitted valueto the critical volume fraction, φ c , increases with the appliedexternal shear force. ume fraction φ c ( F ) (see Fig. 4, inset). Interestingly, thefit of φ c ( F ) from the small φ -value region always overes-timates the volume fraction where we indeed observe thesystem to jam (see Fig. 5), in analogy to glass formingmaterials [13]. Fig. 4 also shows that data recorded atdifferent F collapse one onto the other when plotted as η ( - m ( g / d ) / ) S ( k ) (a)(b) DisorderedOrderedMetastable DisorderedMetastableOrdered φ −1 F DisorderedJammedOrdered (c)
FIG. 5:
Panel (a) and (b) show the dependence of theeffective viscosity η and of the structure factor S ( k ) (for k =(0 , , π/d )) on the inverse volume fraction, φ − , for F =10 . In the low φ regime the system flow is disordered as S ( k ) has typically small values. By increasing φ , there isan initial transient disordered flow ( S ( k ) small) with growing η , followed by a sharp transition to an ordered flow (with amuch higher S ( k ), see text) with a sharp drop in the viscosity(see Fig. 1b). Panel (c):
The system ‘flow diagram’ showingthe long time flow regime as a function of the inverse volumefraction φ − and of the applied force F : circles mark theregion of disordered flow, open squares mark ordered flow,and diamonds the jammed region. Stars indicate the inversecritical volume fraction φ − c defined from the fit in Fig. 4. a function of φ c ( F ) − φ , suggesting that the exponent b ≃ . F . This kind of power law de-pendence is frequently found in colloidal systems undershear [16], with the exponent b varying between 1 and 2.We now consider the regime observed at higher φ where, under shearing, the viscosity exhibits a sharptransition from a disordered to a faster ordered flow asin Fig 1b. The change to the faster flow, with the corre-sponding reduction of the effective viscosity, is associatedto the formation of partially ordered layers parallel to theshearing plate (Fig. 2), as previously observed in experi-ments on granular and colloidal suspensions [8, 14, 15].The connection between ordering and viscosity reduc-tion is apparent in Fig. 5a,b, where we show the depen-dence of both η and S ( k ) on the inverse volume fractionof the system (squares indicate values recorded in thetransient regime, circles the asymptotic ones). In the low φ region, before the ordering transition, S ( k ) has smallvalues typical to disordered arrangements. At higher φ ,the system can flow disordered, with increasing viscosityand S ( k ) having still comparable small values. Eventu-ally, when a sharp drop in the viscosity is observed, it isaccompanied by a clear increase of the structure factor S ( k ), pointing out the presence of order in the system,as seen in the right panel of Fig. 2. In the ordered flowregion, η is almost insensitive to φ in the range here ex-plored (see Fig.5a) and has an approximate Bagnold-likescaling in F , η ∝ F α , with α ≃ . φ regime, if the applied shear stress is not strongenough, the initial transient disordered flow does not gen-erally result in an ordered stationary flow, since a fulldynamical arrest (see Fig. 5c) with v shear = 0 is found.Jamming occurs when the system disordered configura-tions are kinetically trapped in states where further shear(at the given value of F ) becomes impossible, the un-derlying mechanism being still unclear [17, 18]. Such atrapping is broken when F is high enough and orderedflow appears.The flow diagram of Fig. 5c summarizes the propertiesof the long time states reached by the system whenthe initial condition is disordered, and suggests ananalogy between the ( φ, F ) flow diagram and the ( φ, T )phase diagram of usual thermal systems, the disordered,ordered and jammed states corresponding respectivelyto the liquid, crystalline and glassy phases. To reinforcethis analogy we have investigated the behavior of asystem initially prepared in an ordered configuration.Interestingly, we have found that in the region of low F and high φ , where jamming was previously found,the ordered flow is not arrested and keeps going witha finite value of η (although the system might jam attimes longer than those we can investigate). In theother regions of the flow diagram, on the contrary, thelong time state of the system is that shown in Fig. 5cregardless of the initial conditions. This dependence onthe initial conditions is analogous to that observed inthermal systems, where the cristalline phase is stable inthe region where the glassy phase is found.This work has been supported by EU Network MRTN- CT-2003-504712. ∗ Electronic address: [email protected][1] J. B. Knight, C. G. Fandrich, C. N. Lau, H. M. Jaeger,and S. R. Nagel, Phys. Rev. E , 3957 (1995).[2] M. Nicodemi, A. Coniglio, H. J. Herrmann, Phys. Rev.E , 3962 (1997).[3] P. Richard, M. Nicodemi, R. Delannay, P. Ribi`ere, andD. Bideau, Nature Materials , 121 (2005).[4] A. J. Liu, S. R. Nagel, Nature , 21 (1998).[5] A. Coniglio, A. Fierro, H. J. Herrmann, and M. Nicodemi(eds), Unifying Concepts in Granular Media and Glasses(Elsevier, Amsterdam, 2004).[6] E. I. Corwin, H. M. Jaeger, and S. R. Nagel, Nature ,1075 (2005).[7] K. E. Daniels and R. P. Behringer, Phys. Rev. Lett. ,168001 (2005).[8] J.-C. Tsai, G. A. Voth, and J. P. Gollub, Phys. Rev. Lett. , 064301 (2003).[9] E. Bertrand, J. Bibette, and V. Schmitt, Phys. Rev. E , 060401(R) (2002).[10] D. J. Pine (talk, KITP Program on Granular Physics,2005).[11] L. E. Silbert, D. Ertas, G. S. Grest, T. C. Halsey,D. Levine, and S. J. Plimpton, Phys. Rev. E , 051302(2001).[12] Distances, times, velocities, forces, elastic constants, andstresses are, respectively, measured in units of d , t = p d/g , v = √ gd , F = mg , k n = mg/d , and σ = mg/d , where g is gravity acceleration.[13] P. G. Debenedetti and F. H. Stillinger, Nature , 259(2001).[14] J.-C. Tsai and J. P. Gollub, Phys. Rev. E , 031303(2004). Phys. Rev. E , 051304 (2005).[15] R. L. Hoffman, J. Colloid Interface Sci. , 491 (1974);J. J. Erpenbeck, Phys. Rev. Lett. , 1333 (1984);R. L. Hoffman, J. Rheol. , 111 (1998); B. J. Acker-son and P.N. Pusey, Phys. Rev. Lett. , 1033 (1988);L. B. Chen et al. , Phys. Rev. Lett. , 688 (1992).[16] I. M. Krieger, Adv. Colloid Interface Sci. , 111 (1972).[17] M. E. Cates, J. P. Wittmer, J.-P. Bouchaud, and P.Claudin, Phys. Rev. Lett. , 1841 (1998).[18] R. S. Farr, J. R. Melrose, and R. C. Ball, Phys. Rev. E , 7203 (1997).[19] Due to the box rough boundaries (and finite size of theinvestigated system) the volume fraction definition has adegree of ambiguity. We use φ = Nv g / ( l x l y l z ), where v g is the volume of a bead. Other definitions, as for instance φ ′ = Nv g / ( l x l y l z − V b ), where V bb