# Particle shape distribution effects on the critical strength of granular materials

Sergio Carrasco, David Cantor, Carlos Ovalle, Frédéric Dubois

► To cite this version:

Sergio Carrasco, David Cantor, Carlos Ovalle, Frédéric Dubois. Particle shape distribution effects on the critical strength of granular materials. Computers and Geotechnics, 2025, 177, pp.106896. 10.1016/j.compgeo.2024.106896. hal-04779354

HAL Id: hal-04779354

https://hal.science/hal-04779354v1

Submitted on 13 Nov 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

Research paper

# Particle shape distribution effects on the critical strength of granular materials

Sergio Carrasco<sup>a,b</sup>, David Cantor<sup>a,b</sup>, Carlos Ovalle<sup>a,b,*</sup>, Frédéric Dubois<sup>c</sup>

## ARTICLE INFO

### Keywords:

- Discrete-element modeling

- Granular materials

- Shear strength

- Particle size distribution

- Particle shape distribution

## ABSTRACT

This study investigates the influence of correlations between particle morphology and gradation on the critical shear strength of three-dimensional granular assemblies via numerical simulations. While grain shape is acknowledged to play a central role in the mechanical behavior of granular media, only a few works have explored the combined effects of grain shapes varying with grain size. Employing three-dimensional discrete-element simulations, we explore the shear behavior of samples with diverse particle size distributions, where grain shapes (ranging from spheres to very angular polyhedra) are assigned based on the relative size of each particle. Using drained triaxial shear tests, we analyze the macroscopic behavior up to large deformation levels. Micro-mechanical analyses are also conducted to understand the underlying mechanisms governing the observed macroscopic behavior, highlighting the role of grain connectivity and force transmission. Surprisingly, well graded samples composed of coarse angular grains and fine rounded ones are not capable of developing higher shear strengths than uniform samples with the same diversity on grain shapes. Conversely, samples where coarser grains are rounded and fines are angular, show constant shear strength as the particle size distribution becomes broader. These findings underscore the importance of simultaneously considering grain size and shape distributions for the assessment of realistic granular soil behavior.

## 1. Introduction

Geotechnical and mining engineering often deal with materials whose mechanical properties are challenging to characterize. This is particularly evident for coarse granular materials, such as crushed rock obtained after blasting and grinding. These materials, composed of grains of varied geological origin and displaying a wide range of sizes, shapes and strengths, are the building blocks for structures such as rockfill dams and mine waste rock piles (Marsal, 1967; Barton and Kjærnsli, 1981; Aubertin et al., 2002; Ovalle et al., 2014). The challenges for the mechanical characterization of coarse grain assemblies lie in the fact that, while laboratory tests can be performed on materials composed of fine gravels, sands or any finer grains, coarse granular materials may contain large cobbles and clasts that are impossible to fit in standard testing devices (Ovalle et al., 2020). Historical attempts to overcome this limitation have consisted of building large laboratory shearing devices (Marachi et al., 1972; Bauer and Zhao, 1993; Verdugo and de La Hoz, 2007; Hu et al., 2011; Linero-Molina et al., 2020), large oedometric cells (Oldecop and Alonso, 2004; Osses et al., 2024) and in situ direct shear boxes (Barton and Kjærnsli, 1981; Matsuoka et al., 2001; Ovalle et al., 2023).

The condition to test a representative mass of material relates the minimal sample dimension with its maximum particle size (Cerato and Lutenegger, 2006; Amirpour Harehdasht et al., 2019; Cantor and Ovalle, 2023; Girumugisha et al., 2024). Even when very large testing methods are used, it is often necessary to reduce the characteristic grain size of rockfill materials by removing the coarsest size fraction to meet standard geotechnical testing requirements (ASTM D7181, 2020). This is usually achieved using downscaling methods, such as parallel grading or scalping, which inevitably alter the particle size distribution (PSD) of the original material. In doing so, other material properties may also change, such as particle shape distribution, if the material presents grain size-shape correlations (Linero et al., 2017; Meng et al., 2023). A number of studies have proven that changes in PSD do not affect the critical shear strength as long as the particle shape remains constant (Muir-Wood and Maeda, 2008; Li et al., 2013; Yang and Luo, 2018; Polanía et al., 2023). However, if the characteristic grain shape is altered by removing coarse particles, the mechanical behavior of the small-scaled sample might not be representative of the original material (Cho et al., 2006; Cavarretta et al., 2010). Scaling methods are

<footnote><sup>a</sup> Department of Civil, Geological and Mining Engineering, Polytechnique Montréal, Québec, Canada</footnote>

<footnote><sup>b</sup> Research Institute of Mining and Environment (RIME), UQAT-Polytechnique, Québec, Canada</footnote>

<footnote><sup>c</sup> LMGC, Université de Montpellier, CNRS, Montpellier, France</footnote>

<footnote>* Corresponding author at: Department of Civil, Geological and Mining Engineering, Polytechnique Montréal, Québec, Canada. E-mail address: carlos.ovalle@polymtl.ca (C. Ovalle).</footnote>
then not appropriate in materials having grain size-shape correlations. This topic, however, remains poorly understood and the coupled effects of varying grain shape as a function of grain size has only been covered by a few works, such as the experimental studies by Yang and Wei (2012) and Ovalle and Dano (2020), or more recently, using two-dimensional discrete-element simulations (Linero et al., 2019; Carrasco et al., 2022, 2023).

The discrete-element method (DEM) (Cundall and Strack, 1979; Moreau, 1988) provides a privileged framework to explore the mechanics of granular media, as it allows for the explicit description of the particle geometry while exploring the multi-scale response of the material. These discrete multi-scale approaches have become key to understanding the complex behavior of disordered granular systems composed of spheres, regular polyhedra, and ellipsoids (Nicot et al., 2005; Azéma et al., 2013; Wu et al., 2021a; Adesina et al., 2023). More realistic grain shapes have also been considered through the use of irregular polyhedra and clusters of spheres (Nie et al., 2019; Xu et al., 2021; Fan et al., 2024), although the number of particles and length of simulations are quickly limited by computational costs. These studies have certainly provided a comprehensive understanding of macro- and micro-mechanical behavior of granular materials, but have neglected the fact that soil often presents grain shape variations across grain sizes.

In this paper, we introduce a systematic DEM study of the effects of the particle shape distribution on the mechanical behavior of granular materials. While previous works have used mono-size mixtures of different grain shapes (bi-shape or mixtures) (Gong et al., 2024; Illana et al., 2024), this work explores a wide range of gradings in samples presenting grain size-shape correlations. Despite the fact that our study motivation focuses on the critical shear strength of coarse granular materials, our outcomes and observations are also valuable to other granular mechanics problems such as determining the runout distance of dry granular flows (Nie et al., 2024) or the analysis of mechanisms behind grain segregation and layering (Makse et al., 1997/03/01).

This manuscript is organized as follows. In Section 2, we introduce the material model and size-shape correlation that allowed us to build a set of 3D numerical samples, from mono-size to highly variable in particle sizes, and containing grain shapes from spheres to very angular polyhedra. In Section 3, we analyze the macro-mechanical behavior of the numerical samples when subjected to triaxial shearing up to relatively large deformations. In Section 4, we use a series of micro-structural descriptors to characterize the effect of size-shape correlations on grain connectivity and force transmission mechanisms. We also identify the contribution each shape class makes to macroscopic shear strength. Finally, we conclude and draw perspectives of this work.

