The influence of fines content and size-ratio on the micro-scale properties of dense bimodal materials

This paper considers factors inﬂuencing the fabric of bimodal or gap-graded soils. Discrete element method simulations were carried out in which the volumetric ﬁnes content and the size ratio between coarse and ﬁne particles were systematically varied. Frictionless particles were used during isotropic compression to create dense samples; the coefﬁcient of friction was then set to match that of spherical glass beads. The particle-scale data generated in the simulations revealed key size ratios and ﬁnes contents at which transitions in soil fabric occur. These transitions are identiﬁed from changes in the contact distributions and stress-transfer characteristics of the soils and by changes in the size of the void space between the coarse particles. The results are broadly in agreement with available experimental data on minimum void ratio and contact distributions. The results have implications for engineering applications includ-ing assessment of the internal stability of gap-graded soils in embankment dams and ﬂood embankments.


Introduction
Most soil mechanics theory for cohesionless soils has been developed from experiments using uniformly graded or "clean" sands. However, the presence of cohesionless fines (resulting in a gap-graded material) is known to influence soil properties. Prior research on this topic has typically empirically related the volumetric proportion of fines in the material (F fine ) to the macro-scale, overall, soil properties using physical tests (e.g. [1]). The current study uses discrete element modelling [2] to extend consideration to the particle-scale. A bimodal material is considered as this is the simplest type of gap-graded material and the interplay between the fine and coarse particles can be isolated from any effects of the gradation of sizes within either the fine or coarse fraction.
Based on experimental observation of soils at low cohesionless fines contents, the contribution of the fine particles can be neglected and therefore the void ratio of the coarse particles only, termed the granular or skeleton void ratio, e sk , can be used in place of void ratio, e, as a state variable [1]: Thevanayagam et al. [3] extended this concept to account for experimental evidence that at low fines contents the fine particles make a reduced but finite contribution towards shear strength. Skempton and Brogan [4] proposed that where fines do not fill the voids which form between the coarse particles in gapgraded materials, these loose fine particles carry a reduced effective stress. They proposed a stress-reduction factor, α, defined as the ratio of effective stress in the fines to the overall effective stress: where σ fine is the effective stress transferred by the fines and σ is effective overburden stress. In experiments α can be inferred from permeameter testing in which the hydraulic gradient at the onset of piping failure (i c ) is observed. Then α is estimated by combining Eq. (1) with the hydraulic gradient obtained from Terzaghi's theory for heave at zero effective stress (i c,heave ): when α = 1, the stress is shared equally between the coarse and fine fractions. Li and Fannin [5] derived a generalized form of Eq. 2 for the case where an external stress is applied. Shire et al. [6] compared α values calculated directly for virtual samples created using discrete element method (DEM) simulations with experimental data for equivalent physical samples; they concluded that Eq. 3 gives reasonable estimates of α and that α depends on F fine and packing density. Skempton and Brogan [4] identified two key fines contents: (i) S*, the critical fines content at which the fines just fill the voids between the coarse particles, and below which α < 1. S* was estimated to lie between F fine = 24 and 29 % for dense and loose samples respectively. (ii) S max , the fines content at which the fines separate the coarse particles from one another, which should be no higher than F fine = 35 % [4].
when F fine < S* the soil has an "underfilled" fabric, when S * < F fine < S max , the fabric is "filled" and F fine > S max gives an "overfilled" fabric, as shown schematically in Fig. 1. The concepts of α and S* have implications for engineering practice and the current research was carried out to complement a broader study of the internal stability of cohesionless soils [7]. Internal stability describes the ability of the coarse fraction of a soil to prevent the erosion of the fines under seepage [8]. Kenney and Lau [9] defined three prerequisites for internal instability: (i) a primary matrix of coarse particles which transfers stresses; (ii) loose finer par-ticles in the voids between the primary matrix, which do not carry effective stress and can be moved by seepage; (iii) the inter-void constrictions within the primary matrix must be large enough to allow the loose finer particles to be transported from void to void by seepage. Skempton and Brogan [4] found that when F fine < S* internal instability can initiate at lower hydraulic gradients than would be expected to cause failure by heave (i.e. α < 1). The effect of F fine on internal stability has also been recently experimentally demonstrated by Sibille et al. [10]. In their DEM study Shire et al. [6] established a link between α and the empirical Kézdi criterion for internal instability [11], which is based upon a size-ratio between the coarse and fine particles. Shire and O'Sullivan [12] also showed that there is a link between the Kézdi criterion and micro-scale parameters such as coordination number and To et al. [13] used DEM to study the effect of packing arrangement on the characteristics of the primary fabric. Other researchers have shown that, for a given void ratio, cohesionless fines contribute less per unit volume to shear strength, stiffness and liquefaction resistance than coarse particles [1,3,14,15].
This paper considers the effect of varying F fine and the size-ratio between coarse and fine particles, χ = D coarse /D fine , on the micro-scale properties of dense, homogeneous, isotropic collections of bimodal spheres using DEM. The paper examines the relationship between F fine , χ and both the void ratio and sizes of voids between the coarse particles, following which the contact and stress distributions within the samples are analysed. The results are verified using available experimental data for bimodal materials.

