Evolution of Internal Granular Structure at the Flow-Arrest Transition
The evolution of the internal granular structure in shear-arrested and shear-flowing states of granular materials is characterized using fabric tensors as descriptors of the internal contact and force networks. When a dilute system of frictional grains is subjected to a constant pressure and shear stress, the internal stress anisotropy is well-predicted from the anisotropy of its contact and force networks even during transient evolution prior to steady shear flow or shear arrest. Although the onset of shear arrest is a highly stochastic process, the fabric tensors upon arrest are distributed around nearly equal contributions of force and contact network anisotropy to the total stress anisotropy. The distribution becomes seemingly narrower with increasing system size. The anisotropy of the contact network in shear-arrested states is reminiscent of the fabric anisotropy observed in shear-jammed packings.
Keywords:Shear Jamming Fabric Tensor Force Network Granular Friction Critical State
Granular materials respond to applied stress in several remarkable ways. When subjected to a critical ratio of shear stress and pressure, they can plastically deform with accompanying dilation to reach a critical state beyond which further quasistatic deformation is volume preserving and rate independent schofield1968critical . When the applied stress is sufficiently large, granular materials exhibit inertial frictional rheology, which is described by a proportionality between the internal shear stress and rate of shear deformation Jop:2006gh ; Salerno:2018en . Therefore, a well-defined friction-dependent critical ratio of shear stress and pressure demarcates the plastic or fluid-like deformation of granular materials from a static solid-like behavior singh2013effect ; Srivastava:2019bc . In addition to stress, the solids volume fraction of granular materials also controls their mechanical behavior: below a friction-dependent critical volume fraction (also known as the jamming transition), granular materials display fluid-like behavior, whereas above this volume fraction, the grains are jammed and collectively exhibit mechanical rigidity Liu:2010jx ; o2003jamming ; Srivastava:2019bc . However, this picture has been complicated by recent observations that granular materials can jam below the so-called jamming volume fraction by subjecting them to external shear strain Bi:2011bt ; vinutha2016disentangling ; Chen:2018cf ; clark2018critical .
The application of shear to a granular material can introduce significant geometrical anisotropy within its microstructure. When sheared at its boundaries, the constituent grains tend to align preferentially in the direction of shear Majmudar:2005bn . Such directional alignment of the grains also results in the partial alignment of the force network giving rise to the formation of directional force chains—as depicted in Figs. 1(b) and (c)—which carry a majority of the shear stress Radjai:1998hn . Consequently, several metrics have been proposed to quantify structural anisotropy in granular materials, of which the most commonly utilized is a second order symmetric fabric tensor Radjai:2012fc ; Rothenburg:1989iq . Early work by Rothenburg and Bathurst Bathurst:1990bw demonstrated that the shear stress within granular materials can be partitioned into contributions from various fabric tensors, thereby providing a direct connection between bulk mechanics and inherent granular structure. Several later studies have developed continuum models for granular plasticity and rheology incorporating the granular fabric tensor as a key internal variable Radjai:2012fc ; SUN:2011ks . The effect of structural anisotropy on the bulk mechanical response is even more pronounced when granular materials are subjected to complex loading paths.
Prior studies have comprehensively explored the evolution of fabric anisotropy during quasi-static deformation of granular materials leading to critical state and failure Guo:2013gs . The structural anisotropies in granular materials during steady inertial flows as a function of applied shear rate have also been demonstrated Azema:2014ej . However, the transient evolution of granular fabric during its transitions between solid-like and fluid-like states remains largely unexplored. Such studies will especially aid the development of advanced constitutive laws that can predict time-dependent state of stress and deformation in granular materials as they transition between fluid-like and solid-like behaviors.
Recent simulations have demonstrated that the stress-induced transient flow-to-arrest transition in frictional granular materials is a highly stochastic process with a very wide distribution of arrest times, along with power-law divergences of arrest statistics near a critical shear stress, all of which are sensitive to interparticle friction Srivastava:2019bc . In this work we elucidate the transient evolution of the granular fabric during such a flow-arrest transition. We provide evidence that stress partitioning into various fabric components—which has been successfully applied to quasistatic deformation Guo:2013gs and steady inertial flows Azema:2014ej —remains valid even during the transient flow-arrest transition. We highlight that despite the stochastic nature of this transition, the structural anisotropy of the granular material upon arrest remains largely invariant to these statistics, indicating the presence of a unique shear-arrested state for a given loading protocol.
Dilute granular systems—composed of slightly dispersed spheres (within % of average diameter ) with a Hookean spring stiffness and interparticle friction —are prepared at an initial volume fraction . These systems are subjected to a constant shear stress and constant pressure , such that the total applied stress follows: (i) , (ii) for and , and (iii) for all other indices . While the pressure is fixed at to simulate hard particle regime at high volume fraction, is modulated to simulate two distinct long-time scenarios: (i) steady shear flow at a constant strain rate , and (ii) shear arrest where . Here is the second invariant of the 3D strain rate tensor . The simulation method is described schematically in Fig. 1(a).
A constant state of stress at the periodic system boundaries is maintained using the modularly-invariant dynamical equations Martyna:1994ki of the Parrinello-Rahman (PR) method Parrinello:1981it within LAMMPS Plimpton:1995fc , which allows the system boundaries to deform in all possible ways, i.e., shear, dilation/compaction, and rotation, in response to external stress. We simulated granular systems with two values of interparticle friction and representing regimes of low and high friction respectively. More than distinct starting states for and more than distinct starting states for were simulated to provide robust statistics for shear-arrested states, following the observations in Ref. Srivastava:2019bc that the time for the onset of arrest in a flowing granular material is highly unpredictable. We also performed simulations for intermediate interparticle friction , and observed similar phenomenology as low and high friction systems.
The internal (Cauchy) stress state of a granular material is tracked during its evolution under external applied stress from the interparticle contact forces as: , where and are the force and branch vectors between two contacting particles, and the summation is over all contacts. The kinetic contribution to the total internal stress is insignificant compared to contact stress contribution. Two relevant invariants of Cauchy stress are measured, (i) pressure and (ii) shear , and the total stress anisotropy of a granular material is defined by the ratio .
3 Bulk Response
Figure 2 demonstrates the temporal evolution of key bulk variables, i.e., strain rate , volume fraction and average coordination number , for the same starting state of particles subjected to different values of stress ratio , which highlight two distinct dynamical scenarios: (i) steady shear flow at a constant mean for high shear for low friction and for high friction, and (b) shear arrest at a distinct time for low shear (blue) and (red) for low friction, and (blue) and (red) for high friction. The shear arrest is followed by for . The onset time for shear arrest is defined as the time at which drops below a minimal value of , as indicated in Fig. 2. The two distinct responses, which characterize the flow-arrest transition, are observed for both low interparticle friction in Fig. 2(a) and high interparticle friction in Fig. 2(b).
At early times the granular system rapidly compacts to a friction-dependent volume fraction, which remains invariant with time for both shear arrest and shear flow. For low interparticle friction, , which is similar to the random close pack density of frictionless spheres, while for high friction, which is similar to the critical volume fraction in the critical state plasticity theory schofield1968critical . Although the bulk strain rate exhibits a bifurcation of several orders of magnitude depending upon the applied shear stress, volume fraction is insensitive to this transition.
Unlike , the coordination number exhibits differences between long-time steady shear flow and shear arrest. The early time compaction of the granular system is accompanied by a rapid rise in , after which it fluctuates around a mean value during a period of transient flow at a constant volume, prior to steady shear flow or arrest. The mean value of during transient flow is dependent on both shear stress and friction, with a moderately larger for lower shear stress and lower friction, consistent with previous observations on steady inertial flows Azema:2014ej . While continues to fluctuate around a mean value for systems entering steady shear flow, the onset of shear arrest is characterized by a sudden change in to a shear-independent value regardless of the onset of arrest, as shown in Fig. 2. For low friction particles , whereas for high friction particles. Hence, upon arrest these values are similar to the coordination number obtained for isotropic jamming of frictional and frictionless particles o2003jamming ; Silbert.2010 ; Otsuki:2011et . These estimates for were averaged from multiple shear arrest simulations, and the errors represent the standard deviation. These results demonstrate the uniqueness of coordination number for a given interparticle friction for all shear-arrested states attained through the same loading protocol regardless of arrest statistics and applied shear stress.
4 Fabric Evolution
Volume fraction and coordination number describe the mechanical behavior of granular materials in the vicinity of isotropic fluid to solid jamming transition, as observed by the scaling of these quantities with system pressure Goodrich:2016ch . However, the anisotropy induced within contact and force networks in granular materials upon the application of shear stress requires additional characterization beyond such isotropic measures. Fabric tensors provide a tensorial representation of the contact and force networks in granular materials, and consequently, the inherent structural anisotropy KenIchi:1984gt ; ODA:1982dt . The geometry of the contact network can be described by the distribution of orientation of contact normals in terms of a fabric tensor:
where is the -th component of a contact normal , is the number of contacts, and is the volume of the system. It is often convenient to express the contact orientation distribution in terms of a second-order Fourier expansion of as Rothenburg:1989iq ; Bathurst:1990bw :
where is a symmetric and deviatoric contact anisotropy tensor whose magnitude quantifies the magnitude of geometrical anisotropy of the contact network.
Similarly, the force-carrying network in a granular material presents a source of mechanical anisotropy, which can be quantified using tensors representing the force fabric. Particularly, a normal force anisotropy tensor and a tangential force anisotropy tensor can be defined such that orientational distribution of normal forces and tangential forces are represented in their second-order Fourier expansion as Bathurst:1990bw :
where is the mean normal contact force in the granular material.
Prior studies have shown that the total internal stress in granular materials can be partitioned into orientational distributions of contact normal and force vectors as follows Rothenburg:1989iq ; Bathurst:1990bw :
where the integration is over the entire angular domain , is the mean center-to-center distance of contacting particles, and is the orientational distribution of the total contact force vector, which can be decomposed into normal and tangential contributions, as shown in Eq. 3. Combining Eqs. 2, 3 and 4 along with the definition of stress anisotropy reveals the following stress-force-fabric relation Guo:2013gs :
where , and are the second invariants of the deviatoric fabric tensors , and respectively. The equality in Eq. 5 holds true only if the three fabric tensors are co-axial with the stress tensor. Regardless of the applied stress or interparticle friction, we observe co-axiality of all fabric tensors with the stress tensor at all times except during early transients when the contact network is still evolving.
Figure 3 demonstrates the evolution of three scaled fabric anisotropies , , and for the same bulk response depicted in Fig. 2 (same color) for low and high friction at different values of applied . After the early transients during which the contact and force networks develop as a result of rapid compaction, the stress-force-fabric relation in Eq. 5 is satisfied at all times even during the transient flowing state before the system either arrests (red and blue curves) or continues flowing indefinitely in a steady flow (black curves), thereby demonstrating co-axiality of the stress tensor with all the fabric tensors. However, the relative contribution of fabric anisotropies to the total internal stress anisotropy varies significantly with time. Early in the transient flowing state, contact anisotropy dominates the total system anisotropy, whereas the force network remains fairly isotropic, indicating that although a large number of contacts begin to align in the direction of shear, rapid breaking and forming of new contacts does not sustain the directionality of the force network. As time progresses, the relative contribution of contact (force) network decreases (increases) towards seemingly steady values, in which and provide nearly equal contribution to the total stress anisotropy. Figs. 1(b) and (c) respectively depict the force network for a representative simulation during early and late times, and the isotropic (anisotropic) nature of the force network at early (late) time is apparent. The anisotropy of the tangential force network remains significantly smaller compared to contact and normal force network, and is expectedly non-existent for the case of lower friction. The near-equality of and for quasistatic flows near the flow-arrest transition is consistent with previous observations of fabric evolution in steady inertial flows Azema:2014ej .
The fabric anisotropies of the steady flowing systems continue to fluctuate around a mean value at long times. For the systems that arrest, the three fabric anisotropies instantaneously become time-invariant, thereby indicating an arrested microstructure. Interestingly, regardless of interparticle friction or applied shear, upon arrest, while is relatively small, and the equality of Eq. 5 holds. This indicates that the shear-arrested states can be uniquely defined by friction-dependent volume fraction, coordination number and fabric anisotropies, with nearly equivalent contribution of contact and force network anisotropies. Additionally, it is expected that the shear-arrested fabric anisotropies depend significantly upon the loading protocol, which is same for all systems in the present simulations.
5 Fabric States at Shear Arrest
Recent simulations have demonstrated that the flow to arrest transition in frictional granular materials is a highly stochastic process. The onset of arrest below critical shear at time exhibits a heavy-tailed distribution, whose moments diverge near the critical shear threshold Srivastava:2019bc . As an example, Fig. 4(a) demonstrates shear arrest for three different starting states of granular systems subjected to the same sub-critical shear stress , and with the same interparticle friction. Although all the three systems eventually arrest, their onset of arrest times vary by more than an order of magnitude. Here we examine whether this stochasticity associated with shear arrest is reflected in the response of the contact and force fabric of the granular material.
The evolution of three scaled fabric anisotropies for three cases of shear arrest in Fig. 4(a) is shown in Fig. 4(b). The response of the fabric after early transients and before shear arrest is qualitatively equivalent in all the three cases: at early times, dominates the total stress anisotropy, and as time progresses, and evolve to nearly equivalent values. However, more crucially, the state of the fabric upon arrest is nearly the same for all the three cases, with , thereby indicating that the fabric is independent of the time to arrest.
We examined the statistics of the state of the fabric upon arrest in greater detail by simulating multiple starting states towards shear arrest for the same value of applied shear. Figs. 5(a) and (b) respectively highlight the scatter of and for all the cases of shear arrest for low friction and high friction respectively. The fabric anisotropies are scattered around (non-scaled) mean values of , and for low friction, and , and for high friction. Expectedly, higher interparticle friction results in greater microstructural anisotropy which also leads to a larger stress anisotropy through the stress-force-fabric relation. For smaller system sizes, the scatter of the fabric anisotropies is larger, as shown in Fig. 5, although they are scattered around a similar mean value. This signifies the importance of finite system size effects, and it also possibly indicates that for an infinitely large system, the distribution will vanish around unique values of fabric anisotropies. A systematic study of system size effects will better elucidate the uniqueness of arrested fabric.
The narrow distribution of the fabric anisotropies indicates the uniqueness of the fabric at shear arrest. Unlike and which are slightly negatively correlated for finite system sizes (see dashed linear correlation fits in Figs.5(a) and (b)) due to the constraints of Eq.5, we find no correlation between arrest times and fabric anisotropies. Figs.5(c) and (d) highlight the scatter of and , and no correlation between the two is observed, thereby implying that no matter how long the granular system takes to arrest from its flowing state, it will arrest to a uniquely defined internal granular structure.
The results presented here have two important implications. Shear jamming of granular materials was first observed in seminal experiments by Bi et al. Bi:2011bt , which were later simulated for frictional Otsuki:2011et and frictionless Chen:2018cf ; clark2018critical particles through various strain-controlled protocols. These experiments observed the presence of fragile states in addition to shear-jammed states, and a subsequent theory was developed to identify the transition between them Sarkar:2013fj . The shear-arrested states simulated in this study through a completely stress (and pressure)-controlled protocol bear resemblance to those observed in previous studies. They are highly anisotropic in their fabric, and importantly, their internal state of stress is uniquely defined by the stress invariants and , unlike strain-controlled methods where the internal stress is not uniquely defined Chen:2018cf . Further analysis is required to highlight similarities and differences between such shear-arrested and shear-jammed states of granular matter. Additionally, it remains to be seen whether these shear-arrested states exhibit signs of fragility Bi:2011bt , i.e., any directional instability to shear. Recent stress-controlled simulations on frictional suspensions at a constant volume provide insight into the differences between shear-jammed and fragile states of granular materials seto:arxiv .
Recent advances in continuum modeling of granular materials have begun to include anisotropic effects in constitutive models Radjai:2017bl . Particularly, advances to the classical critical state theory have included the role of fabric Li:2012kb ; Li:2015ka ; Dafalias:2016gh , in addition to volume fraction and stress ratio, to uniquely define a critical state of granular plasticity at which granular materials deform indefinitely at a constant volume and stress. The uniqueness of critical state granular fabric has also been backed by recent simulations Wang:2017ch . The simulations presented here highlight the opposite phenomenon of the arrest of a flowing granular material. Currently, no constitutive models exist that can predict transient granular flow prior to shear arrest, and stochastic nature of shear arrest itself. However, the preceding results, which elucidate the unique state of granular fabric upon shear arrest, highlight the importance of including fabric tensor as a key parameter in the future development of such constitutive models.
Acknowledgements.This work was performed, in part, at the Center for Integrated Nanotechnologies, an Office of Science User Facility operated for the U.S. Department of Energy (DOE) Office of Science. Sandia National Laboratories is a multi-mission laboratory managed and operated by National Technology and Engineering Solutions of Sandia, LLC., a wholly owned subsidiary of Honeywell International, Inc., for the U.S. DOE’s National Nuclear Security Administration under contract DE-NA-0003525. The views expressed in the article do not necessarily represent the views of the U.S. DOE or the United States Government.
- (1) A. Schofield, P. Wroth, Critical State Soil Mechanics, vol. 310 (McGraw-Hill London, 1968)
- (2) P. Jop, Y. Forterre, O. Pouliquen, Nature 441, 727 (2006)
- (3) K.M. Salerno, D.S. Bolintineanu, G.S. Grest, J.B. Lechman, S.J. Plimpton, I. Srivastava, L.E. Silbert, Phys. Rev. E 98, 050901 (2018)
- (4) A. Singh, V. Magnanimo, S. Luding, AIP Conference Proceedings, vol. 1542 (2013)
- (5) I. Srivastava, L.E. Silbert, G.S. Grest, J.B. Lechman, Phys. Rev. Lett. 122, 048003 (2019)
- (6) A.J. Liu, S.R. Nagel, Ann. Rev. Cond. Matter Phys. 1, 347 (2010)
- (7) C.S. O’Hern, L.E. Silbert, A.J. Liu, S.R. Nagel, Phys. Rev. E 68, 011306 (2003)
- (8) D. Bi, J. Zhang, B. Chakraborty, R.P. Behringer, Nature 480, 355 (2011)
- (9) H. Vinutha, S. Sastry, Nature Phys. 12, 578 (2016)
- (10) S. Chen, T. Bertrand, W. Jin, M.D. Shattuck, C.S. O’Hern, Phys. Rev. E 98, 042906 (2018)
- (11) A.H. Clark, J.D. Thompson, M.D. Shattuck, N.T. Ouellette, C.S. O’Hern, Phys. Rev. E 97, 062901 (2018)
- (12) T.S. Majmudar, R.P. Behringer, Nature 435, 1079 (2005)
- (13) F. Radjaï, D.E. Wolf, M. Jean, J.J. Moreau, Phys. Rev. Lett. 80, 61 (1998)
- (14) F. Radjai, J.Y. Delenne, E. Azéma, S. Roux, Granular Matter 14, 259 (2012)
- (15) L. Rothenburg, R.J. Bathurst, Géotechnique 39, 601 (1989)
- (16) R.J. Bathurst, L. Rothenburg, Mech. Mat. 9, 65 (1990)
- (17) J. Sun, S. Sundaresan, J. Fluid Mech. 682, 590 (2011)
- (18) N. Guo, J. Zhao, Comp. and Geotech. 47, 1 (2013)
- (19) E. Azéma, F. Radjaï, Phys. Rev. Lett. 112, 93 (2014)
- (20) G.J. Martyna, D.J. Tobias, M.L. Klein, J. Chem. Phys. 101, 4177 (1994)
- (21) M. Parrinello, A. Rahman, J. Appl. Phys. 52, 7182 (1981)
- (22) S. Plimpton, J. Comp. Phys. 117, 1 (1995)
- (23) L.E. Silbert, Soft Matter 6, 2918 (2010)
- (24) M. Otsuki, H. Hayakawa, Phys. Rev. E 83, 051301 (2011)
- (25) C.P. Goodrich, A.J. Liu, J.P. Sethna, Proc. Natl. Acad. Sci. U.S.A 113, 9745 (2016)
- (26) K. Ken-Ichi, Int. J. Engg. Sci. 22, 149 (1984)
- (27) M. Oda, Soils and Foundations 22, 96 (1982)
- (28) S. Sarkar, D. Bi, J. Zhang, R.P. Behringer, B. Chakraborty, Phys. Rev. Lett. 111, 068301 (2013)
- (29) R. Seto, A. Singh, B. Chakraborty, M.M. Denn, J.F. Morris, arXiv:1902.04361 (2019)
- (30) F. Radjaï, J.N. Roux, A. Daouadji, J. Engg. Mech. 143, 04017002 (2017)
- (31) X.S. Li, Y.F. Dafalias, J. Engg. Mech. 138, 263 (2012)
- (32) X.S. Li, Y.F. Dafalias, J. Mech. Phys. Solids 78, 141 (2015)
- (33) Y.F. Dafalias, Acta Geotech. 11, 479 (2016)
- (34) R. Wang, P. Fu, J.M. Zhang, Y.F. Dafalias, J. Engg. Mech. 143, 04017117 (2017)