## 2. Numerical model

### 2.1. Dispersity of grain size and shape

The dispersity of grain geometry and properties, also referred as to polydispersity, is common in natural and industrial granular materials (Ovalle and Dano, 2020; Zhao et al., 2023). First, in terms of size variability of grains, the PSD is usually the preferred tool to quantify the size classes present in a material, based on the mass of sample coarser or finer than a certain size. However, the PSD alone does not provide any insight into particle shape properties. The grain shapes can be described by countless factors, being the most common the angularity, elongation, sphericity and roughness (Linero et al., 2017; Angelidakis et al., 2021).

In this work, we consider three-dimensional numerical assemblies with PSD varying from mono-size, uniform samples up to well-graded grain assemblies. The spectrum of particle shapes is varied from regular octahedral to perfect spheres. Fig. 1 displays the six different shapes (P1 to P6) considered in this study. Note that the typical size of grain *d* is related to the circumscribing sphere around the grain.

<caption>Table 1</caption>

Geometrical characteristics of grains P1 to P6. $n_v$ is the number of vertices, $n_f$ is the number of faces, the angularity is the ratio between the volume of particle ($V_v$) and the corresponding circumscribing sphere ($V_c$) (Angelidakis et al., 2021), and $S$ is the sphericity, which denotes the ratio of the surface area of a sphere with equivalent volume to that of the given particle (Wadell, 1933).

<table><thead><tr><td>Particle</td><td>$n_v$</td><td>$n_f$</td><td>Angularity</td><td>Sphericity ($S$)</td></tr></thead><tbody><tr><td>P1</td><td>6</td><td>8</td><td>0.32</td><td>0.84</td></tr><tr><td>P2</td><td>8</td><td>12</td><td>0.37</td><td>0.80</td></tr><tr><td>P3</td><td>12</td><td>20</td><td>0.61</td><td>0.93</td></tr><tr><td>P4</td><td>74</td><td>144</td><td>0.91</td><td>0.98</td></tr><tr><td>P5</td><td>130</td><td>256</td><td>0.95</td><td>0.99</td></tr><tr><td>P6</td><td>Spheres</td><td>—</td><td>1</td><td>1</td></tr></tbody></table>

Table 1 displays a summary of some relevant geometrical properties of each particle, relative to a spherical shape. Note that the shapes selected for this work have different degrees of axial symmetry. For instance, particles P1 (octahedra), P2 (cube) and P3 (icosahedron) are platonic solids and all their faces have the same shape, while polyhedra P4 and P5 are symmetric at least along one of their axis. While we could have chosen other irregular or non-convex polyhedral shapes, there are a few factors motivating the current selection: (1) a gradual increment of vertices or sides of polyhedra is a simple strategy to modify the angularity of the grains, and more complex shapes might have required a set of multiple parameters for their geometrical description; (2) contact detection between polyhedra is not a trivial problem, and venturing into non-convex grain shapes implies using less optimized strategies to contact detection that would increase computational cost. For the selected particle shapes, the elongation index, flatness index, and convexity are equal to 1, leaving sphericity and angularity as the primary shape descriptors that vary among the studied grain shapes (Zhao et al., 2015; Wu et al., 2023).

In terms of the relations between grain size and shape, two correlations were explored:

* Case A: larger grains are rounded while smaller ones are angular.

* Case B: larger grains are angular while smaller ones are rounded.

The choice of these correlations is not arbitrary. Experimental works focused on assessing grain size-shape correlations have observed that both correlations may exist in natural geological materials (Linero et al., 2017; Ovalle and Dano, 2020), due to the configuration of mineral in the rock from which the grains derived. The process of creating correlated size-shape grains can be broken down into the following steps:

* (i) The PSD breadth is defined by selecting a maximum and minimum particle diameter, $d_{max}$ and $d_{min}$, respectively, and the ratio $R_D = d_{max}/d_{min}$ varied from $R_D = 1$ (mono-size) to $R_D = 5.7$.

* (ii) The range $d_{min}$ to $d_{max}$ is divided into six and the particles shapes P1 to P6 are assigned to each sub-size. This produces a certain number of grains per size that ensures that a uniform particle distribution by mass is achieved. In other words, there are more particles in number of the finer class of grains than in the coarser class; see PSDs in Fig. 2.

* (iii) Depending on whether Case A or B is being generated, an appropriate shape is assigned to each specific particle size $d$.

### 2.2. Particle simulation

This numerical study is developed on a discrete-element approach named non-smooth contact dynamics (NSCD). Different to ‘smooth’ DEM, such as that proposed by Cundall and Strack (1979), the NSCD is a discrete method that considers particles as perfectly rigid bodies interacting through mutual exclusion and dry friction (Moreau, 1988; Jean and Moreau, 1992; Dubois et al., 2018). This approach removes the
<figure><img src="image_3.png" /><figcaption>Fig. 1. Set of different particle shapes considered in this work.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Fig. 2. Particle size distribution (PSD) where each particle diameter *d* is normalized by the maximum particle diameter *d*<sub>max</sub>. The different curves correspond to different ratios R<sub>D</sub> = d<sub>max</sub>/d<sub>min</sub>.</figcaption></figure>

necessity for a regularized relationship between particle overlapping and contact force often characterized by contact stiffness and damping coefficients. Instead, the NSCD operates by considering interactions between grains being ruled by Newton's coefficient of restitution and Coulomb's friction law. In effect, the time-stepping strategy in NSCD relies on tracking changes in relative velocities between particles and momentum transfer, which can be integrated to find contact forces and body displacements. Although a complete rigidity of grains cannot be numerically assured, initializing each step with preceding force calculations reduces fluctuations below acceptable levels, generally falling beneath the numerical solution. For more details about the mathematical framework and numerical implementation of the NSCD method, please see Refs. Renouf et al. (2004) and Radjai and Richefeu (2009b).

The NSCD time-stepping is composed of three main ingredients: (i) contact detection, (ii) resolution of contact force and particle relative velocities, and (iii) update of body positions and forces. The contact detection is based on the distance between the particles, and this distance can be computed differently depending on the interaction type. In numerical assemblies composed of bodies of different angularity, from spheres to very angular particles, a variety of interactions can be produced as sphere-sphere, polyhedron-polyhedron, or sphere-polyhedron touch. Note that – still today – DEM simulations very often are developed using spherical grains or clumps of spheres that mimic more complex grain shapes. Real polyhedral grains are less common due to the complexity of interactions, limitations in computational power and management of the data structure whose size can be significant. In this work, we use full polyhedral bodies as well as an implementation that allows us to combine, for the first time in the frame of NSCD, polyhedra with spherical grains. Such a large variety of grain shapes and interactions is rarely seen in DEM simulations.

