Journal Club Theme of July 2013: Predicting granular flows: A new size-dependent constitutive model
Since Coulomb initiated the topic over 200 years ago, continuum modeling of granular materials has remained an infamously difficult subject. Granular matter is common in everyday life (soil, sand, food grains, pharmaceuticals, etc) and second only to water as the most handled type of industrial material. However, a predictive, general constitutive relation for granular flow is lacking. This has become an expensive world-wide setback, since geotechnical and industrial flows are often on space- and time-scales too large for discrete particle simulation. It begs us to ask: Is there any hope for a "Navier-Stokes" equivalent for sand? With a view toward this question, in this journal club entry I'll outline some encouraging recent results based on a new model obtained with collaborators Georg Koval (of INSA Strasbourg) and recently-former postdoc David Henann (just started as Mech E faculty at Brown). See papers  and  for full details.
Throughout, we limit the discussion to well-developed flows, in order to make the problem more tractable. This focus is part of a bigger strategy --- once we understand the well-developed state, it instructs the long-time behavior when developing critical-state-like models for the transients.
Inertial (local) rheology
Our story begins in the mid-2000's, with the advent of the inertial rheology for granular flow, developed by a group of southern French researchers [3,4]. As reported in my prior journal club entry, the premise of the model can be extracted logically from dimensional arguments in a simple shear cell (see notational details below). In the end, for a given quasi-monodisperse granular composition, the result is a local rheology of the form
Of note, the function f is empirically fit and demonstrates a clear yield criterion; f = 0 if τ/P = μ < μs for material parameter μs. The inertial rheology above has shown to work extremely well for steady-state, homogeneous simple shear flows.
To broaden beyond simple shear, the relation can be written tensorially (i.e. presuming codirectionality and the Drucker-Prager yield surface), and instituted in general geometries as a 3D viscoplastic rheology , or as the flow rule in an elasto-viscoplastic rheology . While these extensions give reasonable results for rapid flows, trouble quickly sets in.
Despite how well the inertial rheology may work in homogeneous simple shear, it fails for slow non-uniform flows, even ones where the steady flow profile is locally simple shearing everywhere (say, annular shear flow) . Of note, the model predictions tend to under-predict the size of flow features; rather, it happens that the grain size itself plays the key role in determining the width of such features. Moreover, we see creeping flow occurring in regions where stress levels are well beneath the μs value obtained from uniform simple shearing, as long as a flow gradient is present. These are tell-tale signs of nonlocality, set by a length-scale associated to the particle size, and consequently, well-developed granular rheology must be recast with nonlocal (or gradient-based) constitutive modeling. Since grains are commonly on our size-scale, the consequences of the nonlocality can be appreciated and noticed with the naked eye, in a truly remarkable fashion. Indeed one does not have to try hard to produce a flow where almost the entire profile is due to size-effects, and differs completely from the inertial rheology prediction.
Extending to a nonlocal flow rule
Our model, the Nonlocal Granular Fluidity (NGF) model, aims at resolving these issues, and can be backed out statistically from a micro-mechanism in which "flow aids flow" . That is, a microscopic rearrangement event at some location can emit stress perturbations that affect material some distance away. Hence, in the continuum limit, it is both the mean applied stress and disturbances from neighboring flow (in the form of a Laplacian term scaled by the particle size) that together determine the shear-rate at some point. Mathematically, NGF resembles an implicit gradient theory, e.g. . Omitting mathematical specifics here (see  for details), we note the following important features of the model:
1) NGF couples the stress/strain-rate relation to a scalar field, the "granular fluidity", which satisfies a PDE calling on the grain-size.
2) NGF reduces to the inertial relation in the absence of flow gradients (as it should).
3) NGF necessitates only one new material parameter; an order-one dimensionless constant we call A, the nonlocal amplitude. The other parameters are borrowed directly from the inertial rheology. In total, NGF uses three material parameters.
The combined results of  and  show that the NGF model predicts flow and stress fields with a newfound level of accuracy, as verified over hundreds of geometries.
Our initial demonstrations in  were in the simplest 2D geometries --- flow in an annular cell, gravity-driven flow down a vertical chute, and dragging of a plate over granular bed. Its predictions matched data from discrete simulations of disk flows (c/o Georg Koval) in the analogous geometries for various different values of system size or gravity, over many orders of magnitude in flow speed. A single calibration of the material parameters was used throughout.
With the success of the 2D prototype cases, the model was recast in 3D, and an Abaqus User-Element was created (c/o David Henann) to permit us to solve the NGF system in arbitrary 3D geometries . We re-calibrated the model parameters for quasi-monodisperse 3D glass beads, and compared model predictions against experimental data for bead flows in many geometries. As a stringent test, we first compared NGF predictions in "split-bottom" flow geometries, a family of flow environments made famous over the last decade for having resisted all previous continuum models. It was shown in  that NGF is the first continuum model to quantitatively capture all features of the flows in these geometries. For example, see below (from ).
(a) A schematic of the split-bottom geometry: Grains fill an annular trough split along its bottom at radius Rs. The outter portion is rotated while holding the inner portion fixed. (b,c) Model predictions for the flow field as measured by the revolution rate ω(r,z) = νθ(r,z) / RΩ for different filling heights H. Note the Heaviside profile at z=0 spreading out z increases. (d) Comparing the predicted flow field on the top surface (z=H) to experimental data  for various filling heights H.
We then applied the same model with same parameters in completely different flow environments and found equally high agreement with experimental data on glass bead flows, reflecting the geometric-generality of NGF. For example, see below (from ).
(Left) Comparison of NGF prediction (solid line) to experimental flow data (symbols) in a 3D, gravity-compacted annular shear apparatus . (Right) Comparison of NGF prediction (solid line) to experimental flow data (symbols) in a 3D plate-dragging geometry in the presence of gravity .
More to do
Despite the encouraging results, there remain several questions to answer. For instance, there are lingering issues about the form of fluidity boundary conditions, how exactly to incorporate transient effects and/or anisotropy internal variables within this framework, and the exact connection between grain properties and the nonlocal amplitude A. Please see the end of our paper  for a comprehensive discussion of the major open issues.
Again, much appreciation to my collaborators Koval and Henann!