Fundamental challenges in packing problems: from spherical to non-spherical particles

Abstract

Random packings of objects of a particular shape are ubiquitous in science and engineering. However, such jammed matter states have eluded any systematic theoretical treatment due to the strong positional and orientational correlations involved. In recent years progress on a fundamental description of jammed matter could be made by starting from a constant volume ensemble in the spirit of conventional statistical mechanics. Recent work has shown that this approach, first introduced by S. F. Edwards more than two decades ago, can be cast into a predictive framework to calculate the packing fractions of both spherical and non-spherical particles.

---

Adrian Baule

Dr Adrian Baule is a Lecturer in the School of Mathematical Sciences at Queen Mary University of London. He received his PhD from the University of Leeds with a thesis on exact relations in non-equilibrium statistical mechanics. He subsequently held postdoctoral positions at the Rockefeller University, New York, and at the Benjamin Levich Institute of the City College of New York. His current research interests include granular matter, active matter, and transport processes in complex media. Outside of work he enjoys swimming and art museums.

Hernán A. Makse

Hernán A. Makse currently serves as Professor of Physics at CUNY in New York City, wherein he is responsible for the Complex Networks and Soft Matter lab at the Levich Institute. He holds a PhD degree in Physics from Boston University and he is a Fellow of the American Physical Society. He has been author of numerous publications on the theory of complex systems and the physics of granular materials. He has visited over 50 countries in the world and whenever he travels, he usually attends a concert at the main Arts Center of the city and enjoys the local culture. He loves jazz, classic literature and opera and he is an avid footballer.

---

Introduction

In 1989 Edwards and Oakeshott made the remarkable proposal that the macroscopic properties of static granular matter can be calculated as ensemble averages over equiprobable jammed microstates controlled by the system volume.¹ In other words, granular matter is amenable to a statistical mechanical treatment, where the role of energy is played by the volume. Clearly, there is no a priori reason why such a treatment should be correct. Granular matter is profoundly out of equilibrium, since thermal fluctuations are essentially absent for the macroscopic length scales considered. In particular, there is no equivalent of Liouville's theorem for equilibrium systems due to the strongly dissipative nature of granular assemblies, which are dominated by static frictional forces. Nevertheless, the Edwards' ensemble approach has proven exceedingly useful in characterizing the properties of these athermal states of matter and continues to intrigue both experimentalists and theoreticians alike.

The main statements of this approach are:^1,2 (i) the distribution of jammed microstates is flat and independent of the compaction history leading to a natural definition of a configurational entropy S = ln Ω_Edw, where Ω_Edw is the number of jammed configurations. (ii) There is an equivalence between ensemble averages and time averages, if the system can explore its jammed configurations by some external drive (tapping or slow shearing). (iii) The compactivity X⁻¹ = ∂S/∂V characterizes the packing state analogous to the temperature in equilibrium systems. These strong assumptions have been scrutinized in various studies over recent years in order to obtain insight into the validity of Edwards' approach.³ Soft compaction experiments under continuous tapping have provided evidence for a reversible branch in the packing fraction for a variation of the tapping amplitude, indicating the existence of thermodynamic states.^4–7 Simple models of such a compaction dynamics have confirmed ergodicity and have been connected to a slow relaxation dynamics akin to the relaxation in glasses.^8–12 One signature of such a slow dynamics is the existence of a non-equilibrium fluctuation–dissipation relation.^8,9,13 Indeed, the effective temperature appearing in FDRs under perturbations agrees with the configurational temperature T_eff⁻¹ = ∂S/∂E in Edwards' framework.¹³

Ergodicity has also been demonstrated explicitly in more realistic simulations.^14,15 The compactivity has been measured in simulations and experiments.^16–20 On the other hand, results on the equiprobability of microstates are mixed. By evaluating the probabilities of jammed microstates in small clusters a break down of the flat distribution assumption has been demonstrated,^21–24 which might be traced back to the packing protocol used.¹⁵ Recent studies have investigated the equilibration of granular subsystems in contact,^15,25,26 providing further insight into the thermodynamic nature of granular matter.

