This is the first paper of a series of three, in which we report on numerical simulation studies of geometric and mechanical properties of static assemblies of spherical beads under an isotropic pressure. The influence of various assembling processes on packing microstructures is investigated. It is accurately checked that frictionless systems assemble in the unique random close packing (RCP) state in the low pressure limit if the compression process is fast enough, higher solid fractions corresponding to more ordered configurations with traces of crystallization. Specific properties directly related to isostaticity of the force-carrying structure in the rigid limit are discussed. With frictional grains, different preparation procedures result in quite different inner structures that cannot be classified by the sole density. If partly or completely lubricated they will assemble like frictionless ones, approaching the RCP solid fraction ΦRCP ≃ 0.639 with a high coordination number: z * ≃ 6 on the force-carrying backbone. If compressed with a realistic coefficient of friction µ = 0.3 packings stabilize in a loose state with Φ ≃ 0.593 and z * ≃ 4.5. And, more surprisingly, an idealized "vibration" procedure, which maintains an agitated, collisional régime up to high densities results in equally small values of z * while Φ is close to the maximum value ΦRCP. Low coordination packings have a large proportion (>10%) of rattlers -grains carrying no force -the effect of which should be accounted for on studying position correlations, and also contain a small proportion of localized "floppy modes" associated with divalent grains. Low pressure states of frictional packings retain a finite level of force indeterminacy even when assembled with the slowest compression rates simulated, except in the case when the friction coefficient tends to infinity. Different microstructures are characterized in terms of near neighbor correlations on various scales, and some comparisons with available laboratory data are reported, although values of contact coordination numbers apparently remain experimentally inaccessible.

In this third and final paper of a series, elastic properties of numerically simulated isotropic packings of spherical beads assembled by different procedures, as described in the first companion paper, and then subjected to a varying confining pressure, as reported in the second companion paper, are investigated. In addition to the pressure, which determines the stiffness of contacts because of Hertz's law, elastic moduli are chiefly sensitive to the coordination number z , which should not be regarded as a function of the packing density. Comparisons of numerical and experimental results for glass beads in the 10 kPa-10 MPa pressure range reveal similar differences between dry samples prepared in a dense state by vibrations and lubricated packings, so that the greater stiffness of the latter, in spite of their lower density, can be attributed to a larger coordination number. Effective medium type approaches, or Voigt and Reuss bounds, provide good estimates of bulk modulus B , which can be accurately bracketed, but badly fail for shear modulus G , especially in low z configurations under low pressure. This is due to the different response of tenuous, fragile networks to changes in load direction, as compared to load intensity. In poorly coordinated packings, the shear modulus, normalized by the average contact stiffness, tends to vary proportionally to the degree of force indeterminacy per unit volume, even though this quantity does not vanish in the rigid limit. The elastic range extends to small strain intervals and compares well with experimental observations on sands. The origins of nonelastic response are discussed. We conclude that elastic moduli provide access to mechanically important information about coordination numbers, which escape direct measurement techniques, and indicate further perspectives.

This is the second paper of a series of three investigating, by numerical means, the geometric and mechanical properties of spherical bead packings under isotropic stresses. We study the effects of varying the applied pressure P (from 1 or 10 kPa up to 100 MPa in the case of glass beads) on several types of configurations assembled by different procedures, as reported in the preceding paper [1]. As functions of P , we monitor changes in solid fraction Φ, coordination number z, proportion of rattlers (grains carrying no force) x0, the distribution of normal forces, the level of friction mobilization, and the distribution of near neighbor distances. Assuming that the contact law does not involve material plasticity or damage, Φ is found to vary very nearly reversibly with P in an isotropic compression cycle, but all other quantities, due to the frictional hysteresis of contact forces, change irreversibly. In particular, initial low P states with high coordination numbers lose many contacts in a compression cycle, and end up with values of z and x0 close to those of the most poorly coordinated initial configurations. Proportional load variations which do not entail notable configuration changes can therefore nevertheless significantly affect contact networks of granular packings in quasistatic conditions.

a b s t r a c tIn micromechanics of the elastic behaviour of granular materials, the macro-scale continuum elastic moduli are expressed in terms of micro-scale parameters, such as coordination number (the average number of contacts per particle) and interparticle contact stiffnesses in normal and tangential directions. It is well-known that mean-field theory gives inaccurate micromechanical predictions of the elastic moduli, especially for loose systems with low coordination number. Improved predictions of the moduli are obtained here for loose two-dimensional, isotropic assemblies. This is achieved by determining approximate displacement and rotation fields from the force and moment equilibrium conditions for small sub-assemblies of various sizes. It is assumed that the outer particles of these sub-assemblies move according to the mean field. From the particle displacement and rotation fields thus obtained, approximate elastic moduli are determined. The resulting predictions are compared with the true moduli, as determined from the discrete element method simulations for low coordination numbers and for various values of the tangential stiffness (at fixed value of the normal stiffness). Using this approach, accurate predictions of the moduli are obtained, especially when larger sub-assemblies are considered. As a step towards an analytical formulation of the present approach, it is investigated whether it is possible to replace the local contact stiffness matrices by a suitable average stiffness matrix. It is found that this generally leads to a deterioration of the accuracy of the predictions. Many micromechanical studies predict that the macroscopic bulk modulus is hardly influenced by the value of the tangential stiffness. It is shown here from the discrete element method simulations of hydrostatic compression that for loose systems, the bulk modulus strongly depends on the stiffness ratio for small stiffness ratios.

The elastic moduli of four numerical random isotropic packings of Hertzian spheres are studied. The four samples are assembled with different preparation procedures, two of which aim to reproduce experimental compaction by vibration and lubrication. The mechanical properties of the samples are found to change with the preparation history, and to depend much more on coordination number than on density.Secondly, the fluctuations in the particle displacements from the average strain are analyzed, and the way they affect the macroscopic behavior analyzed. It is found that only the average over equally oriented contacts of the relative displacement these fluctuations induce is relevant at the macroscopic scale. This average depends on coordination number, average geometry of the contact network and average contact stiffness. As far as the separate contributions from particle displacements and rotations are concerned, the former is found to counteract the average strain along the contact normal, while the latter do in the tangential plane. Conversely, the tangential components of the center displacements mainly arise to enforce local equilibrium, and have a small and generally stiffening effect at the macro-scale.Finally, the fluctuations and the shear modulus that result from two approaches available in the literature are estimated numerically. These approaches are both based on the equilibrium of a small-sized representative assembly. The improvement of these estimate with respect to the average strain assumption indicates that the fluctuations relevant to the macroscopic behavior occur with short correlation length.