Eight different contact situations between particles can arise when including polyhedra and spheres in an assembly. Among sphere-sphere contacts we can find (a) point contacts; while between polyhedra we can find (b) vertex-to-face, (c) edge-to-face, (d) face-to-face and (e) edge-to-edge. Finally, between spheres and polyhedra we can find (f) point-to-vertex, (g) point-to-edge and (h) point-to-face. See Fig. 3 for a

schematic representation of each contact type. Numerically, the contact geometry can be represented by a contact point, a contact line or a contact surface. A contact point represents the contact geometry of all sphere-particle, vertex-face, and edge-edge contacts between polyhedra. To model contact transmitted by a line (edge-face contacts) or surface (face-face contacts), it is assumed that the contact is transmitted by a finite number of points; a contact line is defined by two points and a contact surface by a minimum of three points. Although the location for these contact points may vary along the interaction, only their resultant force and point of application have a physical meaning (Cantor et al., 2020).

In terms of contact detection, Fig. 4(a) shows the simplest set of interactions that can occur in our simulations referring to touching spheres, where the distance between two spheres is calculated using the center-to-center distance and respective radii. More challenging is the contact detection between polyhedra, where the distance between a pair of bodies is determined by calculating a separation plane using Cundall’s “common plane” technique (Cundall, 1988), as shown in Figs. 4(b)–4(e). This method iteratively defines the plane that minimizes the distances from the vertices of each polyhedra to the separating plane, based on the perturbation of the orientation of the normal vector of the plane (Dubois, 2011). Finally, Fig. 4(f) shows that the distance between a sphere and a polyhedron is calculated by projecting the sphere center onto the polyhedra surface and determining the distance between the projected point and the sphere’s center. For a more detailed explanation of the contact detection mechanisms used in our simulations, we refer the reader to Radjai and Dubois (2011).

The contact force resolution in NSCD is an iterative process that stabilizes contact force and particle relative velocities fluctuations within a certain tolerance. The implicit nature of the time integrator makes it unconditionally stable. Only the time step ($\Delta t$) needs to be adapted for each simulation or set of boundary conditions, so it is a fraction of the collision time between particles $\tau_c$. In practice, this is achieved by prescribing the relationship $\Delta t \ll \tau_c \approx d/v_{max}$, where *d* is the typical particle size and $v_{max}$ is the maximum velocity a particle can reach at any given time.

## 2.3. Sample generation

Numerical samples for this study were created with a number of particles ($N_p$) that varied depending on the particle size distribution. $N_p$ ranged between 10 000 for $R_D$ = 1 to 25 000 for $R_D$ = 5.7. Assemblies of grains were built by sequentially placing grains in a box of square base W and height D, using geometrical protocols (Jullien et al., 1992; Voivret et al., 2007). The box dimensions and scale followed the laboratory testing standard ASTM D7181, thus, keeping a consistent sample ratio $H/W = 2$ and a sample size that ensured a ratio $W/d_{max} \ge 10$. Fig. 5(a) shows a screenshot of a sample created following such a protocol.

In order to prepare the samples for triaxial shearing, they were previously compressed by applying controlled forces to rigid walls on all faces of the sample box up to an isotropic pressure $P = 10$ kPa. This compression is continued until the void ratio $e = V_v/V_s$ ($V_v$ denotes the volume of voids and $V_s$ the volume of solid grains) is stabilized. Fig. 5(b) to (e) display screenshots of the samples after isotropic compression for the different size-shape correlations and $R_D$.
<figure><img src="image_5.png" /><figcaption>Fig. 3. Contact types between spheres, polyhedra, and spheres in contact with polyhedra.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Fig. 4. Contact detection scheme between spheres, polyhedra, and spheres in contact with polyhedra.</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>Fig. 5. Screenshot of a sample after construction protocol for R<sub>D</sub> = 1 (a); samples with particle size span R<sub>D</sub> = 1.2 (b) and (c), and R<sub>D</sub> = 5.7 (d) and (e) for different size-shape correlations.</figcaption></figure>

During isotropic compression, gravity and friction between particles and walls were set to zero, so no internal pressure gradients are created and the applied stresses on the box correspond to principal stresses.

Drained triaxial shear tests were conducted by applying a constant vertical velocity ($v_s$) to the upper and lower walls, while lateral walls maintain constant confinement pressure $P$; see schematic representation of boundary conditions in Fig. 6. We used rigid boundaries in our simulations to simplify the model and avoid strain localization. However, it is well known that the critical shear strength is not significantly influenced by the boundary conditions (Binesh et al., 2018; Wu et al., 2021b). Upon shearing, the interaction between grains was modeled using a Coulomb's dry frictional coefficient $\mu = 0.4$ and Newton's coefficient of restitution $e_R = 0$. The shearing rate was determined based on the inertial number $I = d\dot{\gamma}\sqrt{\rho/P}$, where $\rho = 2600 \text{ kg/m}^3$ is the particle's density, $\dot{\gamma} = v_s/h_0$ is the shear rate and $h_0$ is the height of the sample after isotropic compression. In order to ensure a quasi-static shearing regime, the inertial number was set to $1 \times 10^{-3}$ (GDR-MiDi, 2004). The axial strain $\epsilon_a = \Delta h/h_0$, where $\Delta h$ is the change in sample height along the direction of shearing, reached 0.4 in all the simulations.

All the simulations were performed on the modeling platform LMGC90 (Dubois et al., 2011, 2023), a free, open-source software capable of simulating discrete mechanical systems using the NSCD method.

### 3. Macroscopic mechanical behavior

The mechanical behavior was studied through measurements for stress, strain and deformation quantities during shearing. First, the volumetric strain was computed using its relationship with the sample's void ratio as $\epsilon_v = (e - e_0)/(1 + e_0)$, where $e_0$ is void ratio after isotropic consolidation. Fig. 7 illustrates the evolution of volumetric strain as a function of the applied vertical strain $\epsilon_a$. At lower levels of axial deformation, samples with low $R_D$ display an initial phase of contraction from $\epsilon_a = 0$ to 0.15, followed by a period of dilatation. All samples exhibited a steady volumetric strain beyond $\epsilon_a \approx 0.3$, except samples with $R_D = 1.2$ and $R_D = 1.5$ of Case A, which nonetheless exhibit a tendency to stabilize. This is a classical limitation of triaxial

<figure><img src="image_8.png" /><figcaption>Fig. 6. Schematic representation of triaxial boundary conditions employing rigid lateral walls around the sample.</figcaption></figure>

shearing on dense samples, which usually require larger strains to fully reach the critical void ratio (Azéma et al., 2013; Cantor et al., 2024).

Secondly, the macroscopic mechanical response was characterized through the normalized shear stress $q/p$, where $q$ is the deviatoric stress and $p$ is the average pressure within the assembly. These quantities are measured using the granular stress tensor $\sigma$, defined as (Potyondy and Cundall, 2004; Andreotti et al., 2013; Nicot et al., 2013):

$$ \sigma_{ij} = \frac{1}{V} \sum_{vc} f_i^c \ell_j^c, \qquad (1) $$