Ultimately, the success of any statistical mechanical theory needs to be measured by the comparison with experiments. One key problem in Edwards' approach is to identify a suitable volume function, which parametrizes the total volume of the packing as a function of the particle configurations (positions and orientations), replacing the role of the Hamiltonian.² Here, different conventions can be employed to partition the total volume into cells associated with each particle,^27,28 the simplest of which is the Voronoi tessellation.^29,30 In 3D these exact volume functions are difficult to handle analytically, so that reduced representations are sought. The thermodynamic nature of packings suggests to use a coarse-grained description of the volume function in terms of observables such as the average number of contacts z (coordination number).^31–34 In turn, z is determined by the force transmission in the contact network leading to a mechanically stable packing in which forces and torques on each particle balance.³⁵ For the force network ensemble approaches similar to the volume ensembles have been introduced in order to explain the observed force distribution³⁶ from an entropy maximization.^37–41 The stress tensor has also been considered as a conserved quantity leading to a different class of ensembles.^41,42

Recently, the force transmission has been treated on a random graph under local mechanical stability constraints resulting in quantitative predictions for the force distribution and the value of z using a cavity method.⁴³ The problem of finding the densest random packing can be similarly formulated as a constraint optimization problem: random close packings appear as the ground state of the volume ensemble restricted to disordered packings as X → 0 for a given z.^31,32 This picture highlights that jamming falls into the class of NP optimization problems,⁴⁴ which can be tackled successfully with the methods of statistical mechanics such as cavity methods.⁴³ A full solution needs to combine the two approaches for the force and volume ensembles, where the Hamiltonian that enforces jamming is a function of both the particle configurations and the contact forces on a random contact network. These ensemble approaches are thus similar in spirit to other recent works that consider jamming as the infinite pressure limit of metastable glass phases.^44–51 Here, one considers instead of the Edwards entropy S, the “glass complexity” in order to obtain the statistics of the metastable basins as the pressure diverges. Treatments of this problem based on the random-first order transition picture and replica theory have been performed.⁵⁰

From spheres to non-spheres

Random packings of hard objects appear in a broad range of scientific and engineering fields like self-assembly of nanoparticles, liquid crystals, glass formation and bio-materials.⁵² In fact, the question of how densely objects of a particular shape can fill a given volume is probably one of the most ancient scientific problems and still of great practical importance for all industries involved in granular processing. The densest random packing has been extensively studied in experiments and simulations for spheres, which typically reach a maximal volume fraction of ϕ ≈ 0.64 in monodisperse assemblies.^53,54 This value is quite robustly reproducible and commonly referred to as random close packing (RCP) density. However, much less is known about anisotropic shapes, despite the fact that all shapes in nature deviate from the ideal sphere. A theoretical investigation of the packing problem has proven notoriously difficult due to the strong positional and orientational correlations of dense packings. In fact, a mathematically rigorous treatment of random sphere packings has been the outstanding component in T.C. Hales' proof of Kepler's conjecture on the densest packing of spheres.⁵⁵

Recent empirical work has focused on packings of anisotropic shapes like ellipsoids, spherocylinders, and tetrahedra, which can achieve considerably denser volume fractions than the spherical RCP^56–68 (see Table 1). In fact, a conjecture attributed to Ulam in the context of regular packings and recently also formulated for random packings⁶⁷ states that the sphere is, indeed, the worst packing object among all convex shapes. This suggests to improve packing fractions by searching in the space of object shapes, but in the absence of any general theory, this search could so far only be performed on a case-by-case basis using experiments and simulations. A caveat of such empirical studies is the strong protocol dependence of the final close packed state even for the same shape. While the range of achieved volume fractions is relatively small for spheres,⁵⁴ recent studies of spherocylinder packings, e.g., exhibit a much greater variance depending on the algorithm used.^{56,57,60–62,65–68} Further theoretical insight is needed, which can be obtained by considering a coarse-grained distribution for the Voronoi volumes in the packing, as discussed next.

Table 1 Overview of maximal packing fractions found for disordered packings of a selection of regular shapes. We observe that spheres are the worst-packing objects among all shapes in 3d as conjectured by Ulam,⁸⁸ while tetrahedra achieve the densest disordered packing. We note that the tetrahedron is the only shape known that packs in a disordered arrangement denser than spheres in the FCC crystal (ϕ_FCC = 0.7405). Ellipsoids and lens-shaped particles pack very close to this value

