Independence of shear strength with particle size dispersity still valid in polyhedral particle assemblies

A very staggering result that has been constantly highlighted in granular media is that the shear strength of granular assemblies is independent of the particle size dispersity. In other words, a packing composed of monodisperse particles has similar strength properties to those of polydisperse systems. This has been shown numerically for the simplified case of disc and polygon assemblies in 2D and spheres in 3D. In this paper, we use three-dimensional contact dynamics simulations to revisit these results for the more complex case of assemblies composed of highly polydisperse rigid polyhedra. Although non-spherical shapes induce more intricated spatial correlations than spherical shapes because of the multiple contact types (i.e., vertex-face, edge-edge, edge-face, face-face), our numerical data provide evidence that the shear strength independence as the particle size dispersity increases still holds up for assemblies of polyhedra. We explain this finding from compensation mechanisms at the micro-scale between geometrical and mechanical anisotropies developed within the assemblies.


Introduction
Most granular materials in their natural state or in industrial processes are characterized by a very wide range of grain shapes and sizes. Polydispersity can result from various natural processes by which grains are formed (by fragmentation or mineral precipitation). Polydispersity may also be a target in the production, transformation or use of grains. Although the space-filling properties of highly polydisperse systems and their elastic properties under uniaxial compression tests have widely been studied in the past [1][2][3][4], it is remarkable to note that very few studies are devoted to their shear behavior.
In this case, the few existing works have constantly highlighted that the shear strength of dry granular assemblies is independent of the particle size polydispersity. Such surprising result has been evidenced in 2D simulations (via discrete-element modeling) using disks [5][6][7], using regular and irregular polygons [8], and in 3D simulations with spheres [9][10][11]. The present paper is not intended to review on the practical or industrial implications of phenomenon of the shear strength independence on grain size dispersion. It is instead a generalization of such counterintuitive observation to general angular threedimensional grains.
In the following, we present the numerical procedures to build and test the shear strength of samples composed ⇤ e-mail: emilien.azema@umontpellier.fr ⇤⇤ e-mail: david.cantor@polymtl.ca ⇤⇤⇤ e-mail: itthichai_p20@hotmail.com of polyhedra (Sec. 2). In Sec. 3, we focus on the evolution of shear strength and packing fraction with particle size polydispersity. The microstructure, described in terms of force and fabric anisotropies, is discussed in Sec. 4. We conclude, in Sec. 5 with a summary and perspectives of this work.

Numerical procedures
We use the Contact Dynamics (CD) [12,13] to simulate assemblies of angular rigid grains in interaction. The CD method is based on implicit time integration of the equations of motion and a nonsmooth formulation of mutual exclusion and dry friction between particles. The contact forces and body velocities are found by means of an iterative nonlinear Gauss-Seidel algorithm. This method requires no elastic repulsive potential and no smoothing of the Coulomb friction law to determine interaction forces. For this reason, simulations can be performed with larger time steps compared to molecular dynamics simulations. We used LMGC90 [14], a multipurpose software developed in Montpellier, capable of modeling a collection of deformable or rigid particles of various shapes by di↵erent algorithms.
We consider grains of octahedral shape whose size, defined by the diameter d of the circumscribed sphere, is varied in the range of [d min , d max ] with a uniform distribution of particle volume fractions. We define the size span S of In this work, we systematically varied S in the range [0, 0.8] in steps of 0.1. For each value of S , 19000 frictionless grains where geometrically deposited in a cubic box and then isotropically compressed to reach a dense state. Figure 1 shows a screenshots for S = 0.2 and 0.8 after the isotropic compression. Then, we removed the four lateral walls, we set periodic boundary conditions and the friction coe cient between particles is set to 0.4. The shear tests were performed by fixing the bottom wall and applying on the upper wall a pressure P in the vertical direction and a constant velocity v x in the x-direction so the inertial number I is below 10 3 for all the simulations. All the particles in contact either with the upper or the lower wall are bonded and follow the boundary constrains. Gravity was set to 0 to avoid a force gradient within the samples and to prevent particle size segregation. The shearing was undertaken up to reaching the so-called steady flow state in which the volume and shear strength remain constant or present only small fluctuations.