where $f_i^c$ is the i component of the force at contact c, and $\ell_j^c$ is the j component of the branch vector, the vector joining the center of mass of touching particles at contacts c; these definitions are illustrated in Fig. 8. We use the principal stresses of $\sigma$ (i.e., $\sigma_1, \sigma_2, \sigma_3$) to calculate $q = (\sigma_1 - \sigma_3)$ and $p = (\sigma_1 + \sigma_2 + \sigma_3)/3$. The evolution of $q/p$ as a function
<figure><img src="image_9.png" /><figcaption>Fig. 7. Evolution volumetric deformation $\epsilon_v$ for Cases A (a) and B (b) as a function of the axial deformation $\epsilon_a$.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Fig. 8. Schematic representation of two interacting grains, displaying the associated force and branch vectors.</figcaption></figure>

of the axial strain is presented Fig. 9. We can observe a gradual increase in shear resistance in all cases up to $\epsilon_a \approx 0.1$. Beyond this point, and despite some fluctuations, $q/p$ reaches a steady value. The stress–strain behavior lets us suggest that the set of granular samples consistently reach a critical strength after an axial strain $\epsilon_a \approx 0.3$, even if $\epsilon_v$ may still not fully stabilize in some samples.

To compare the size-shape correlation effects on the mechanical behavior, we plot in Fig. 10(a) the evolution of the mean values for $\langle e \rangle$ and $\langle q/p \rangle$ at critical state, as a function of $R_D$. Note that $\langle e \rangle$ and $\langle q/p \rangle$ are computed using the last 10% of axial deformation $\epsilon_a$. For $\langle e \rangle$, Case A exhibits a gradual decrease as function of $R_D$, while Case B displays a parabolic evolution, first increasing up to $R_D \approx 1.5$ and then gradually decreasing for larger grain size spans. This behavior is consistent with the results of previous studies, which showed that the void ratio of granular materials decreases as the $R_D$ increases (Voivret et al., 2007; Azéma et al., 2017). This effect is clearly more pronounced in Case A since the smaller particles have angular shapes, allowing them to fit more closely together.

In terms of the normalized shear stress in the critical state $\langle q/p \rangle$, Fig. 10(b) shows that shear strength is practically independent of $R_D$ for Case A (around $\langle q/p \rangle \approx 1.2$), while Case B exhibits a drop in $\langle q/p \rangle$ as $R_D$ increases. The minor discrepancies between Case A and Case B at $R_D = 1.0$ can be attributed to variations in the random seed used during sample generation and a slight deviation in the $R_D$ value (1.005 instead of 1.0 to avoid crystallization). However, these differences are statistically insignificant, as the standard deviations closely align with the mean values. For the grain dispersion studied, Case A shows a similar behavior to that observed in mono-shape grain assemblies where critical shear strength remains unchanged for differences PSDs (Muir-Wood and Maeda, 2008; Cantor et al., 2020). Nevertheless, Case B shows a counterintuitive behavior, where samples with larger angular grains and fine rounded ones are incapable of

developing any larger strength when the grain size dispersion increases. While the behavior of correlations B were first observed numerically on 2D simulations (Carrasco et al., 2023), it is remarkable that these 3D simulations strongly differ from the 2D studies regarding correlations A.

In order to better understand the mechanisms at the origin of these macroscopic observations, the following sections explore the microstructural and micromechanical characteristics of the different samples during shear steady flow.

# 4. Characterization at the microscale

## 4.1. Connectivity

Particle connectivity refers to the degree to which solid grains are in contact with each other. The first indicator of connectivity we analyze is the number of load-bearing grains ($N_p^*$). Fig. 11 illustrates screenshots of the specimens in the critical state for Cases A and B, with $R_D$ values of 1.5 and 5.7. Load-bearing grains are visualized in gray, while no-force-bearing (floating) particles are displayed in black.

For both correlations A and B, we remark that fewer particles belong to the force-bearing structure as PSD becomes broader as small particles are able to fit in the void spaces between larger ones. Fig. 12(a) shows the evolution of the ratio of force-bearing grains at critical state, with respect to the total number of grains in the sample ($N_p^*/N_p$) as a function of $R_D$. As expected, the ratio $N_p^*/N_p$ decreases as $R_D$ increases. There is only a slight increment of the proportion of grains engaged in force-bearing for Case B up to $R_D \approx 1.5$, presumably due to the capacity of angular large grains to reach more distant particles with their sharp edges. Such gain in connectivity is then lost as soon as the smaller grains start to rattle in the poral space and no longer contribute to force transmission.

Regarding mean number of contact per load-bearing grain or me- chanical coordination number ($Z_m = 2N_c/N_p^*$, where $N_c$ is the number of force-carrying contacts) (Thornton, 2000), Fig. 12(b) presents the evolution of this parameter at critical state as $R_D$ increases. For cor- relation A, $Z_m$ seems to be not affected by the grain size distribution, presenting a value of $Z_m \approx 3.5$. In contrast, in Case B $Z_m$ shows an asymptotic trend, increasing from $Z_m = 3.5$ for $R_D = 1$ to $Z_m = 4.5$ for $R_D = 2$. It is again surprising that Case B presents a larger average number of contacts per grain at critical state, while the shear strength of such samples decreases with $R_D$. This behavior reveals that other microstructural elements, instead of the relatively simple scalar measure $Z_m$, must be acting together to explain the drop of macroscopic shear strength.

## 4.2. Force network

At the contact scale, a local framework ($\hat{n}$, $\hat{i}$) is formed by $\hat{n}$, the unit vector perpendicular to the contact plane between grains, and $\hat{i}$, an
<figure><img src="image_11.png" /><figcaption>Fig. 9. Evolution of the normalized shear stress $q/p$ for size-shape correlations A and B, respectively, as a function of the axial deformation $\epsilon_a$.</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>Fig. 10. Critical state parameters: (a) evolution of mean void ratio $\langle e \rangle$ and (b) normalized mean shear strength in the critical state $\langle q/p \rangle$ as a function of grains size span $R_D$. Error bars display the standard deviation of the data.</figcaption></figure>

orthonormal unit vector oriented along the tangential force direction, as displayed in Fig. 8. In such frame, interparticle forces can be defined as $\hat{f} = f_n\hat{n} + f_t\hat{t}$, where $f_n$ and $f_t$ are the normal and tangential contact forces, respectively.

Fig. 13 illustrates a map of the normal forces at critical state, displayed as columns between the mass centers of touching grains. The thickness of the columns are related to the magnitude of the corresponding normal contact force. Visual inspection reveals the anisotropic and heterogeneous nature of these force networks. In Case A, where smaller particles are angular and larger ones are rounded, the heterogeneity of contact forces appears to be amplified at higher values of $R_D$. However, in Case B there are no major differences in force network architecture as $R_D$ changes (see Figs. 13(b) and 13(d)). This finding, together with the connectivity analysis in the previous section, indicates that a class of fine rounded particles improves connectivity at critical state and contributes to a more uniform force distribution in our granular samples.

Force distribution within a granular assembly is commonly studied through the probability density function (PDF) of contact force magnitudes, generally exhibiting two features: (i) the PDF generally follows a decreasing exponential pattern of forces exceeding the average; (ii) for forces below the mean, the probability of having small-magnitude forces does not decline to zero (Richefeu et al., 2009; Cantor et al., 2018).