Shape	ϕ max simulation	ϕ max experiment	ϕ max theory
Sphere	0.645 (ref. 54)	0.64 (ref. 53)	0.634 (ref. 31)
M&M candy		0.665 (ref. 58)
Dimer	0.703 (ref. 91)		0.707 (ref. 34)
Oblate ellipsoid	0.707 (ref. 58)
Prolate ellipsoid	0.716 (ref. 58)
Spherocylinder	0.722 (ref. 68)		0.731 (ref. 34)
Lens-shaped particle			0.736 (ref. 34)
Octahedron	0.697 (ref. 67)
Icosahedron	0.707 (ref. 67)
Dodecahedron	0.716 (ref. 67)
General ellipsoid	0.735 (ref. 58)	0.74 (ref. 59)
Tetrahedron	0.7858 (ref. 63)	0.76 (ref. 64)

---

A mean-field theory for random close packings

In the Voronoi convention one associates with each particle the fraction of space that is closer to this particle than to any other one. This defines the Voronoi volume W_i of a particle i, which depends on the configurations of all remaining particles x_j = (r_j, t̂_j), (including position r_j and orientation t̂_j). The total volume V occupied by N particles is and the packing fraction of monodisperse particles of volume V_α is ϕ = NV_α/V. In order to determine W_i one has to know the Voronoi boundary (VB) between two particles i and j, which is the hypersurface that contains all points equidistant to the surfaces of both particles and thus depends on the particle shape and their relative configuration. The boundary of W_i then follows from a global minimization procedure over all pairwise VB.³² In order to take into account multi-particle correlations in the packing, we use a statistical treatment where the overall volume is expressed in terms of an average Voronoi volume: V = NW̄(z), so that ϕ = V_α/W̄(z). Instead of an exact description in terms of all configurations {x₁,…, x_N}, the average Voronoi volume is characterized by the coordination number z, which denotes the average number of contacting neighbours in the packing. We derive a self-consistent equation for the coarse-grained volume function W̄(z) of monodisperse particles:^31,32,34Here, V_α is the volume of a single particle and σ(z) is the average free-surface of particles at contact, which can be estimated from local configurations of z contacting particles. Formally, the integrand on the right hand side can be considered as the cumulative distribution function P(c) containing the probability to find the boundary of the Voronoi volume in the direction ĉ at a value larger than c. This quantity can be interpreted geometrically as the probability to find all N − 1 particles outside a volume Ω centred at c from the reference particle (see Fig. 1a). The particular form of P(c) results from a factorization into bulk and contact terms, which are motivated from the dominant contributions in the radial distribution function.^31,32,34,69

Fig. 1 The Voronoi excluded volume and surface for spherocylinders. (a) The volume Ω (red) is excluded for the remaining N − 1 particles in the packing because otherwise the Voronoi boundary would be found at a value smaller than c in the direction ĉ. We draw the usual hard-core excluded volume V_ex⁷⁰ in blue. (b) The overlap of Ω and V_ex defines the Voronoi excluded volume V* (red) and the Voronoi excluded surface S* (green). Figure taken from ref. 34.

The quantities V* and S* are the Voronoi excluded volume and surface, which extend the usual hard-core excluded volume of equilibrium systems V_ex⁷⁰ to packings. The volume V* is the volume excluded by Ω for bulk particles and takes into account the overlap between Ω and V_ex: where the bar denotes an orientational average. Likewise, S* denotes the surface excluded by Ω for contacting particles: Plots of V* and S* for spherocylinders are shown in Fig. 1b. Analytical expressions for V* and S* can be derived in the spherical limit in closed form.^31,32 For non-spherical shapes analytic expressions for the VB can be derived using a suitable decomposition of the shape into overlapping and/or intersecting spheres. This leads to exact expressions for V* and S*, which can be evaluated numerically.⁷¹ Interestingly, in the limit α → 1, eqn (1) admits an exact solution for spheres: As a consequence, we obtain an equation of state for spherical packings^31,32

which predicts the limiting values ϕ = 0.536 and ϕ = 0.634 under the isostatic conditions z = 4 and z = 6 for infinitely rough and frictionless spheres, respectively. Using the thermodynamic framework one can show that these two values are reached in the limits of infinite and zero compactivity, respectively.³¹ Therefore, the spherical equation of state leads to a statistical interpretation of RCP as the ground state of disordered sphere packings. The predictions for the limiting values are in good agreement with the values found in experiments and simulations for both random loose packings and RCP of spheres.

