Dense , collisional , shearing flows of compliant spheres

We outline the development of theory to describe, dense, collisional shearing flows of identical compliant spheres. We begin with two simple theories: one for rigid, nearly elastic spheres that interact through instantaneous, binary collisions; the other for compliant spheres that interact through multiple, enduring contacts. We then join the two extremes by adding compliance to the collisions and collisions to the spheres in enduring contact. Finally, we compare the predictions of the resulting theory with the results of discrete numerical simulations of steady, homogeneous shearing of compliant frictional spheres.


Introduction
In steady, homogeneous shearing flows of identical spheres at volume fractions less than about 0.49, spheres interact through collisions that can be regarded as instantaneous, binary, and uncorrelated. In this case, methods from the kinetic theory of dense gases that take into account the energy loss in a collision [1][2][3] can be used to predict the relationship between the forces necessary to maintain the flow and the rate of shearing. [4][5][6] Above a volume fraction of 0.49, at which a firstorder phase transition between disordered and ordered states in an equilibrated system of colliding hard spheres is first possible [7,8], molecular dynamics simulations show that correlations between collisions begin to influence the relationship between the components of the stress and the shear rate [9][10][11]. The introduction of an additional length scale in the relation for the rate of collisional dissipation of fluctuation energy associated with the size of clusters of interacting spheres modifies the stress relations in an appropriate way [12][13][14][15]. The length scale is determined by the competition between the orienting influence of the flow and the randomizing influence of the collisions, using a local balance between the rates of production and dissipation of fluctuation energy. This approach has been tested against discrete element simulations of steady flows in a variety of flow configurations [16][17][18][19].
These extended stress relations apply until the mean separation distance between the edges of the spheres vanishes, at least along the direction of principal compression, at which point, the stresses for hard spheres become singular [19]. The volume fraction at which this occurs is less than that for random close packing, φ c = 0.64, and is seen in numerical simulation to depend on the coefficient of sliding friction [20].
Replacing a rigid contact with a compliant contact has three consequences: it permits the introduction of a time associated with the duration of a collision before the hard-sphere singularity [21]; it relaxes the singularity of the stresses; and it allows prediction of the stresses at volume fractions greater than that at the singularity. The stresses at volume fractions beyond the hard-sphere singularity have parts that depend on the shear rate and parts that depend on the deformation of the contact. We refer to these as the rate-dependent and rate-independent parts of the stress [22,23].
Here, we focus on continuous, steady, homogeneous shearing of a dense aggregate of identical spheres over a range of solid volume fraction, , φ above 0.49. We first consider theories that are appropriate in limits on either side of, but distant from, the critical volume fraction c φ at which chains of spheres first span the flow. For φ < φ c , we outline the simplest theory for rigid spheres that interact through instantaneous, frictionless binary collisions [24]. Then, we indicate how this theory must be modified to apply to realistic situations. For , > c φ φ we treat deformable spheres that interact through a constantly changing network of longer lasting contacts with multiple neighbors [25]. Finally, we join the two extremes over the entire range of dense volume fractions by incorporating contact deformation in the collisional regime and impulsive momentum transfer in the deformational regime [26].
We treat spheres with a diameter, ,