Fig. 14 shows the probability density functions of normal forces $P(f_n)$ in log-linear scales for cases A and B at critical state, normalized by the mean normal force $\langle f_n \rangle$. For Case A, we note an important increase in the capability of the grain assembly to transmit larger forces as $R_D$ increases. Samples with a broad particle size distribution can present forces up to three times the maximum forces observed

in a mono-size sample. On the other hand, Case B presents only a slight increase in force magnitude as $R_D$ increases. These figures are in agreement with the previous observations in grain connectivity and force network at critical state.

### 4.3. Fabric and contact and force anisotropies

To understand the drop of critical strength for Case B as a function of $R_D$ and the apparent independence of strength on $R_D$ for case A, it is useful to use the decomposition of the granular stress tensor in terms of micromechanical descriptors, following the works of Oda (1982) and Rothenburg and Bathurst (1989). This approach is based on the fact that contacts, forces and branch vectors have distributions in space that can be described using harmonic functions.

First, the contact anisotropy ($a_c$) represents the excess of contacts oriented in a certain direction and the lack of contacts along the perpendicular direction. This anisotropy is linked to the Fabric tensor ($F_{ij}$), which is defined based on the contact normal vectors $\hat{n}$ as:

$$F_{ij} = \frac{1}{N_c} \sum_{c \in V} n_i^c n_j^c, \qquad (2)$$

so $a_c = 5/6(F_{33} - F_{11})$, where $F_{33}$ is the component of the fabric tensor in the vertical direction (i.e., along axis $z$) and $F_{11}$ is the horizontal component (along axis $x$).

Secondly, the geometric arrangement of grains of different shapes and sizes leads to vector branch lengths that can vary in space and present preferred orientations. Consider the branch vector (i.e., the vector connecting the centers of two touching particles) represented in the contact framework as $\hat{\ell} = \ell_n\hat{n} + \ell_t\hat{t}$, where $\ell_n$ is the branch projection on the normal contact orientation, and $\ell_t$ is the branch
<figure><img src="image_13.png" /><figcaption>Fig. 11. Screenshots of floating particles (black) and bearing-force particles (gray) for Case A (left) and Case B (b) for different grain size dispersion at critical state.</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>Fig. 12. Evolution of (a) ratio of particles engaged in force transmission respect to the total number of particles in the sample and (b) mechanical coordination number as a function of grain size span R<sub>D</sub>, both measured at critical state.</figcaption></figure>

projection on the tangential orientation to the contact. In this case, the anisotropy of branch orientations needs to be considered for its normal and tangential components. In other words, this anisotropy is associated to two different tensors defined as:

$$L_{ij}^{n} = \frac{1}{\langle \ell_n \rangle} \sum_{c \in V} \ell_{n}^{c} n_{i}^{c} n_{j}^{c}, \qquad (3)$$

and

$$L_{ij}^{t} = \frac{1}{\langle \ell_n \rangle} \sum_{c \in V} \ell_{i}^{c} n_{i}^{c} t_{j}^{c}, \qquad (4)$$

where $\langle \ell_n \rangle$ is the average branch length projected onto the normal contact orientation. Using these definitions, the associated anisotropy
<figure><img src="image_15.png" /><figcaption>Fig. 13. Screenshots of force network at critical state for Case A (left) and Case B (b) for different grain size ratios $R_D$. The red connections between the centers of touching grains represent the normal forces whose thickness is proportional to magnitude of normal forces.</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>Fig. 14. Probability density functions of normal contact forces at critical state in log-linear scales for Case A (a) and Case B (b).</figcaption></figure>

in branch length orientations is simply defined as $a_{\epsilon_n} = \frac{5}{6} \frac{L_{11}^n - L_{33}^n}{L_{11}^n + L_{22}^n + L_{33}^n} - a_c$ and $a_{\epsilon_t} = \frac{5}{6} \frac{L_{11} - L_{33}}{L_{11} + L_{22} + L_{33}} - a_c - a_{\epsilon_n}$, where $L_{ij}^n$ corresponds to the ij component of the tensor $L^n$, and $L_{ij}$ corresponds to the ij component of the tensor $L_{ij} = L_{ij}^n + L_{ij}^t$.

Finally, anisotropies in the orientation of contact forces can also be represented using similar tensors, defined as

$$ \xi_{ij}^n = \frac{1}{\langle f_n \rangle} \sum_{c \in V} f_n^c n_i^c n_j^c \quad (5) $$

and

$$ \xi_{ij}^t = \frac{1}{\langle f_n \rangle} \sum_{c \in V} f_t^c n_i^c t_j^c. \quad (6) $$

The corresponding anisotropies of normal and tangential force magnitudes can be found as $a_{f_n} = \frac{\xi_{11}^n - \xi_{33}^n}{6(\xi_{11}^n + \xi_{22}^n + \xi_{33}^n)} - a_c$ and $a_{f_t} = \frac{\xi_{11}^t - \xi_{33}^t}{6(\xi_{11}^t + \xi_{22}^t + \xi_{33}^t)} - a_c - a_{f_n}$, where $\xi_{ij}^n$ corresponds to the ij component of the tensor $\xi^n$ and $\xi_{ij}^t$ corresponds to the ij component of the tensor $\xi_{ij} = \xi_{ij}^n + \xi_{ij}^t$.
<figure><img src="image_17.png" /><figcaption>Fig. 15. Microstructural descriptors at critical state: (a) evolution of geometrical and mechanical anisotropies as a function of the grain size ratio $R_D$; (b) comparison between $\langle q/p \rangle$ computed using the granular stress tensor (Eq. (1)) and its decomposition in microscopic anisotropies (Eq. (7)).</figcaption></figure>

The variables $a_c$, $a_{\ell_n}$, $a_{\ell_t}$, $a_{fn}$ and $a_{ft}$ not only provide indicators for the degree of anisotropy in spacial orientations for each of their corresponding parameters. They also have a deeper connection with the stress state. As shown analytically and validated in the literature (Rothenburg and Bathurst, 1989; Radjai and Richefeu, 2009a; Voivret et al., 2009; Guo and Zhao, 2013; Carrasco et al., 2022), their relationship with $q/p$ reads as follows:

$$ \frac{q}{p} = \frac{6}{5}(a_c + a_{\ell_n} + a_{\ell_t} + a_{fn} + a_{ft}). \qquad (7) $$

This remarkable property allows us to understand the underlying geometrical and mechanical mechanisms at the origin of macroscopic critical shear strength variations discussed before.

Fig. 15(a) shows the evolution of $a_c$, $a_{\ell_n}$, $a_{\ell_t}$, $a_{fn}$ and $a_{ft}$ as a function of $R_D$ for both Cases A and B at critical state. In terms of geometrical anisotropies, we observe that $a_c$ decreases as the grain size dispersion increases, which suggests that adding small particles to the assemblies contributes to creating contacts more homogeneously in all orientations. This phenomenon occurs despite the fact that many of the small particles may rattle the poral space during shearing. In any case, the drop in contact orientation anisotropy is relatively similar in both correlations A and B. For branch length anisotropies, similar trends are observed for both cases A and B, with a slightly increment in their values as the grain size ratio $R_D$ increases. This behavior is simply due to the fact that branch lengths have all the same length for mono-size grains, but, as the grain size dispersion broadens, different branch lengths can be created as particles of different classes interact.