Simulation approach and samples analysed
The DEM simulations were carried out using a modified version of the open-source code Granular LAMMPS [16]. The samples were created by generating a cloud of non-contacting frictionless spheres at random locations within a periodic cell, which avoids boundary effects associated with walls. A  Hertz-Mindlin contact model was used and the simulation input parameters are presented in Table 1. The samples were compressed isotropically and monotonically to a mean normal stress of p = 50 kPa. A servo-controlled algorithm was used to adjust the strain rate until this target stress was reached [17]. This approach generates samples with the densest possible packing at this stress level. These samples effectively have a relative density of 100 %, enabling comparison with the experimental data on dense bimodal assemblies [18,19], as discussed below. Following isotropic compression the coefficient of friction was set to μ = 0.3, which is approximately equal to the experimental value reported for spherical glass beads [20]. To ensure a stable state was achieved the simulations were terminated once p and the coordination number (the average number of contacts per particle in the system), Z, remained unchanged for 20,000 timesteps. All the results presented here correspond to this end state. All the simulations were carried out in a gravity-free environment to allow the use of periodic boundaries, thus removing boundary effects [21], and allowing easy identification of those particles which participate in effective stress transfer through the sample. Further details of the simulation methodology are given in Shire [7].
In selecting the samples for analysis, consideration was given to the range of F fine and χ values that merited consideration. Following the contribution of Skempton and Brogan, discussed above, the research specifically considered F fine values about the critical fines content, i.e., F fine = 20, 25, 30 and 35 %. The definitions of α and S* are based on an assumption the fines can fit within the voids formed between the coarse particles. However, there is a limit to the size ratio, χ = D coarse /D fine , at which a fine particle (with diameter D fine ) can fit between coarse particles (with diameter D coarse ), and below this limiting ratio granular materials cannot be considered to be gap-graded. Based on a consideration of mutually touching uniform circles, Lade et al. [18] suggested that at χ ≈ 6.5 a single fine particle can fit within the smallest possible constriction formed between three coarse particles. Ratios of χ ≈ 6 have been adopted in the definition of "gap-graded" materials in studies considering the effect of non-plastic fines on soil behaviour (e.g. [3,22,23]). However, the smallest circle which can be inscribed between four mutually contacting circles occurs at χ ≈ 2.4, and this was taken as a lower limit χ for gap-graded soils by Choo and Burns [24]. Both these limits are based on two-dimensional considerations of inter-void constrictions. In three dimensions the largest sphere which can fit within the void body of the densest face centred cubic packing of uniform spheres (e = 0.35) occurs at χ ≈ 4.45 and for the looser orthorhombic packing (e = 0.65) this occurs at χ ≈ 2. It is clear from these theoretical considerations that for χ ≤ 2 the material cannot be considered gap-graded. Between χ = 2 and χ = 6 there could exist an intermediate type of gap-graded behaviour. This has been shown experimentally for the variation of void ratio by Yerazunis et al. [25]. For χ > 6, when the size ratio increases, the local increase in void ratio of the fines close to coarse particle surfaces becomes less significant [26] and the fine packings between the coarse particles will progressively densify so that with increasing χ the critical fines content, S*, increases. Consideration of these earlier results motivated use of χ values of 2, 4, 6, 8 and 10 in the analyses.
The 22 bimodal or binary samples analysed are summarised in Table 2. In each sample there were at least 100 coarse particles, and the total number of particles in the simulations ranged from 307 to 54033. One simulation was carried out with F fine = 50 % and χ = 4 to allow comparison with the experimental data of Pinson et al. [19], and one simulation was carried out with F fine = 0, i.e. with a monodisperse (uniform) size distribution in order to provide baseline data. A parametric study was carried out in order to show that a representative element volume (REV) had been obtained, meaning there was no statistically-significant change in the results (e.g. in the connectivity distribution) with an increase in the number of particles within the periodic cell [7].