Under deformation from the sphere, higher packing fractions are typically reached, where the spherical limit appears as a singular point in the ϕ(α) plane. Moreover, smooth shapes close to the sphere are not isostatic but hypostatic with z < 2d_f due to redundancies in the force and torque balance equations.^72,73 The variation z(α) is obtained by considering the average effective number of degrees of freedom d̃_f defined as the number of linearly independent force and torque balance equations: z = 〈d̃_f(α)〉.³⁴ Here, the probability of redundant configurations with d̃_f < d_f can be estimated by re-weighting all configurations by rotating into states of maximal redundancy.

The existence of redundant configurations explains the observed convergence in z(α) to values close to 8 for spherocylinders with large aspect ratios:^62,68 for long spherocylinders the contacts are predominantly on the cylindrical part so that all normal forces are coplanar. As a consequence, the effective number of degrees of freedom is reduced by one, leading to z = 8.³⁴ The requirement of local force and torque balance can also be formulated as a constraint optimization problem on a factor graph, which describes the force transmission on a single particle.⁴³ Solving this problem with standard methods such as the cavity method predicts values of z in frictional packings and also allows for the computation of the distribution of contact forces in good agreement with experimental results.⁴³

Phase-diagram of jammed isotropic and anisotropic particles

The combination of the results for W̄(z) and z(α) leads to a complete theoretical prediction for the packing density ϕ(α) = V_α/W̄(z(α)) of non-spherical particles without any adjustable parameters.³⁴ We estimate the maximum density of spherocylinders at α = 1.3 with a density ϕ_max = 0.731 in good agreement with empirical data. The theory also reproduces well the density of dimers, estimating a maximum at α = 1.3 with ϕ_max ≈ 0.707. We have also calculated the packing fraction of lens-shaped particles, which can serve as approximations for oblate ellipsoidal shapes. Our theory yields ϕ_max = 0.736 for α = 0.8. This shape represents the densest random packing of an axisymmetric shape known so far. The appearance of a maximum in ϕ for non-spherical shapes close to the sphere has been explained in a simple qualitative picture on the basis of the excluded volume V_ex.⁵⁶ For α close to 1, the ratio V_ex/V_α changes only slightly from the spherical value and a density increase results due to the additional orientational degrees of freedom, whereby the particles can fit into gaps by rotating, similar to the increase in packing efficiency due to polydispersity.⁷⁴ For larger α, V_ex exceeds V_α while z remains constant, so that the packing is dominated by the excluded volume and the packing fraction decreases. This argument can explain qualitatively the observed larger packing fraction of spherocylinders compared with dimers. The ratio V_ex/V_α is approximately equal for both shapes up to α ≈ 1.2, but for larger α the ratio for dimers increases beyond that of spherocylinders. The packing densities derived in our framework are interpreted as upper bounds to the empirically obtained densities and correspond to maximally random jammed states⁷⁵ by construction, since the distribution of contact angles in the first coordination shell is imposed to be uniform, avoiding any partial order.

By plotting z(α) against ϕ(α) parametrically as a function of α, we obtain a phase diagram in the z–ϕ plane (Fig. 2). Surprisingly, we find that both dimer and spherocylinder packings appear as smooth continuations of spherical packings. The analytic form of this continuation from the spherical random branch can be derived (blue dashed line in Fig. 2).³⁴ A comparison of our theoretical results with empirical data for a large variety of shapes highlights that the analytic continuation provides a boundary line in the z–ϕ phase diagram. Maximally dense disordered packings appear to the left of this boundary, while the packings to the right of it are partially ordered. The spherical ordered branch provides another boundary, which separates tetrahedra from all other shapes: tetrahedra are the only shape that pack in a disordered way denser than spheres in a FCC crystal. We observe that the maximally dense packings of dimers, spherocylinders, lens-shaped particles and tetrahedra all lie surprisingly close to the analytic continuation of RCP. Whether there is any deeper meaning to this remains an open question.