Shear strength and packing fraction
The shear strength is computed from the well known formula of the granular stress tensor , defined as i j = ( P c f c i`c j )/V where, f is the contact force,`is the branch vector (i.e., the vector joining the centers of particles interacting at contacts c) and the sum run over all contacts in a control volume V. Due to the geometry of the simple shear test, the stress state is considered invariant along the y direction. Thus, we may only consider the stresses of the tensor on the shear plane xz.
According to the Mohr-Coulomb model, the e↵ective macroscopic friction coe cient ' is given by sin ' = q/p, where p = ( 1 + 2 )/2 is the mean stress, q = ( 1 2 )/2 is the deviatoric stress, and 1 and 2 are the principal values of the stress tensor on the shear plane [15]. The sample's density is quantified by the packing fraction ⌫ = P p V p /V, where V p is the volume of a particle p. As a general observation (see Fig.2a), sin ' increases rapidly at the beginning of the shearing, and reaches a maximum value before decreasing towards a constant value at large strains (steady state). Along the same line, starting with an initially dense state, all packings dilate during shear and, hence, the volume increases (i.e., ⌫ declines) and tend to a constant value presenting only minor fluctuations. Figure 2 shows sin ' ⇤ and ⌫ ⇤ , the averaged steady state values of sin ' and ⌫, as a function of S . For the set of samples tested, we see that the average shear strength remains almost constant and around sin ' ⇤ ' 0.5. In contrast the packing fraction ⌫ ⇤ increases with S . This finding thus extends to polyhedral particles assembly a similar behavior previously observed in the case of circular/polygonal particles in 2D and spherical particles in 3D.
The fact that the shear strength remains constant while the particle size polydispersity is increased reflects the organization of the microstructure and the features of force transmission that we discuss below.

Micro-scale compensation mechanisms
At the contact scale, the local geometry and the mechanical equilibrium associated with two touching particles can be statistically characterized via the analysis of the angular distributions of the branch vector h`i(⌦) and the contact force h f i(⌦), where ⌦ is the solid angle. For coherence with the previous section, the branch and force vectors are projected on the shear plane 'xy'and can be thus split into their normal contributions h`ni(✓), h f n i(✓), respectively, and tangential contributions h`ti(✓), h f t i(✓), respectively, where ✓ is the orientation of n 0 , the projection of the real normal contact vector on the shear plane. Along with these distributions, we also define the angular distribution of contacts along the shear plane as P(✓). As shown in a number of previous studies, the analysis of these distributions should reveal the key microscopic elements in play as the particle size polydispersity increases. Under simple shear, these distribution take relatively simple shapes. Figure 3 displays polar representations of the five angular distributions P c (✓), h`ni(✓), h`ti(✓), h f n i(✓) and h f t i(✓) as a function of the angle ✓, and averaged in the steady state for three di↵erent values of S . For all S values, we observe an anisotropic behavior in P c (✓) and h f n i(✓), where the peak value occurs, in average, along the major principal stress direction ✓ ' 3⇡/4. Although less pronounced, an anisotropic orientation along the major principal stress is also seen for h`ni(✓). In contrast, the peak values for h`ti(✓) and h f t i(✓) occurs at ✓ f t ' ⇡/2. These distributions are well approximated by their lowest-order Fourier expansions [16,17] 8 > > > > > > > > < > > > > > > > > : where a c is the contact orientation anisotropy, a ln is the normal branch anisotropy, a lt is the tangential branch anisotropy, a f n is the normal force anisotropy, and a f t is tangential force anisotropy. From Fig. 3, we get the corresponding privileged directions ✓ c , ✓ ln , ✓ lt , ✓ f n , and ✓ f t , that are nearly equal and aligned respect to the shear direction. Figure 4 displays the averaged contact and branch anisotropies in the steady state, as function of S . We see that a c declines from 0.55 to 0.4 as S increases. The decrease of a c traduce the fact that the mean connectivity of the particles increases. In fact, a detailed analysis of the connectivity as a function of particle size evidence that large particles are increasingly surrounded by a number of small particles as the size span is increased [9,18]. In same time, the increases of normal branch anisotropies, from 0 to ' 0.09, while the tangential branch anisotropy remains close to 0, evidence that longer branch appears preferentially in the direction of shear. This suggest that biggest particles self-organise along the strong force chains.  The variations of the normal and tangential force anisotropies with S are shown in Fig. 5. Basically, a f n slightly increases with S , from 0.4 to 0.5. In contrast the tangential force anisotropy remains independent with S and close to 0.2. In fact, the above anisotropic parameters are the key elements behind the macroscopic friction coe cient '. Indeed, it can be shown that the general expression of the stress tensor leads to the following simple relation [17,19]: The predicted values of sin ' ⇤ by Eq. (3) are shown in Fig. 2 together with the measured values as a function of S . We see that the approximation of our data by Eq. Figure 5. Evolution of the normal and tangential force anisotropies (a f n , a f t ) as a function of the particle size span S . Bars represent the standard deviation of the data.
(3) is outstanding for all values of S . By virtue of this equation, the independence of sin ' ⇤ with respect to S results from a compensation mechanism between geometrical and mechanical source of anisotropy. More precisely, the decreases of a c is compensated by the increases of a`n + a f n together with the fact that a f t is independent of S . Such compensation mechanism was also observed in previous studies with assemblies of discs, polygons and spheres but was unexpected for polyhedra particles assemblies. Indeed, when dealing with three-dimensional regular grains, a rich contact configuration appears given the di↵erent contact types that can occur: vertex-face, edgeedge, edge-face, face-face and punctually, vertex-vertex. In particular, it is now well known that these contact contribute di↵erently to the shear strength [17]. For example, face-face contacts are more likely to mobilize friction while vertex-face contacts, as more unstable interactions, are supporting relatively lighter stresses. But a surprising observation, detailed in [18], is that, in fact, the contributions of each family of contacts to the shear strength and anisotropies remain independent of the grain size span too.

Conclusion
In this paper, we developed a set of numerical experiments using the three dimensional Contact Dynamics method to generalize -or the opposite -the shear strength independence on particle size using very angular polyhedral grains. We built a set of samples composed of octahedra, where the grain size span varied from mono-size particles, to highly polydisperse particles with a size ratio (maximum/minimum) equivalent to ⇠ 10. The systems were sheared until a steady state, under quasi-static conditions, was reached.
Our numerical data clearly allow us to generalize to three-dimensional angular grains the observation that grain size span does not a↵ect the shear strength although the solid fraction increases as a broader grain size distribution is considered. This unexpected result, consistently evidenced with more simple particle geometry, appear thus to be a robust observation. By means of an additive decomposition of the shear strength in terms of anisotropic parameters related to particle connectivity, force transmission and friction mobilization, we show that subtle compensation mechanisms between anisotropies explains the invariance of shear strength as the particle size polydispersity is increased.
This long series of numerical experiments that include more and more realistic properties (from 2D discs assembly to 3D polyhedra systems), suggest that it is now time to develop controlled and systematic physical experiments to definitively conclude on this phenomenon.