Calculation of α
To calculate α from the DEM data, the average stress tensor within a particle,σ p i j , was determined using the approach described in [27]; the mean particle stress is then For the isotropic stress state considered here, α is defined in terms of the mean (effective) normal stress, p , for the whole sample: where V p is particle volume, V is the sample volume and N p is the number of particles. The mean normal stress for the fines is: where: n: sample porosity; N p,fine : number of fine particles and α is: 3 Results

Analysis of void space
As noted above, previous researchers have carefully considered how to calculate e in gap-graded materials. The relationship between F fine and e is presented for each χ value in Fig. 2. For the F fine values presented here there is little overall variation of e with F fine for a given size ratio. For each χ a minimum void ratio can be identified (e min ) and the F fine at which e min is obtained decreases with increasing χ , in agreement with experimental observations [18].
McGeary [28] carried out experiments on bimodal mixtures of glass beads with F fine varying between 10 and 50 % to investigate the maximum densities that could be obtained; the data from these experiments was considered by Lade et al. [18]. The e min values for each χ value observed in Fig. 2 are compared with the experimental data in Fig. 3. There is good agreement between the simulations and experiments, although the DEM data give a slightly denser packing at a given χ . This observed difference is to be expected as the DEM particles are perfectly spherical and frictionless (during the isotropic compression stage); physical glass ballotini deviate from this ideal [29]. As χ increases, the rate at which e min decreases reduces. Lade et al. [18] concluded the relationship had a bi-linear shape, with a distinct change in the gradient of the experimental data at χ ≈ 6.5. The DEM data show a reduction in the gradient of the e min − χ plot with increasing χ , but the relationship does not replicate the distinct bi-linear shape attributed to the experimental data.
Lade et al. [18] hypothesized that the change in gradient of the e min − χ plot at χ ≈ 6.5 occurs because the fines pack more efficiently between the coarse particles when the size ratio increases, as illustrated schematically in Fig. 1. The void sizes between the coarse particles for the bimodal materials with F fine > 0 were determined by considering only the coarse particles in void partitioning. Voids were identified using the Delaunay method proposed by Reboul et al. [30], which considers a void to be defined by the tetrahedra formed by a Delaunay tessellation of particle centers, as shown schematically in Fig. 4. The diameter (size) of a void, D void , is defined by the largest sphere which can be inscribed between the particles forming the tetrahedra. The resultant coarse particle void size distributions (VSDs) for four bimodal samples are given on Fig. 5. Figure 5 also includes the VSD for the monodisperse sample and these data agree with that of Bryant et al. [31] for a similar analysis of monodisperse spheres. The binary samples contain only 100 coarse particles and therefore have a less smooth distribution than the monodisperse sample, which contains 600 particles. Analytically known values of D void for regular packings of monodisperse spheres (close-packed cubical/hexagonal, body-centred cubical/tetragonal and orthorhombic) are also included for reference.
The smallest voids in the monodisperse random packing are D void = 0.2245 D coarse , which is equal to the minimum void between the densest possible regular packings (closepacked cubical and hexagonal). However, the majority of the voids are larger than this, with most (∼70 %) falling within the range of D void = 0.3−0.5 D coarse . 15 % of the voids have D void > 0.5 D coarse . When χ = 2, D fine = 0.5D coarse , and therefore D fine > D void for the majority of the voids, meaning that the fines will not be able to sit between the coarse particles under reduced stress, confirming that materials with χ = 2 should not be considered to be gap-graded. As χ increases, D fine < D void meaning single fines and groups of fines are able to fit more efficiently within voids. When the gap-ratio is large (χ = 10) and F fine ≈ S * (F fine = 25 %) the void size distribution is similar to the sample containing no fines indicating that the coarse particles form a dense network very similar to that if there were no fines present.