Simple kinetic theory
The simplest kinetic theory [24] assumes instantaneous binary collisions, between a typical pair of frictionless, rigid spheres, with an impulse J exerted by sphere 2 upon sphere 1, and the change in velocities of the two spheres in a collision given by where the prime denotes values after the collision. The component of the relative velocity ′ c 12 ≡ ′ c 1 − ′ c 2 , along the unit vector k directed from the center of sphere 1 to that of 2 is related to that before the collision by a coefficient of restitution e, with 0 1 ≤ < e : The change in kinetic energy is important when calculating the rate of collisional dissipation. Averages over the exchanges of momentum and energy in collisions are carried out by integrating over a complete pair particle distribution that governs the likelihood of collisions between pairs of spheres. In the simplest dense theory, the positions and velocities of colliding pairs are assumed to be completely uncorrelated -the assumption of Molecular Chaos -and the distribution is assumed to be expressed as the product of single particle velocity distributions and a correction to the collision frequency that incorporates the distance between edges, rather than centers.
The single particle velocity distribution is ( ) where the velocity fluctuation ( ) ( ) ≡ − C r c u r is given in terms of the local average velocity u(r), 2 , ≡ ⋅ C C C at the point of contact r and the measure 2 / 3 ≡ T C of the strength of the velocity fluctuations is often referred to as the granular temperature. The frequency correction is given in terms of the volume fraction dependence of the radial distribution function for a contacting pair in a dense aggregate [8]: With these, the averages are carried out using the distribution function in which care is taken to distinguish between the locations of the enters of the colliding pairs and, in the simplest theory, spatial gradients of the average fields enter the theory through series expansions of these about the point of contact. Integration yields expressions for the pressure (10) and the rate of collisional dissipation The last of these may be written as to emphasizes that it is the product of the number of spheres per volume, the energy lost per collision, and the collision frequency. The inverse of the collision frequency is the average time of flight, , f τ between collisions: In a steady, homogeneous shearing flow, the granular temperature is determined through the energy balance: When the expressions for shear stress and the rate of dissipation are employed in this, it may be solved for T: Then, the shear stress and the pressure can be expressed in terms of the shear rate: ( ) ( ) and ( ) 2 2 4(1 ) . 15 (1 ) We note for later comparison that the ratio s/p is independent of φ in the simplest theory.

Simple deformation theory
Pairs of spheres interact through contact forces that have normal and tangential components, P and T, parallel and perpendicular to their line of centers. In the simplest linear theory these are, respectively, proportional to the normal and tangential displacements δ and s of the contact. The normal force exerted by sphere 2 on sphere 1 is while the tangential force is , if 4 , otherwise.
The distribution of the magnitudes of the contact forces is roughly exponential [27]; so, for example, (20) where P is the average normal force. Discrete numerical simulations and physical experiments on photo-elastic disks indicate both the anisotropy and inhomogeneity of the distribution of contact forces in a sheared aggregate. The orientational distribution of contacts is described by a function A(k) of orientation k defined so that ( ) Ω k A d is the average number of contacts per particle within an element dΩ of solid angle centered at k [25]. The average total number of contacts per particle Z is called the coordination number and is given in terms of A(k) by ( ) .
When the contact distribution is isotropic, Anisotropic distributions of contacts are often described using a contact fabric tensor, : a ( ) The stress, , t is defined in terms of the force transmitted across area element ΔΣ with unit normal N by contacting pairs with orientations k within solid angle dΩ [28]: Similarly, for 0 ⋅ k N < . Taking half the sum and integrating over all orientations k gives For example, for an isotropic contact distribution, an average volume strain, Δ, positive in compression, and the contact displacement assumed to be given by the average strain, δ = σΔ / 3, [25] the pressure is In applying this result to continued shearing beyond the critical volume fraction, it's natural to relate the average volume strain to the excess of the volume fraction above the critical: ( ) Numerical simulations of very slow, continuous, homogeneous shearing [29] indicate that the stress is given by Because only a xy is nonzero, , xy s t p and η yet to be determined.

Corrections to simple kinetic theory
In order to make realistic, quantitative predictions for dense, dissipative, shearing flows the simple kinetic theory whose derivation was outlined above has to be modified in a number of ways. First, a more realistic complete pair distribution function must be derived and employed. This has been done for elastic spheres [30]. The resulting modification to the theory for steady, homogeneous shearing outlined above for dissipative spheres [2] requires that in the shear stress (10), the factor 1+e be replaced by the scalar factor J: ( ) Second, we must take into account the rotational degrees of freedom, and the additional loss of translational fluctuation energy to sliding, tangential restitution, and rotational fluctuation energy. This can be done by employing an effective normal coefficient of translation restitution, ε, that is derived using the balance of rotational fluctuation energy to determine the partition between rotational and translation energy in a steady, homogeneous flow. [31,32] The effective coefficient of restitution is simply expressed in terms of the coefficients of sliding friction and tangential restitution in two limits: when all contacts slide in a collision and μ is about 0. 1, 2 9 , 2 2 = − + e π ε μ μ (33) and when all contacts come to rest at some time during a collision and μ is about 1, Third, we must include the dependence seen in discrete numerical simulations of the critical volume fraction on friction. To do this, in g 0 of Eq.  Finally, we must incorporate the influence of velocity correlations that occur above 0.49. = φ At greater volume fractions, the assumption of molecular chaos breaks down and correlated velocities occur in clusters of size L, oriented along the axis of compression. The need for making such a correction is shown in Fig. 3, in which predictions and measurements of the ratio of shear stress to pressure versus volume fraction are shown for several typical systems.
A simple means of doing this is to introduce the cluster size in place of the particle diameter in the rate of collisional dissipation [12][13][14] ( ) 2 3 /2 1/ 2 12 1 The size of the cluster is determined by a balance between the orienting influence of the flow and the randomizing influence of the collisions: [19]. That is, it is singular at random close-packing, not at the critical volume fraction. Then, ( ) where, in contrast to the ratio given by Eq. (16) and (17), there is a dependence on the volume fraction This relation is one of those plotted in Fig. 3, where it is referred to as that of the extended kinetic theory..