In terms of mechanical anisotropies, the tangential force anisotropy $a_{ft}$ surprisingly remains constant for both cases, with a value of $a_{ft} \simeq 0.15$. This suggests that the mobilization of friction forces is not significantly affected by the size span and varied grain shapes. In contrast, the anisotropy related to the normal contact force $a_{fn}$ shows rather different trends for both cases of size-shape correlations. For Case A, $a_{fn}$ gradually increases as $R_D$ increases, while the variations of this parameter in Case B are substantially smaller. This behavior is linked to the force distributions displayed in Fig. 14, where the force network shows a more heterogeneous spacial distribution in Case A, while the architecture of the networks remains quite similar across different PSDs in Case B.

In Fig. 15(b), the macroscopic critical shear strength is presented using the granular stress tensor (Eq. (1)) and its approximation based on the addition of anisotropies (Eq. (7)). As observed in the figure, the microscopic description of anisotropies exhibit a very good match with

the critical shear strength measures via $\langle q/p \rangle$. Only small differences are observed as the grain size distribution becomes broader, still the values lie within the error bars for the data in the steady flow state.

We can now conclude that, for Case A, the increments in $a_{fn}$ are offset by the drop in $a_c$, which explains why the macroscopic shear strength remains constant with $R_D$. A similar behavior has been observed in mono-shape samples, where $\langle q/p \rangle$ does not depend on PSD because the increase in $a_{fn}$ is matched by the decrease in $a_c$ (Cantor et al., 2018, 2020). For Case B, such a compensation between geometric and mechanic anisotropies only partially occurs for certain grain size ratios $R_D$. Before $R_D \approx 3$, we observe that both the contact orientation and normal force anisotropies gradually decrease as grain size dispersion increases. Only once $R_D \ge 3$, $a_{fn}$ is able to increase and can compensate further reductions in $a_c$. This behavior highlights the underrated role that rounded fine grains were thought to play in polydisperse granular assemblies up to date.

## 5. Conclusions

This research investigated the mechanical behavior of granular assemblies presenting correlations between the size and shape of grains. Two distinct correlations, labeled Case A and Case B, were examined. In Case A, fine grains are more angular than coarse ones, whereas in Case B, fine grains are rounder. The study also considers granular assemblies ranging from mono-sized to polydisperse with a maximum size ratio $d_{max}/d_{min} \simeq 6$. Utilizing the discrete-element method known as non-smooth contact dynamics, three-dimensional granular assemblies were built and tested under triaxial shearing conditions, reaching an accumulated axial deformation of $\epsilon_a = 0.4$.

The investigation of the shear behavior in steady flow state revealed that size-shape correlations significantly influence the mechanical response of granular materials. In Case A, shear strength remains practically unchanged as grain size dispersion increases, whereas Case B exhibits a decline in shear strength with increasing grain size span.

This behavior is accompanied by complex microstructural mechanisms linked to contact organization and force transmission mechanisms. For instance, the loss of macroscopic shear strength in Case B occurs together with particles having more contacts on average. Conversely, Case A showed a rather constant mechanical coordination number as grain size dispersion increased.

Further analysis involved decomposing the granular stress tensor in terms of microstructural anisotropies that take into account the heterogeneous distribution in space of contacts and forces between grains.
This approach has been systematically employed in the literature and was successfully validated in this work. In particular, it was found that the anisotropy linked to contact orientations and normal force magnitudes primarily drive variations in macroscopic shear strength. This detailed microscopical analysis helped us reveal that smaller grains may have a larger role in mechanical strength than previously considered. While for mono-shape granular assemblies the smaller grains are mainly considered to belong to the class of grains that likely to rattle in the poral space, this simplified vision seems no longer sufficient when dealing with assemblies of grains with varied particle shapes. In particular, it seems that small rounded grains have a lubricating effect that restrain assemblies with coarse angular grains (Case B) from developing larger shear strengths. In turn, assemblies of large rounded grains with a class of angular fine grains (Case A) are also limited from developing larger shear strengths, although they do not show any relevant change in macroscopic shear strength as the grain size dispersion is broadened.

The results of this investigation highlight the importance of analyzing both particle size and shape distributions to establish the adequacy of downscaling methods in coarse granular materials. Future studies in the vein of this work should explore more realistic mechanisms that may modify the mechanical behavior of granular matter. In particular, grains with different sizes and shapes are prone to wear and break as the material is sheared. For this, different fragmentation models can be implemented along the grain size-shape correlations. Also, further studies in size-shape correlated materials are needed on the effects of different stress paths and boundary conditions, such as triaxial with flexible walls and direct shear tests. Experimental validation is certainly also critical and highly wanted in order to highlight the relevance of a discrete-element approach. Finally, the poral space is another element to explore in these assemblies where the geometry and size of pores can be strongly affected by variable grain size and shape composition. The analysis of this void space may be relevant for problems related to hydraulic conductivity or internal erosion.

### CRediT authorship contribution statement

**Sergio Carrasco:** Writing – review & editing, Writing – original draft, Visualization, Software, Methodology, Investigation, Formal analysis, Data curation. **David Cantor:** Writing – review & editing, Supervision, Software, Formal analysis, Conceptualization. **Carlos Ovalle:** Writing – review & editing, Supervision, Resources, Funding acquisition, Formal analysis, Conceptualization. **Frédéric Dubois:** Writing – review & editing, Software.

### Declaration of competing interest

The authors declare the following financial interests/personal relationships which may be considered as potential competing interests: Carlos Ovalle reports financial support was provided by Natural Sciences and Engineering Research Council of Canada. Carlos Ovalle reports financial support was provided by Quebec Research Fund Nature and Technology. If there are other authors, they declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

### Acknowledgments

This research work benefited from the financial support of the Natural Sciences and Engineering Research Council of Canada (NSERC) [Ref. RGPIN-2019-06118], the Fonds de recherche du Québec – Nature et technologies (FRQNT), Canada through the ‘Programme de recherche en partenariat sur le développement durable du secteur minier-II’ [Ref. 2020-MN-281267], the industrial partners of the Research Institute on Mines and the Environment (RIME) UQAT-Polytechnique (irme.ca/en). This research was also partially enabled by Calcul Quebec (calcul Quebec.ca) and Compute Canada (computecanada.ca) support.

### Data availability

The data supporting this study's findings are available from the corresponding author upon reasonable request.

### References

Adesina, P., O'Sullivan, C., Wang, T., 2023. DEM study on the effect of particle shape on the shear behaviour of granular materials. Comput. Part. Mech. 1–20.

Amirpour Harehdasht, S., Hussien, M.N., Karray, M., Roubtsova, V., Chekired, M., 2019. Influence of particle size and gradation on shear strength-dilation relation of granular materials. Can. Geotech. J. 56 (2), 208–227.

Andreotti, B., Forterre, Y., Pouliquen, O., 2013. Granular Media: Between Fluid and Solid. CUP.

Angelidakis, V., Nadimi, S., Utili, S., 2021. Shape analyser for particle engineering (SHAPE): Seamless characterisation and simplification of particle morphology from imaging data. Comput. Phys. Comm. 265, 107983.