Fig. 2 Phase diagram of jammed matter. We plot our results for dimers and spherocylinders in the z–ϕ plane together with results from the literature for frictionless disordered packings of a selection of regular shapes. We have selected those shapes for which the z and ϕ values have been determined in the same simulation. The predicted spherical random branch eqn (2) (ref. 31) (thick black line) and a conjectured first order disorder–order transition at RCP for spheres⁷⁷ (dotted and thin black lines) are also indicated. We observe that the analytic continuation of RCP under deformation into dimers and spherocylinders provides an empirical bound to disordered packings in the phase diagram. The symmetry of the shape indicates the possible values of the coordination number z: (i) spheres have z between 4 (infinitely rough) and 6 (frictionless). (ii) Axisymmetric particles have z between 6 and 10. (iii) Fully aspherical particles have z between 10 and 12. Note that for polyhedra, z is associated with the total degrees of freedom blocked by the different types of contacts (face–face, face–vertex, vertex–vertex, face–edge).⁶⁷ The data point for lens-shaped particles is a theoretical prediction.³⁴

The picture that emerges is that spherical packings can be generated between the RLP and RCP limits by variation of the inter-particle friction, since this leads to an increase in the coordination number under the isostatic condition from z = 4 to z = 6. Beyond RCP, the spherical equation of state can be continued smoothly by deforming the sphere into elongated shapes. Moreover, the spherical RCP is interpreted as the freezing point of disordered sphere packings, associated with a melting point at ϕ = 0.68.^76,77 The signature of this disorder–order transition is a discontinuity in the entropy density of jammed configurations as a function of the compactivity. This highlights the fact that beyond RCP, denser packing fractions of spheres can only be reached by partial crystallization up to the homogeneous FCC crystal phase.⁷⁵

Conclusions and outlook

The first-order transition of jammed spheres identified within Edwards' thermodynamics⁷⁷ is reminiscent of the entropy induced phase transition of equilibrium hard spheres, which is found at ϕ = 0.494 and ϕ = 0.545, respectively. However, it should be emphasized that the physical origin of these two transitions is fundamentally different: the equilibrium phase transition is a consequence of the maximization of the conventional entropy, while the transition at RCP of jammed spheres is driven by the competition between volume minimization and maximization of the entropy of jammed configurations S. For anisotropic particles at equilibrium, a disorder–order phase transition appears, e.g., between isotropic and nematic phases of elongated shapes: for large α, Onsager's theory of equilibrium hard rods predicts a first order isotropic–nematic transition with freezing point at the rescaled density ϕα = 3.29 and melting point at ϕα = 4.19.⁷⁰ For colloidal suspensions of more complex shapes like polyhedra, both liquid crystal as well as plastic crystal and even quasicrystal phases have been found.^78–80 By analogy with the case of jammed spheres, one might wonder whether packings of non-spherical particles exhibit similar transitions that might be characterized in the z–ϕ phase diagram. Packings of hard thin rods indeed satisfy a scaling law, where the RCP has been experimentally identified at ϕα ≈ 5.4.⁸¹

The Edwards' approach thus helps to elucidate how macroscopic properties of granular matter arise from the anisotropy of the constituents – one of the central questions in present day materials science.^82,83 A better understanding of this problem will facilitate, e.g., the engineering of new functional materials with particular mechanical responses by tuning the shape of the building blocks. A search in the space of object shapes for optimization can be performed by considering a small number of spheres and systematically exploring the different possible configurations.⁸⁴

Our approach eqn (1) can be applied to a large variety of both convex and non-convex shapes. The key is to parametrize the Voronoi boundary between two such shapes, which allows for the calculation of the Voronoi excluded volume and surface. In fact, analytical expressions for the Voronoi boundary can be derived following an exact algorithm for arbitrary shapes by decomposing the shape into overlapping and intersecting spheres (see Fig. 3). Therefore, a systematic search for maximally dense packings in the space of given object shapes can be performed using our framework. Extensions to mixtures and polydisperse packings can also be formulated. So far, exhaustive searches for dense packings have only been performed for ordered packings using computer simulations^85,86 and a combination of analytic and simulation techniques.⁸⁷ This has elucidated in particular the validity of Ulam's conjecture that the sphere is the worst packing object in 3d.⁸⁸ Analytical progress to prove this conjecture locally, that is, for shapes deformed from the sphere, has recently been made.⁸⁹