Deformation in collisions and collisions in deformation
We next incorporate the deformation of the contacts during a collision into the collisional theory. [21] This provides a way to regularize the singularity in the collisional constitutive relations that without deformation occurs at the critical volume fraction.
When the deformation of a contact is linear, the duration, τ c , of a contact is while we have seen that the time between collisions is Consequently, if we take into account the duration of a collision, the time between collisions increases from τ f to τ f + τ f . In this case, Similarly, 2 15 The need to incorporate dynamics in the range of volume fraction above the critical can be seen in the results of experiments carried out on photo-elastic disks. Here, a dramatic difference in the strength of the fluctuations in stress is seen in Fig. 4 below and above an area fraction that we identify with the critical. Consequently, for , ≥ c φ φ there is still a collis-ional component of the exchange of momentum associated with the energy released by the breaking of chains. We assume that the frequency of these collisions is inversely proportional to the duration of a contact and evaluate the momentum transfer at the critical volume fraction.
We take the pressure in this regime to be the sum of contributions from collisions and deformation -that is, the sum of rate-dependent and rate-independent parts: We treat the shear stress in a similar way, taking into account the relation (30) between the shear stress and pressure when deformations alone occur: In this regime, we assume that the anisotropy of the collisions is the same as that of the deformation and take η to be the ratio of the collisional shear stress to the collisional pressure at φ c : Finally, in order to enable the determination of the strength of the velocity fluctuations at volume fractions greater than the critical, we adopt similar assumptions for the frequency of collision and the volume fraction in the rate of collisional dissipation: Then, in simple shearing flows, both below and above the critical volume fraction, T may be obtained as a solution to the energy balance (14). In this energy balance, the elastic part of the shear stress above the critical volume fraction is assumed to be recoverable and not to contribute to the rate of working.

Results
With T determined on either side of the critical volume fraction, the stress may be calculated and compared to the results of discrete numerical simulations [20,33,34]. The comparison for the pressure and shear stress, A more stringent test results from plotting the ratio of the shear stress to the pressure. In this case, because he stiffness is assumed to enter into both quantities in the same way, the ratio is independent of stiffness. In contrast. the measured data is sensitive to the stiffness, especially for the more compliant spheres. The agreement for the stiffer spheres is relatively good. Predictions of the normalized pressure versus volume fraction for e = 0.7, μ = 0.5, values of the dimensionless contact stiffness ranging from 10 to 10 11 (dashed curves) and data from the discrete numerical simulations [20,33,34].

Fig 6:
Predictions of the normalized shear stress versus volume fraction for e = 0.7, μ = 0.5, values of the dimensionless contact stiffness ranging from 10 to 10 7 (dashed curves) and data from the discrete numerical simulations [20,33,34]. Predictions of the ratio of shear stress to pressure versus volume fraction (solid curve) for e = 0.7 and μ = 0.5 and data from the discrete numerical simulations [20,33,34] for values of the dimensionless contact stiffness ranging from 10 to 10 7 (symbols).

Conclusions
We have outlined the derivation of a theory for steady, homogeneous shearing flows of dense aggregates of deformable, dissipative spheres. The incorporation of contact deformation in the collisional regime and impulsive momentum transfer in the deformational regime permits a smooth transition through the critical volume fraction and relatively good descriptions of the stresses on either side of it. The predictions compare well over the range of volume fraction before and after the hard sphere singularity and over nine orders of magnitude of the contact stiffness. The model requires no parameters other than the contact stiffness, and the coefficients of restitution and sliding friction of the spheres.