ASTM D7181, 2020. Method for Consolidated Drained Triaxial Compression Test for Soils. Standard, American Society for Testing and Materials.

Aubertin, M., Bussière, B., Bernier, L., 2002. Environnement et Gestion des Rejets Miniers. Presses Internationales Polytechniques.

Azéma, E., Linero, S., Estrada, N., Lizcano, A., 2017. Shear strength and microstructure of polydisperse packings: The effect of size span and shape of particle size distribution. Phys. Rev. E (ISSN: 24700053) 96.

Azéma, E., Radjai, F., Dubois, F., 2013. Packings of irregular polyhedral particles: Strength, structure, and effects of angularity. Phys. Rev. E 87 (6), 062203.

Barton, N., Kjaernsli, B., 1981. Shear strength of rockfill. J. Geotech. Eng. 107 (7), 873–891.

Bauer, G.E., Zhao, Y., 1993. Shear strength tests for coarse granular backfill and reinforced soils. Geotech. Test. J. 16 (1), 115–121.

Binesh, S.M., Eslami-Feizabad, E., Rahmani, R., 2018. Discrete element modeling of drained triaxial test: flexible and rigid lateral boundaries. Int. J. Civ. Eng. 16, 1463–1474.

Cantor, D., Azéma, E., Ovalle, C., 2024. Failure of an effective stress approach in polydisperse wet granular materials. Phys. Rev. Res. 6 (2), L022008.

Cantor, D., Azéma, E., Preechawuttipong, I., 2020. Microstructural analysis of sheared polydisperse polyhedral grains. Phys. Rev. E 101, 062901.

Cantor, D., Azéma, E., Sornay, P., Radjai, F., 2018. Rheology and structure of polydisperse three-dimensional packings of spheres. Phys. Rev. E 98 (5), 052910.

Cantor, D., Ovalle, C., 2023. Sample size effects on the critical state shear strength of granular materials with varied gradation and the role of column-like local structures. Géotechnique 1–12. http://dx.doi.org/10.1680/jgeot.23.00032.

Carrasco, S., Cantor, D., Ovalle, C., 2022. Effects of particle size-shape correlations on steady shear strength of granular materials: The case of particle elongation. Int. J. Numer. Anal. Methods. Geomech. 46 (5), 979–1000.

Carrasco, S., Cantor, D., Ovalle, C., Quiroz-Rojo, P., 2023. Shear strength of angular granular materials with size and shape polydispersity. Open Geomech. 4, 1–14.

Cavarretta, I., Coop, M., O'sullivan, C., 2010. The influence of particle characteristics on the behaviour of coarse grained soils. Géotechnique 60 (6), 413–423.

Cerato, A.B., Lutenegger, A.J., 2006. Specimen size and scale effects of direct shear box tests of sands. Geotech. Test. J. 29 (6), 507–516.

Cho, G.-C., Dodds, J., Santamarina, J.C., 2006. Particle shape effects on packing density, stiffness, and strength: Natural and crushed sands. J. Geotech. Geoenviron. 132 (5), 591–602.

Cundall, P., 1988. Formulation of a three-dimensional distinct element model—Part I. A scheme to detect and represent contacts in a system composed of many polyhedral blocks. In: Int J Rock Mech Min Sci Geomech Abst. Vol. 25, pp. 107–116.

Cundall, P., Strack, O., 1979. A discrete numerical model for granular assemblies. geotechnique 29 (1), 47–65.

Dubois, F., 2011. Numerical modeling of granular media made of polyhedral particles. In: Discrete-Element Modeling of Granular Materials. ISTE - Wiley, pp. 233–262, URL https://hal.science/hal-00807348.

Dubois, F., Acary, V., Jean, M., 2018. The contact dynamics method: A nonsmooth story. C. R. Méc. 346 (3), 247–262.

Dubois, F., Jean, M., Renouf, M., Mozul, R., Martin, A., Bagnéris, M., 2011. LMGC90. In: 10e colloque national en calcul des structures. p. 8 p, 10e colloque national en calcul des structures.

Dubois, F., Jean, M., et al., 2023. LMGC90 wiki page. https://git-xen.lmgc.univ-montp2.fr/lmgc90/lmgc90_user/wikis/home. (Online Accessed 2 July 2024).

Fan, M., Su, D., Chen, X., 2024. Framework for incorporating multi-level morphology of particles in DEM simulations: independent control of polydisperse distributions of roundness and roughness while preserving form distributions in granular materials. Acta Geotech. 1–24.

GDR-MiDi, 2004. On dense granular flows. Eur. Phys. J. E 14, 341–365.

Girumugisha, G., Ovalle, C., Ouellet, S., 2024. Grading scalping and sample size effects on critical shear strength of mine waste rock through laboratory and in-situ testing. Int. J. Rock Mech. Min. Sci. 183, 105915.
Gong, J., Pang, X., Tang, Y., Yang, Z., Jiang, J., Ou, X., 2024. Effects of angularity and content of coarse particles on the mechanical behaviour of granular mixtures: a DEM study. Granul. Matter 26, 17.

Guo, N., Zhao, J., 2013. The signature of shear-induced anisotropy in granular media. Comput. Geotech. 47, 1–15.

Hu, W., Dano, C., Hicher, P.-Y., Le Touzo, J.-Y., Derkx, F., Merliot, E., 2011. Effect of sample size on the behavior of granular materials. Geotech. Test. J. 34 (3), 186–197.

Illana, E., Qyteti, K., Scharnowski, M., Brömmer, M., Wirtz, S., Scherer, V., 2024. Shape-changing particles for locally resolved particle geometry in DEM simulations. Particuology 89, 185–190.

Jean, M., Moreau, J.J., 1992. Unilaterality and dry friction in the dynamics of rigid body collections. In: 1st Contact Mechanics International Symposium. pp. 31–48, 1st Contact Mechanics International Symposium.

Jullien, R., Pavlovitch, A., Meakin, P., 1992. Random packings of spheres built with sequential models. J. Phys. A Math. Theor. 25 (15), 4103.

Li, G., Ovalle, C., Dano, C., Hicher, P.-Y., 2013. Influence of grain size distribution on critical state of granular materials. In: Constitutive Modeling of Geomaterials. Springer, pp. 207–210.

Linero, S., Azéma, E., Estrada, N., Fityus, S., Simmons, J., Lizcano, A., 2019. Impact of grading on steady-state strength. Geotech. Lett. 9 (4), 328–333.

Linero, S., Fityus, S., Simmons, J., Lizcano, A., Cassidy, J., 2017. Trends in the evolution of particle morphology with size in colluvial deposits overlying channel iron deposits. In: EPJ Web of Conferences. Vol. 140, EPJ Web of Conferences, p. 14005.

Linero-Molina, S., Bradfield, L., Fityus, S.G., Simmons, J.V., Lizcano, A., 2020. Design of a 720-mm square direct shear box and investigation of the impact of boundary conditions on large-scale measured strength. Geotech. Test. J. 43, 1463.

Makse, H.A., Havlin, S., King, P.R., Stanley, H.E., 1997/03/01. Spontaneous stratification in granular mixtures. Nature 386 (6623), 379–382. http://dx.doi.org/10.1038/386379a0,