Fig. 3 Decomposition of various shapes to calculate the Voronoi boundary. The Voronoi boundary (VB) between two particles is defined as the hypersurface that contains all points equidistant to their surfaces. This implies that the VB between two equal spheres, e.g., is that between two points at the centers of the spheres, so that the VB is generated effectively by the interaction of two points (a). Likewise, the VB between two spherocylinders is due to the effective interaction of two lines, since spherocylinders can be represented as dense overlaps of spheres (d). Arbitrary shapes can be decomposed into dense overlaps of spheres following certain design principles.⁹⁰ The VB between two such shapes can then be calculated following an exact algorithm that considers the effective Voronoi interactions between points and lines (a–d).³⁴ For shapes that are not naturally given as overlapping spheres (e–h), we propose alternatively an approximation in terms of a small number of intersecting spheres. In this way, two intersecting spheres (a lens-shaped particle) approximate an oblate ellipsoid and four intersecting spheres approximate a tetrahedron. The effective Voronoi interactions are then between points, lines, and anti-points (indicated by crosses).³⁴ Anti-points arise from the inversion of the effective interaction between the spheres in the decomposition. This is evident in the case of lens-shaped particles (e), where the interaction between the spheres is inverted compared to the case of dimers (b). The VB between two tetrahedra is then due to the interaction between the vertices (leading to four point interactions), the edges (leading to six line interactions), and the faces (leading to four anti-point interactions). This approach can be generalized to arbitrary polyhedra. Figure taken from ref. 34.

A more systematic investigation of disordered packings can shed light on the validity of a random variant of Ulam's conjecture, which so far has only been investigated in simulations.⁶⁷ Our analytic continuation from RCP highlights that this conjecture might hold more generally than previously assumed, containing not only convex shapes, but also a significant class of non-convex ones. Ultimately, our approach might lead to a more exhaustive theoretical investigation of Ulam's conjecture. Along the way one might be able to answer important questions such as if a shape that packs denser in a random configuration than in a regular one exists.⁷² Such objects could represent optimal glass formers with far reaching consequences for materials science.

References