Contact density
The extent to which the finer particles carry a reduced stress, i.e. the α value, is influenced by the contact network, and the connectivity (i.e. number of contacts per particle) of the finer particles. Pinson et al. [19] identified contacts between coarse and fine particles in bimodal packings of spheres with χ = 2 and 4 using a liquid bridge technique. The resultant connectivity data can be compared with the DEM data generated in this study. The DEM samples are somewhat denser than the experimental samples. However, while the experimental void ratio was measured for the whole sample, connectivity was measured away from the sides of the container in order to avoid wall effects and therefore void ratio is probably overestimated in the experiments.
The distributions of connectivity for χ = 2 and 4 are given in Fig. 6. Figure 6a, b give the connectivity for fine to fine (C fine−fine ) and fine to coarse (C fine−coarse ) contacts respectively for χ = 2 and F fine = 25, 30 and 35 % for the DEM simulations. Equivalent data for χ = 4 are presented in Fig. 6c, d. For both χ values experimental data are included in the figures; for χ = 2 the experimental data considers F fine = 28 %, while experimental data for F fine = 28 and 50 % were available for χ = 4. In all cases the experimental and DEM data show the same upper limits to the distribution of connectivities and the proportions of particles with 0 contacts are broadly similar for equivalent F fine values. For χ = 2 the experimental C fine−fine distribution with F fine = 28 % is similar to the DEM dis- Experimental data, F fine = 28% Experimental data, F fine = 50% Fig. 6 Cumulative distribution of connectivity for size ratio: χ = 2: a fine to fine contacts; b fine to coarse contacts. χ = 4: c fine to fine contacts; d fine to coarse contacts. Experimental data taken from Pinson et al. [19] tributions for F fine = 30 and 35 % (Fig. 6a), although the experimental data has a greater proportion of particles with C fine−fine > 4 and fewer with C fine−fine = 0. The DEM sample with F fine = 25 % shows fewer fine to fine contacts per particle, specifically there are many more particles with C fine−fine = 0 in this sample. The particles with C fine−fine = 0 are likely to be either trapped between two coarse particles or isolated within the voids between the coarse particles. As shown in Fig. 6b, the experimental and DEM distributions of C fine−coarse show good agreement for all three DEM samples despite the difference in void ratio.
As shown in Fig. 6c, for χ = 4, C fine−fine increases as F fine attains and then exceeds the critical fines content at which the fines fill the voids. For F fine = 50 %, the C fine−fine distributions are very similar for both experimental and DEM data. While the DEM data for F fine = 25 % and F fine = 30 % show far fewer fine to fine contacts per particle than the experimental data for F fine = 28 %, there is a close agreement between the experimental data for F fine = 28 % and the DEM data for F fine = 35 %.
As well as looking at the connectivity data, it is useful to consider the overall coordination number, given as : where N c is the number of interparticle contacts in the system and N p is the number of particles in the system. Figure 7a and Table 2 show the variation of Z with χ and F fine for each of the samples. For F fine = 20 % there is a steep reduction from Z = 5.04 at χ = 2 to Z = 0.17 at χ = 6. In the DEM data Z < 1 is possible as gravity has been neglected. In a bimodal material the number of fine particles far exceeds the number of coarse particles, even when accounting for only a fraction of the total volume (when χ = 10, one coarse particle has the same volume as 1000 fine particles). This means that the coordination number is largely determined by contacts involving the fines. The reduction of Z between χ = 2 and χ = 6 for the F fine = 20 % material is evidence that the fines transition from being well connected to being predominantly unconnected and non-stress transmitting. This supports the argument of Lade et al. [18] that as χ increases, fines are better able to fit within voids and so play a reduced role in stress transfer. For F fine ≥ 25 % there is little variation of Z with χ . This is because the fines completely fill the voids between the coarse particles and so have many interparticle contacts, regardless of F fine or χ . With reference to Skempton and Brogan [4], for materials at their highest relative density, the critical fines content S* at which fines just fill the voids between coarse particles occurs at F fine ≈ 24 % and can be identified by an increase in Z with F fine . Table 2 show the variation of the stressreduction factor, α, with χ and F fine . In samples with F fine ≥ 30 %, α ≈ 1, indicating that the coarse and fine particles contribute approximately equally to stress transfer. The same is true of samples with χ = 2 regardless of F fine , as the fines are unable to completely fit within the voids. For samples with F fine = 20 % and χ ≥ 6, α ≈ 0 and the fines are completely loose within the voids and play almost no role in stress transfer, hence F fine < S * .