Marachi, N., Chan, C., Seed, H., 1972. Evaluation of properties of rockfill materials. J. Soil Mech. Found. Div. 98 (1), 95–114.

Marsal, R., 1967. Large scale testing of rockfill materials. J. Soil Mech. Found. Div. 93 (2), 27–43.

Matsuoka, H., Liu, S., Sun, D., Nishikata, U., 2001. Development of a new in-situ direct shear test. Geotech. Test. J. 24 (1), 92–102.

Meng, M., Duan, X., Shi, J., Jiang, X., Cheng, L., Fan, H., 2023. Influence of particle gradation and morphology on the deformation and crushing properties of coarse-grained soils under impact loading. Acta Geotech. 18 (11), 5701–5719.

Moreau, J.J., 1988. Unilateral contact and dry friction in finite freedom dynamics. In: Nonsmooth Mechanics and Applications. Springer, pp. 1–82.

Muir-Wood, D., Maeda, K., 2008. Changing grading of soil: Effect on critical states. Acta Geotech. 3, 3–14.

Nicot, F., Darve, F., Group, R., 2005. A multi-scale approach to granular materials. Mech. Mater. 37 (9), 980–1006.

Nicot, F., Hadda, N., Guessasma, M., Fortin, J., Millet, O., 2013. On the definition of the stress tensor in granular media. Int. J. Solids Struct. 50 (14–15), 2508–2517.

Nie, J.-Y., Cui, Y., Wu, Z., Zhang, L., Fang, J., 2024. DEM study on role of fines in mobility of dry granular flows considering particle size-shape correlation. Comput. Geotech. 166, 105980.

Nie, Z., Zhu, Y., Wang, X., Gong, J., 2019. Investigating the effects of Fourier-based particle shape on the shear behaviors of rockfill material via DEM. Granul. Matter 21, 1–15.

Oda, M., 1982. Fabric tensor for discontinuous geological materials. Soils Found. 22 (4), 96–108.

Oldecop, L.A., Alonso, E.E., 2004. Testing rockfill under relative humidity control. Geotech. Test. J. 27 (3), 269–278.

Osses, R., Pineda, J., Ovalle, C., Linero, S., Sáez, E., 2024. Scale and suction effects on compressibility and time-dependent deformation of mine waste rock material. Eng. Geol. 340, 107668.

Ovalle, C., Dano, C., 2020. Effects of particle size-strength and size-shape correlations on parallel grading scaling. Geotech. Lett. 10 (2), 191–197.

Ovalle, C., Frossard, E., Dano, C., Hu, W., Maiolino, S., Hicher, P.Y., 2014. The effect of size on the strength of coarse rock aggregates and large rockfill samples through experimental data. Acta Mech. 225, 2199–2216.

Ovalle, C., Girumugisha, G., Cantor, D., Ouellet, S., 2023. Size effects assessment of mine waste-rock shear strength combining numerical, laboratory and in situ approaches. In: SSIM 2023: Third International Slope Stability in Mining Conference. Australian Centre for Geomechanics, pp. 291–300.

Ovalle, C., Linero, S., Dano, C., Bard, E., Hicher, P.-Y., Osses, R., 2020. Data compilation from large drained compression triaxial tests on coarse crushable rockfill materials. J. Geotech. Geoenviron. (ISSN: 1090-0241) 146, 06020013.

Polanfa, O., Cabrera, M., Renouf, M., Azéma, E., Estrada, N., 2023. Grain size distribution does not affect the residual shear strength of granular materials: An experimental proof. Phys. Rev. E 107 (5), L052901.

Potyondy, D.O., Cundall, P.A., 2004. A bonded-particle model for rock. Int. J. Rock Mech. Min. Sci. 41 (8), 1329–1364.

Radjai, F., Dubois, F., 2011. Discrete-Element Modeling of Granular Materials. Wiley-Iste.

Radjai, F., Richefeu, V., 2009a. Bond anisotropy and cohesion of wet granular materials. Phil. Trans. R. Soc. A 367 (1909), 5123–5138.

Radjai, F., Richefeu, V., 2009b. Contact dynamics as a nonsmooth discrete element method. Mech. Mater. (ISSN: 0167-6636) 41 (6), 715–728, Advances in the Dynamics of Granular Materials.

Renouf, M., Dubois, F., Alart, P., 2004. A parallel version of the non smooth contact dynamics algorithm applied to the simulation of granular media. J. Comput. Appl. Math. 168 (1–2), 375–382.

Richefeu, V., El Youssoufi, M.S., Azéma, E., Radjai, F., 2009. Force transmission in dry and wet granular media. Powder Technol. 190 (1–2), 258–263.

Rothenburg, L., Bathurst, R.J., 1989. Analytical study of induced anisotropy in idealized granular materials. Geotechnique 39 (4), 601–614.

Thornton, C., 2000. Numerical simulations of deviatoric shear deformation of granular media. Géotechnique 50 (1), 43–53.

Verdugo, R., de La Hoz, K., 2007. Strength and stiffness of coarse granular soils. In: Soil Stress-Strain Behavior: Measurement, Modeling and Analysis. Springer, pp. 243–252.

Voivret, C., Radjai, F., Delenne, J.-Y., El Youssoufi, M.S., 2007. Space-filling properties of polydisperse granular media. Phys. Rev. E 76 (2), 021301.

Voivret, C., Radjai, F., Delenne, J.Y., El Youssoufi, M.S., 2009. Multiscale force networks in highly polydisperse granular media. Phys. Rev. Lett. 102, 178001.

Wadell, H., 1933. Sphericity and roundness of rock particles. J. Geol. 41 (3), 310–331.

Wu, M., Lu, J., Li, X., Pan, S., Wang, J., Yin, Z., 2023. Quantifying the morphology of crushed sand particles using X-ray micro-tomography. Granul. Matter 25 (4), 80.

Wu, K., Sun, W., Liu, S., Cai, G., 2021a. Influence of particle shape on the shear behavior of superellipsoids by discrete element method in 3D. Adv. Powder Technol. 32 (11), 4017–4029.

Wu, M., Wang, J., Russell, A., Cheng, Z., 2021b. DEM modelling of mini-triaxial test based on one-to-one mapping of sand particles. Géotechnique 71 (8), 714–727.

Xu, M.Q., Guo, N., Yang, Z.X., 2021. Particle shape effects on the shear behaviors of granular assemblies: irregularity and elongation. Granul. Matter 23 (2), 25.

Yang, J., Luo, X.D., 2018. The critical state friction angle of granular materials: does it depend on grading? Acta Geotech. 13, 535–547.

Yang, J., Wei, L., 2012. Collapse of loose sand with the addition of fines: the role of particle shape. Géotechnique 62 (12), 1111–1125.

Zhao, B., Wang, J., Coop, M.R., Viggiani, G., Jiang, M., 2015. An investigation of single sand particle fracture using X-ray micro-tomography. Géotechnique 65 (8), 625–641.

Zhao, J., Zhao, S., Luding, S., 2023. The role of particle shape in computational modelling of granular matter. Nat. Rev. Phys. 5 (9), 505–525.