1. S. F. Edwards and R. B. S. Oakeshott, Physica A, 1989, 157, 1080.
2. S. F. Edwards, in Granular matter: an interdisciplinary approach, ed. A. Mehta, Springer, New York, 1994, pp. 121–140.
3. M. P. Ciamarra, P. Richard, M. Schroter and B. P. Tighe, Soft Matter, 2012, 8, 9731.
4. E. R. Nowak, J. B. Knight, E. Ben-Naim, H. M. Jaeger and S. R. Nagel, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 1998, 57, 1971.
5. P. Philippe and D. Bideau, Europhys. Lett., 2002, 60, 677.
6. J. Brujić, P. Wang, C. Song, D. L. Johnson, O. Sindt and H. A. Makse, Phys. Rev. Lett., 2005, 95, 128001.
7. P. Richard, M. Nicodemi, R. Delannay, P. Ribiére and D. Bideau, Nat. Mater., 2005, 4, 121.
8. M. Nicodemi, Phys. Rev. Lett., 1999, 82, 3734.
9. A. Barrat, J. Kurchan, V. Loreto and M. Sellitto, Phys. Rev. Lett., 2000, 85, 5034.
10. J. Brey, A. Prados and B. Sanchez-Rey, Phys. A, 2000, 275, 310 , ISSN 0378-4371.
11. D. S. Dean and A. Lefèvre, Phys. Rev. Lett., 2001, 86, 5639.
12. G. D'Anna and G. Gremaud, Nature, 2001, 413, 407.
13. H. A. Makse and J. Kurchan, Nature, 2002, 415, 614.
14. M. Pica Ciamarra, A. Coniglio and M. Nicodemi, Phys. Rev. Lett., 2006, 97, 158001.
15. K. Wang, C. Song, P. Wang and H. A. Makse, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 2012, 86, 011305.
16. M. Schröter, D. I. Goldman and H. L. Swinney, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 2005, 71, 030301.
17. F. Lechenault, F. da Cruz, O. Dauchot and E. Bertin, J. Stat. Mech.: Theory Exp., 2006, 2006, P07009.
18. P. Ribiére, P. Richard, P. Philippe, D. Bideau and R. Delannay, Eur. Phys. J. E, 2007, 22, 249.
19. C. Briscoe, C. Song, P. Wang and H. A. Makse, Phys. Rev. Lett., 2008, 101, 188001.
20. S. McNamara, P. Richard, S. K. de Richter, G. Le Caër and R. Delannay, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 2009, 80, 031301.
21. N. Xu, J. Blawzdziewicz and C. S. O'Hern, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 2005, 71, 061306.
22. G.-J. Gao, J. Bławzdziewicz and C. S. O'Hern, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 2006, 74, 061304.
23. G.-J. Gao, J. Blawzdziewicz, C. S. O'Hern and M. Shattuck, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 2009, 80, 061304.
24. N. Xu, D. Frenkel and A. J. Liu, Phys. Rev. Lett., 2011, 106, 245502.
25. F. Lechenault and K. E. Daniels, Soft Matter, 2010, 6, 3074.
26. J. G. Puckett and K. E. Daniels, Phys. Rev. Lett., 2013, 110, 058001.
27. R. C. Ball and R. Blumenfeld, Phys. Rev. Lett., 2002, 88, 115505.
28. R. Blumenfeld and S. F. Edwards, Phys. Rev. Lett., 2003, 90, 114303.
29. A. Okabe, B. Boots, K. Sugihara and S. Nok Chiu, Spatial Tessellations: Concepts and Applications of Voronoi Diagrams, Wiley-Blackwell, 2000.
30. H. A. Makse, J. Brujić and S. F. Edwards, in The Physics of Granular Media, ed. H. Hinrichsen and D. E. Wolf, Wiley-VCH, 2004.
31. C. Song, P. Wang and H. A. Makse, Nature, 2008, 453, 629.
32. C. Song, P. Wang, Y. Jin and H. A. Makse, Physica A, 2010, 389, 4497.
33. P. Wang, C. Song, Y. Jin and H. A. Makse, Physica A, 2011, 390, 427.
34. A. Baule, R. Mari, L. Bo, L. Portal and H. A. Makse, Nat. Commun., 2013, 4, 2194.
35. S. Alexander, Phys. Rep., 1998, 296, 65.
36. G. W. Delaney, T. Di Matteo and T. Aste, Soft Matter, 2010, 6, 2992.
37. N. Kruyt and L. Rothenburg, J. Solids Struct., 2002, 39, 571.
38. K. Bagi, Granular Matter, 2003, 5, 45.
39. J. Goddard, Int. J. Solids Struct., 2004, 41, 5851.
40. J. H. Snoeijer, T. J. H. Vlugt, M. van Hecke and W. van Saarloos, Phys. Rev. Lett., 2004, 92, 054302.
41. S. Henkes and B. Chakraborty, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 2009, 79, 061301.
42. S. Henkes, C. S. O'Hern and B. Chakraborty, Phys. Rev. Lett., 2007, 99, 038002.
43. L. Bo, R. Mari, C. Song and H. A. Makse, 2013, arXiv preprint p. 1310.7973.
44. F. Krzakala and J. Kurchan, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 2007, 76, 021122.
45. T. Aste and A. Coniglio, Europhys. Lett., 2004, 67, 165.
46. R. D. Kamien and A. J. Liu, Phys. Rev. Lett., 2007, 99, 155501.
47. R. Mari, F. Krzakala and J. Kurchan, Phys. Rev. Lett., 2009, 103, 025701.
48. I. Biazzo, F. Caltagirone, G. Parisi and F. Zamponi, Phys. Rev. Lett., 2009, 102, 195701.
49. M. Hermes and M. Dijkstra, EPL, 2010, 89, 38005.
50. G. Parisi and F. Zamponi, Rev. Mod. Phys., 2010, 82, 789.
51. G. Parisi and F. Zamponi, J. Chem. Phys., 2005, 123, 144501.
52. S. Torquato and F. H. Stillinger, Rev. Mod. Phys., 2010, 82, 2633.
53. J. D. Bernal and J. Mason, Nature, 1960, 910–911.
54. M. Skoge, A. Donev, F. H. Stillinger and S. Torquato, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 2006, 74, 041127.
55. T. C. Hales, Ann. Math., 2005, 162, 1065.
56. S. Williams and A. Philipse, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 2003, 67, 051301.
57. C. Abreu, F. Tavares and M. Castier, Powder Technol., 2003, 134, 167.
58. A. Donev, I. Cisse, D. Sachs, E. Variano, F. Stillinger, R. Connelly, S. Torquato and P. Chaikin, Science, 2004, 303, 990.
59. W. Man, A. Donev, F. H. Stillinger, M. T. Sullivan, W. B. Russel, D. Heeger, S. Inati, S. Torquato and P. M. Chaikin, Phys. Rev. Lett., 2005, 94, 198001.
60. X. Jia, M. Gan, R. A. Williams and D. Rhodes, Powder Technol., 2007, 174, 10.
61. M. Bargiel, Computational Science-ICCS2008, 2008, vol. 5102, p. 126.
62. A. Wouterse, S. Luding and A. P. Philipse, Granular Matter, 2009, 11, 169.
63. A. Haji-Akbari, M. Engel, A. S. Keys, X. Zheng, R. G. Petschek, P. Palffy-Muhoray and S. C. Glotzer, Nature, 2009, 462, 773.
64. A. Jaoshvili, A. Esakia, M. Porrati and P. M. Chaikin, Phys. Rev. Lett., 2010, 104, 185501.
65. P. Lu, S. Li, J. Zhao and L. Meng, Science, 2010, 53, 2284.
66. A. V. Kyrylyuk, M. A. van de Haar, L. Rossi, A. Wouterse and A. P. Philipse, Soft Matter, 2011, 7, 1671.
67. Y. Jiao and S. Torquato, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 2011, 84, 041309.
68. J. Zhao, S. Li, R. Zou and A. Yu, Soft Matter, 2012, 8, 1003.
69. Y. Jin, P. Charbonneau, S. Meyer, C. Song and F. Zamponi, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 2010, 82, 051126.
70. L. Onsager, Ann. N. Y. Acad. Sci., 1949, 51, 627.
71. L. Portal, M. Danisch, A. Baule, R. Mari and H. A. Makse, J. Stat. Mech.: Theory Exp., 2013, P11009.
72. P. M. Chaikin, A. Donev, W. Man, F. H. Stillinger and S. Torquato, Ind. Eng. Chem. Res., 2006, 45, 6960.
73. A. Donev, R. Connelly, F. H. Stillinger and S. Torquato, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 2007, 75, 051304.
74. M. Danisch, Y. Jin and H. A. Makse, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 2010, 81, 051303.
75. S. Torquato, T. Truskett and P. Debenedetti, Phys. Rev. Lett., 2000, 84, 2064.
76. C. Radin, J. Stat. Phys., 2008, 131, 567.
77. Y. Jin and H. A. Makse, Physica A, 2010, 389, 5362.
78. U. Agarwal and F. A. Escobedo, Nat. Mater., 2011, 10, 230.
79. P. F. Damasceno, M. Engel and S. C. Glotzer, Science, 2012, 337, 453.
80. M. Marechal and H. Löwen, Phys. Rev. Lett., 2013, 110, 137801.
81. A. Philipse, Langmuir, 1996, 12, 1127.
82. S. C. Glotzer and M. Solomon, Nat. Mater., 2007, 6, 557.
83. T. Borzsonyi and R. Stannarius, Soft Matter, 2013, 9, 7401.
84. M. Z. Miskin and H. M. Jaeger, Nat. Mater., 2013, 12, 326.
85. J. de Graaf, R. van Roij and M. Dijkstra, Phys. Rev. Lett., 2011, 107, 155501.
86. J. de Graaf, L. Filion, M. Marechal, R. van Roij and M. Dijkstra, J. Chem. Phys., 2012, 137, 214101.
87. E. R. Chen, D. Klotsa, M. Engel, P. F. Damasceno and S. C. Glotzer, 2013, arXiv preprint p. 1309.2662.
88. M. Gardner, The Colossal Book of Mathematics: Classic Puzzles, Paradoxes, and Problems, Norton, 2001.
89. Y. Kallus, 2013, arXiv preprint p. 1305.0289.
90. C. L. Phillips, J. A. Anderson, G. Huber and S. C. Glotzer, Phys. Rev. Lett., 2012, 108, 198304.
91. S. Faure, A. Lefebvre-Lepot and B. Semin, in ESAIM: Proceedings, ed. M. Ismail, B. Maury and J.-F. Gerbeau, 2009, vol. 28, pp. 13–32.

---