Figure 7b and
For the samples at F fine = 25 %, which is just above the critical fines content S*, α reduces steeply with increasing χ . Comparing Fig. 7a, b, it can be seen that for F fine = 25 %, α falls rapidly with χ , whereas Z does not. This supports the hypothesis that although fines are not carrying the same proportion of stress as the coarse particles, they are still in contact with the coarse particles and therefore forming weak lateral force chains to support the main strong force chains [6,32]. Although these fines are under low stress, they still perform an important supporting role to the more highly stressed particles and their removal, for example due to internal instability under seepage, could lead to the collapse of the stress-transferring matrix [6,32]. Figure 8a shows the relationship between the stressreduction factor, α, and the coarse to coarse coordination number, Z coarse−coarse [33]: where N c,coarse−coarse is the number of contacts between coarse particles, and N p,coarse is the number of coarse particles. The value Z coarse−coarse < 4 has been highlighted on Figure 8a, as when Z coarse−coarse < 4 coarse particles cannot be considered to be forming a mechanically stable matrix on    [34] their own and must therefore be separated from one another by fines (i.e. be overfilled). Figure 8a shows that the fine and coarse particles carry approximately equal stress (i.e. α ≈ 1) in every sample with Z coarse−coarse < 4, confirming the hypothesis that overfilled samples must be internally stable. The stress-reduction α values drop rapidly as Z coarse−coarse increases beyond 4 and the coarse particles are able to form a stress-transmitting primary matrix, leaving fines transmitting little stress. Figure 8b, considers the relationship between α and the skeleton void ratio, e sk (Eq. 1). When e sk ≥ e coarse,max , the experimental maximum void ratio for monodispersed spheres alone [34], the coarse particles must be separated from one another by fine particles and α ≈ 1. The e sk = e coarse,max condition was termed the limit void ratio by Salgado et al. [1]. who found a distinct change in the stress-strain response of silty sands when this was exceeded. Considering effective stress transfer, when e sk < e coarse,max the coarse particles must be in contact with one another and therefore dominate stress transfer, as shown by α < 1. This effect becomes more prominent as e sk reduces further below e coarse,max . Considering both Fig. 8a and b, it is clear that even at low fines contents (20 %) samples with χ = 2 cannot form a stable fabric comprised of coarse particles alone, whereas samples with χ = 4 form an intermediate fabric in which the coarse particles transfer more stress than fine.
As noted above the void size distributions (Fig. 5) suggest that for large χ values and when F fine ≈ S * (i.e. χ = 10 and F fine = 25 %), the coarse matrix is similar to a monodisperse material. This is confirmed by the very low stress in the fines (α = 0.09). The high coordination number (Z = 5.35), suggests that this is just at the point where the fines fill the voids. As F fine increases to 30 % the void diameters between the coarse particles increase noticeably as fines begin to separate them and the stress in the fines also increases to α = 1.02. For the samples with χ = 2 the coarse voids are significantly larger because, as discussed in Sect. 3.2, the fines are larger than the majority of the voids in the monodisperse sample. For these samples α ≈ 1.
Thornton and Antony [17] showed that columns or chains of "strong" contact forces (i.e. forces of above average magnitude) transfer almost all the deviatoric load through granular materials. It is therefore reasonable to take the probability of a fine particle forming part of a strong force chain, P(strong) fine as an indicator of the extent to which the finer particles form part of the load-transferring matrix. For internal instability to occur the load-transferring matrix should be formed predominantly of coarse particles. Therefore the greater the extent to which the fines contribute to the loadtransferring matrix, the higher the internal stability of the soil will be.
The variation of P(strong) fine with χ and F fine is shown in Fig. 9. The pattern is similar to the relationship between α, χ and F fine presented in Fig. 7b, where for F fine ≥ 30 %, P(strong) fine > 0.5 and therefore fines play a significant role in supporting the fabric of the samples. For samples with F fine ≤ 25 %, P(strong) fine reduces with increasing χ as  Fig. 9 Variation of probability of a fine particle forming part of strong force chain with size ratio and fines content the fines are able to fit more efficiently in the voids and so are less likely to interact with the coarse particles, which dominate the strong force chains (for each sample the probability of a coarse particle forming part of a strong force chain is greater than 80 %).
This supports the hypothesis of Rahman et al. [22] that the role which cohesionless fines play in stress transfer diminishes with both χ and F fine when F fine < S * , (they refer to a threshold fines content equivalent to S*). For soils with χ = 2 the concept of a threshold fines content has little meaning. For soils with χ ≥ 4 care must be taken in defining this threshold content-when F fine < S * the role of the fines is primarily dependent on χ , in particular for the range 4 < χ < 6. However, for S * < F fine < S max the fines play a lesser role in stress transfer and this is dependent on both χ and F fine as shown in Fig. 9 for F fine = 25 %. An added complication is also that the values of S* and S max are density-dependent [6].

Conclusions
This paper has considered results from DEM simulations that investigated the effect of varying the cohesionless fines content and the size ratio between coarse and fine particles on the properties of dense bimodal materials. In particular the role of the fines in stress transfer has been studied, which has practical significance, for example when considering the internal stability of embankment dams. The stress-transfer and contact distributions in bimodal soils were shown to be controlled by χ and F fine and role of fines in stress transfer was quantified by the stress-reduction αfactor [4]. The following conclusions can be drawn from the work: (a) The DEM data are physically reasonable. It was found that the minimum void ratio which can be obtained by bimodal samples falls as the size-ratio, χ = D coarse /D fine increases, in agreement with experimental results presented by Lade et al. [18]. The connectivity distributions for fine-fine and fine-coarse contacts obtained from DEM simulations are broadly similar to the experimental results of Pinson et al. [19]. (b) The size ratios (χ ) for which a soil can be considered gap-graded was considered. Bimodal soils with a size ratio of χ ≤ 2 have stress transfer and coordination number characteristics that are similar to those for uniform soils, as the fines cannot completely fit between the coarse particles. Fine and coarse particles play an approximately equal role in stress transfer (α ≈ 1) regardless of F fine . Such soils should not be considered gap-graded in terms of stress-transmission.
(c) For χ ≥ 6 particle-scale evidence is shown for two critical fines contents at which transitions in soil fabric occur: (i) S*, at which the fines just fill the voids between coarse particles can be seen by an increase in the coordination number, Z; (ii) S max , where the fines begin to separate the coarse particles from one another is distinguished by a reduction in the coarse to coarse coordination number to Z coarse−coarse < 4. At the macro-scale this transition can be demonstrated by the skeleton void ratio, e sk , increasing above e coarse,max . For χ ≥ 6 stress transfer is dependent on both F fine and χ . In particular: (i) when F fine < S* the fines play only a minor role in stress transfer and α ≈ 0; (ii) when F fine > S max the coarse and fine play approximately equal roles and α ≈ 1; (iii) when S * < F fine < S max an increase in F fine leads to an increase in α, and an increase in χ leads to a reduction in α. (d) For χ = 4 at low fines contents behaviour which is intermediate between χ = 2 and χ = 6 is observed, with fines playing a reduced but significant role in stress transfer. (e) When F fine = 25 % the fines have a similar coordination number to the samples with F fine ≥ 30 % but α is lower. This shows that the fines play a supporting role to the coarse primary, stress-transmitting, matrix. (f) The influence of F fine and χ on α and on P(strong) fine , the probability of a fine particle forming part of a strong force chain, are similar. (g) Measurement of the void size distribution between the coarse particles shows that as more fine particles form part of the stress-transmitting matrix the size of voids between coarse particles increases. When α ≈ 0 the void size distribution reduces to that for a monodisperse sample. This occurs when χ ≥ 6 and F fine ≤ 25 %.