Basic Understanding of Condensed Phases of Matter via Packing Models

Basic Understanding of Condensed Phases of Matter via Packing Models

S. Torquato Department of Chemistry, Department of Physics, Princeton Center for Theoretical Science,Princeton Institute for the Science and Technology of Materials, and Program in Applied and Computational Mathematics, Princeton University, Princeton New Jersey, 08544 USA

Packing problems have been a source of fascination for millenia and their study has produced a rich literature that spans numerous disciplines. Investigations of hard-particle packing models have provided basic insights into the structure and bulk properties of condensed phases of matter, including low-temperature states (e.g., molecular and colloidal liquids, crystals and glasses), multiphase heterogeneous media, granular media, and biological systems. The densest packings are of great interest in pure mathematics, including discrete geometry and number theory. This perspective reviews pertinent theoretical and computational literature concerning the equilibrium, metastable and nonequilibrium packings of hard-particle packings in various Euclidean space dimensions. In the case of jammed packings, emphasis will be placed on the “geometric-structure” approach, which provides a powerful and unified means to quantitatively characterize individual packings via jamming categories and “order” maps. It incorporates extremal jammed states, including the densest packings, maximally random jammed states, and lowest-density jammed structures. Packings of identical spheres, spheres with a size distribution, and nonspherical particles are also surveyed. We close this review by identifying challenges and open questions for future research.

I Introduction

We will call a packing a large collection of nonoverlapping (i.e., hard) particles in either a finite-sized container or in -dimensional Euclidean space . Exclusion-volume effects that often arise in dense many-particle systems in physical and biological contexts are naturally captured by simple packing models. They have been studied to help understand the symmetry, structure and physical properties of condensed matter phases, including liquids, glasses and crystals, as well as the associated phase transitions. Mayer and Mayer (1940); Bernal (1960, 1965); Stillinger et al. (1964); Stillinger and Salsburg (1969); Weeks et al. (1971); Ashcroft and Mermin (1976); Woodcock and Angell (1981); Hansen and McDonald (1986); Speedy (1994); Dijkstra and Frenkel (1994); Chaikin and Lubensky (1995); Rintoul and Torquato (1996) Packings also serve as excellent models of the structure of multiphase heterogeneous materials, Felderhof (1982); Torquato (2002); Zohdi (2006); Liasneuski et al. (2014) colloids, Russel et al. (1989); Chaikin and Lubensky (1995); Torquato (2009) suspensions,Kim and Torquato (1991); Spangenberg et al. (2014) and granular media, Edwards (1994); Zohdi (2004) which enables predictions of their effective transport, mechanical and electromagnetic properties.Torquato (2002)

Packing phenomena are also ubiquitous in biological contexts and occur in systems across a wide spectrum of length scales. This includes DNA packing within the nucleus of a cell, Cremer and Cremer (2001)“crowding” of macromolecules within living cells,  Ellis (2001) the packing of cells to form tissue, Torquato (2002); Gevertz and Torquato (2008); Jiao et al. (2014); Chen et al. (2016) the fascinating spiral patterns seen in plant shoots and flowers (phyllotaxis), Prusinkiewicz and Lindenmayer (1990); Nisoli et al. (2010) and the competitive settlement of territories. Tanemura and Hasegawa (1980); Torquato (2002)

Understanding the symmetries and other mathematical properties of the densest sphere packings in various spaces and dimensions is a challenging area of long-standing interest in discrete geometry and number theory Conway and Sloane (1998); Cohn and Elkies (2003); Rogers (1964) as well as coding theory. Shannon (1948); Conway and Sloane (1998); Cohn and Kumar (2007) Packing problems are mathematically easy to pose, but they are notoriously difficult to solve rigorously. For example, in 1611, Kepler was asked: What is the densest way to stack equal-sized cannon balls? His solution, known as “Kepler’s conjecture,” was the face-centered-cubic (fcc) arrangement (the way your green grocer stacks oranges). Gauss Gauss (1831) proved that this is the densest Bravais lattice packing (defined below). Remarkably, nearly four centuries passed before Hales proved the general conjecture that there is no other sphere packing in three-dimensional (3D) Euclidean space whose density can exceed that of the fcc packing.Hales (2005) Even the proof that the densest packing of congruent (identical) circles in the plane is the triangular-lattice packing appeared only 80 years ago; see Refs. Rogers, 1964; Conway and Sloane, 1998 for the history of this problem.

One objective of this perspective is to survey recent developments concerning the simplest but venerable packing model: identical frictionless spheres in the absence of gravity subject only to a nonoverlap constraint, i.e., the spheres do not interact for any nonoverlapping configuration. This “stripped-down” hard-sphere model can be viewed as the particle-system analog of the idealized Ising model for magnetic systems,Torquato and Stillinger (2010) which is regarded as one of the pillars of statistical mechanics. Onsager (1944); Domb (1960); Gallavotti (1999) More complex packing models can include interactions that extend beyond hard-core distances, but our main concern here is the aforementioned classic version. This model enables one to uncover unifying principles that describe a broad range of phenomena, including the nature of equilibrium states of matter (e.g., liquids and crystals), metastable and nonequilibrium states of matter (e.g., supercooled liquids and structural glasses), and jamming phenomena.

Jammed sphere packings represent an important subset of hard-sphere configurations and have attracted great attention because they capture salient characteristics of crystals, glasses, granular media and biological systems Edwards (1994); Torquato et al. (2000); O’Hern et al. (2003); Parisi and Zamponi (2010); Torquato and Stillinger (2010); Jiao et al. (2014); Jaeger (2015) and naturally arise in pure mathematical contexts.Conway and Sloane (1998); Martinet (2003); Cohn et al. (2011) Jammed packings are those in which each particle has a requisite number contacting particles in order for the packing to achieve a certain level of mechanical stability. This review focuses on the so-called “geometric-structure” approach to jammed particle packings, which provides a powerful and unified means to quantitatively understand such many-particle systems. It incorporates not only the maximally dense packings and disordered jammed packings, but a myriad of other significant jammed states, including maximally random jammed states,Torquato et al. (2000); Torquato and Stillinger (2010) which can be regarded to be prototypical structural glasses, as well as the lowest-density jammed packings.Torquato and Stillinger (2007)

Importantly, the simplified hard-sphere model embodies the primary structural attributes of dense many-particle systems in which steep repulsive interparticle interactions are dominant. For example, the densest sphere packings are intimately related to the (zero-temperature) ground states of such molecular systemsStillinger (2001) and high-pressure crystal phases. Indeed, the equilibrium hard-sphere model Hansen and McDonald (1986) also serves as a natural reference system in the thermodynamic perturbation theory of liquids characterized by steep isotropic repulsive interparticle interactions at short distances as well short-range attractive interactions. Weeks et al. (1971); Verlet and Weis (1972) Moreover, the classic hard-sphere model provides a good description of certain classes of colloidal systems. Russel et al. (1989); Rutgers et al. (1996); Dinsmore et al. (1998); Dijkstra et al. (1999); Dijkstra (2014) Note that the hard-core constraint does not uniquely specify the hard-sphere model; there are an infinite number of nonequilibrium hard-sphere ensembles, some of which will be surveyed.

We will also review work that describes generalizations of this simplified hard-sphere model to include its natural extension to packings of hard spheres with a size distribution. This topic has relevance to understanding dispersions of technological importance, (e.g., solid propellants combustionKerstein (1987) flow in packed beds, Scheidegger (1974) and sintering of powdersRahaman (1995)), packings of biological cells, and phase behavior of various molecular systems. For example, the densest packings of hard-sphere mixtures are intimately related to high-pressure phases of compounds for a range of temperatures.Demchyna et al. (2006); Cazorla et al. (2009) Another natural extension of the hard-sphere model that will be surveyed are hard nonspherical particles in two and three dimensions. Asphericity in particle shape is capable of capturing salient features of phases of molecular systems with anisotropic pair interactions (e.g., liquid crystals) and is also a more realistic characteristic of real granular media. In addition, we will review aspects of dense sphere packings in high-dimensional Euclidean spaces, which provide useful physical insights and are relevant to error-correcting codes and information theory.Shannon (1948); Conway and Sloane (1998)

We begin this perspective by introducing relevant definitions and background (Sec. II). This is followed by a survey of work on equilibrium, metastable and nonequilibrium packings of identical spheres and polydisperse spheres in one, two and three dimensions (Sec. III). The geometric-structure approach to jammed and unjammed packings, including their classification via order maps, is emphasized. Subsequently, corresponding results for sphere packings in high Euclidean and non-Euclidean space dimensions, (Sec. IV) and packings of nonspherical particles in low-dimensional Euclidean spaces (Sec. V) are reviewed. Finally, we describe some challenges and open questions for future research (Sec. VI).

Ii Definitions and Background

A packing is a collection of nonoverlapping solid objects or particles of general shape in -dimensional Euclidean space . Packings can be defined in other spaces (e.g., hyperbolic spaces and compact spaces, such as the surface of a -dimensional sphere), but our primary focus in this review is . A saturated packing is one in which there is no space available to add another particle of the same kind to the packing. The packing fraction is the fraction of space covered by the particles.

A -dimensional particle is centrally symmetric if it has a center that bisects every chord through connecting any two boundary points of the particle, i.e., the center is a point of inversion symmetry. Examples of centrally symmetric particles in include spheres, ellipsoids and superballs (defined in Section V). A triangle and tetrahedron are examples of non-centrally symmetric 2D and 3D particles, respectively. A -dimensional centrally symmetric particle for is said to possess equivalent principal axes (orthogonal directions) associated with the moment of inertia tensor if those directions are two-fold rotational symmetry axes. Whereas a -dimensional superball has equivalent directions, a -dimensional ellipsoid generally does not. The reader is referred to Ref. Torquato and Stillinger, 2010 for further discussion concerning particle symmetries.

A lattice in is a subgroup consisting of the integer linear combinations of vectors that constitute a basis for . In the physical sciences, this is referred to as a Bravais lattice. The term “lattice” will refer here to a Bravais lattice only. Every lattice has a dual (or reciprocal) lattice in which the lattice sites are specified by the dual (reciprocal) lattice vectors. Conway and Sloane (1998) A lattice packing is one in which the centroids of the nonoverlapping identical particles are located at the points of , and all particles have a common orientation. The set of lattice packings is a subset of all possible packings in . In a lattice packing, the space can be geometrically divided into identical regions called fundamental cells, each of which has volume and contains the centroid of just one particle of volume . Thus, the lattice packing fraction is


Common -dimensional lattices include the hypercubic , checkerboard , and root lattices; see Ref. Conway and Sloane, 1998. Following Conway and SloaneConway and Sloane (1998), we say two lattices are equivalent or similar if one becomes identical to the other possibly by a rotation, reflection, and change of scale, for which we use the symbol . The and lattices are -dimensional generalizations of the face-centered-cubic (FCC) lattice defined by ; however, for , they are no longer equivalent. In two dimensions, is the triangular lattice. In three dimensions, is the body-centered-cubic (BCC) lattice. In four dimensions, the checkerboard lattice and its dual are equivalent, i.e., . The hypercubic lattice and its dual lattice are equivalent for all .

A periodic packing of congruent particles is obtained by placing a fixed configuration of particles (where ) with arbitrary orientations subject to the nonoverlapping condition in one fundamental cell of a lattice , which is then periodically replicated. Thus, the packing is still periodic under translations by , but the particles can occur anywhere in the chosen fundamental cell subject to the overall nonoverlap condition. The packing fraction of a periodic packing is given by


where is the number density, i.e., the number of particles per unit volume.

For simplicity, consider a packing of identical -dimensional spheres of diameter centered at the positions } in a region of volume in -dimensional Euclidean space . Ultimately, we will pass to the thermodynamic limit, i.e., , such that the number density is a fixed positive constant and its corresponding packing fraction is given by


where is the number density,


is the volume of a -dimensional sphere of radius , and is the gamma function. For an individual sphere, the kissing or contact number is the number of spheres that may simultaneously touch this sphere. In a sphere packing, the mean kissing or contact number per particle, , is the average of over all spheres. For a lattice packings, . For statistically homogeneous sphere packings in , the quantity is proportional to the probability density associated with simultaneously finding sphere centers at locations } in ; see Ref. Hansen and McDonald, 1986 and references therein. With this convention, each -particle correlation function approaches unity when all particles become widely separated from one another for any system without long-range order. Statistical homogeneity implies that is translationally invariant and therefore only depends on the relative displacements of the positions with respect to some arbitrarily chosen origin of the system, i.e., , where .

The pair correlation function is of basic interest in this review. If the system is also rotationally invariant (statistically isotropic), then depends on the radial distance only, i.e., . The total correlation function is defined by . Importantly, we focus in this review on disordered packings in which decays to zero for large sufficiently rapidly so that its volume integral over all space exists.Torquato and Stillinger (2006a)

As usual, we define the nonnegative structure factor for a statistically homogeneous packing as


where is the Fourier transform of and is the wave vector. The nonnegativity of for all follows physically from the fact that it is proportional to the intensity of the scattering of incident radiation on a many-particle system. Hansen and McDonald (1986) The structure factor provides a measure of the density fluctuations in particle configurations at a particular wave vector . For any point configuration in which the minimal pair distance is some positive number, such as a sphere packing, and , the structure factor obeys the following sum rule:Torquato (2018)


For a single-component system in thermal equilibrium at number density and absolute temperature , the structure factor at the origin is directly related to the isothermal compressibility via the relation Hansen and McDonald (1986)


where is the Boltzmann constant.

It is well known from Fourier transform theory that if a real-space radial function in decreases sufficiently rapidly to zero for large such that all of its even moments exist, then its Fourier transform , where is a wavenumber, is an even and analytic function at . Hence, , defined by (5), has an expansion about in any space dimension of the general form


where and are the -dependent constants defined by






is the surface area of -dimensional sphere of radius . This behavior is to be contrasted with that of maximally random jammed sphere packings,Torquato et al. (2000) which possess a structure factor that is nonanalytic at (Ref. Donev et al., 2005a), as discussed in greater detail in Sec. III.5.

A hyperuniform Torquato and Stillinger (2003) point configuration in is one in which the structure factor tends to zero as the wavenumber tends to zero, i.e.,


implying that single scattering of incident radiation at infinite wavelengths is completely suppressed. Equivalently, a hyperuniform system is one in which the number variance of particles within a spherical observation window of radius grows more slowly than the window volume, i.e., , in the large- limit. This class of point configurations includes perfect crystals, many perfect quasicrystals and special disordered many-particle systems. Torquato and Stillinger (2003); Zachary and Torquato (2009); Oğuz et al. (2017) Note that structure-factor definition (5) and the hyperuniformity requirement (12) dictate that the following sum rule involving that a hyperuniform point process must obey:


This sum rules implies that must exhibit negative correlations, i.e., anticorrelations, for some values of .

The hyperuniformity concept was generalized to incorporate two-phase heterogeneous media (e.g., composites, porous media and dispersions). Zachary and Torquato (2009); Torquato (2016) Here the phase volume fraction fluctuates within a spherical window of radius and hence can be characterized by the volume-fraction variance . For typical disordered two-phase media, the variance for large goes to zero like . However, for hyperuniform disordered two-phase media, goes to zero asymptotically more rapidly than the inverse of the window volume, i.e., faster than , which is equivalent to the following condition on the spectral density :Torquato (2016)


The spectral density is proportional to scattered intensity associated with “mass” (volume) content of the phases.Debye et al. (1957)

Iii Sphere Packings in Low Dimensions

The classical statistical mechanics of hard-sphere systems has generated a huge collection of scientific publications, dating back at least to Boltzmann; Boltzmann (1898) see also Refs. Alder and Wainwright, 1957; Reiss et al., 1959; Salsburg and Wood, 1962; Hansen and McDonald, 1986. Here we focus on packings of frictionless, congruent spheres of diameter in one, two and three dimensions in the absence of gravity.

It is important to observe that the impenetrability constraint alone of this idealized hard-sphere model does not uniquely specify the statistical ensemble. Hard-sphere systems can be in thermal equilibrium (as discussed in Sec. III.1) or derived from one of an infinite number of nonequilibrium ensembles Torquato (2002) (see Sec. III.2 for examples).

iii.1 Equilibrium and Metastable Phase Behavior

The phase behavior of hard spheres provides powerful insights into the nature of liquid, crystal and metastable states as well as the associated phase transitions in molecular and colloidal systems. It is well known that the pressure of a stable thermodynamic phase in at packing fraction and temperature is simply related to the contact value of the pair correlation function, : Torquato (2002)


Away from jammed states, it has been proved that the mean nearest-neighbor distance between spheres, , is bounded from above by the pressure,Torquato (1995a) namely, .

Figure 1 schematically shows the 3D phase behavior in the - plane. At sufficiently low densities, an infinitesimally slow compression of the system at constant temperature defines a thermodynamically stable liquid branch for packing fractions up to the “freezing” point (). Increasing the density beyond the freezing point putatively results in an entropy-driven first-order phase transition Alder and Wainwright (1957); Frenkel (1999) to a crystal branch that begins at the melting point (). While there is no rigorous proof that such a first-order freezing transition occurs in three dimensions, there is overwhelming simulational evidence for its existence, beginning with the pioneering work of Alder and Wainwright.Alder and Wainwright (1957) Slow compression of the system along the crystal branch must end at one of the optimal (maximally dense) sphere packings with , each of which is a jammed packing in which each particle contacts 12 others (see Sec. III.4). This equilibrium state has an infinite pressure and is putatively entropically favored by the fcc packing.Mau and Huse (1999)

Figure 1: The isothermal phase behavior of the 3D hard-sphere model in the pressure-packing fraction plane, adapted from Ref. Torquato, 2002. Three different isothermal densification paths by which a hard-sphere liquid may jam are shown. An infinitesimal compression rate of the liquid traces out the thermodynamic equilibrium path (shown in green), including a discontinuity resulting from the first-order freezing transition to a crystal branch that ends at a maximally dense, infinite-pressure jammed state. Rapid compressions of the liquid while suppressing some degree of local order (blue curves) can avoid crystal nucleation (on short time scales) and produce a range of amorphous metastable extensions of the liquid branch that jam only at the their density maxima.

However, compressing a hard-sphere liquid rapidly, under the constraint that significant crystal nucleation is suppressed, can produce a range of metastable branches whose density end points are disordered “jammed” packings,Rintoul and Torquato (1996); Torquato (2002); Jadrich and Schweizer (2013) which can be regarded to be glasses. A rapid compression leads to a lower random jammed density than that for a slow compression. The most rapid compression ending in mechanically stable packing is presumably the maximally random jammed (MRJ) stateTorquato et al. (2000) with . Accurate approximate formulas for the pressure along such metastable extensions up to the jamming points have been obtained. Torquato (1995b, 2002) This ideal amorphous state is described in greater detail in Sec. III.5. Note that rapid compression a hard-sphere system is analogous to supercooling a molecular liquid.

Pair statistics are exactly known only in the case of 1D “hard-rod” systems.Zernike and Prins (1927) For , approximate formulas for are known along liquid branches.Hansen and McDonald (1986) Approximate closures of the Ornstein-Zernike integral equation linking the direct correlation function to the total correlation function Ornstein and Zernike (1914) such as the Percus-Yevick (PY) and hypernetted chain schemes,Percus and Yevick (1958); Stell (1963); Hansen and McDonald (1986) provide reasonably accurate estimates of for hard-spheres liquids. Because decays to unity exponentially fast for liquid states, we can conclude that it must have a corresponding structure factor that is an even function and analytic at ; see Eq. (8). Also, since the leading-order term must be positive because the isothermal compressibility is positive [cf. (7)], hard-sphere liquids are not hyperuniform. Figure 2 shows and the corresponding structure factor in three dimensions at as obtained from the PY approximation.

Figure 2: Pair statistics for an equilibrium hard-sphere fluid in three dimensions at as obtained from the PY approximation. Percus and Yevick (1958); Hansen and McDonald (1986) Top panel: Pair correlation function versus , where is the sphere diameter. Bottom panel: The corresponding structure factor as a function of the dimensionless wavenumber , which clearly shows nonhyperuniformity.

iii.2 Nonequilibrium Disordered Packings

Here we briefly discuss three different nonequilibrium sphere packings: random sequential addition, “ghost” random sequential addition, and random close packing.

iii.2.1 Random Sequential Addition

Perhaps one of the most well-known nonequilibrium packing models is the random sequential addition (RSA) packing procedure, which is a time-dependent process that generates disordered sphere packings in ; see Refs. Reńyi, 1963; Widom, 1966; Finegold and Donnell, 1979; Feder, 1980; Swendsen, 1981; Cooper, 1988; Talbot et al., 1991; Torquato et al., 2006; Zhang and Torquato, 2013. The RSA packing process in the first three space dimensions has been used to model a variety of different condensed phases, including protein adsorption, Feder (1980) polymer oxidation, Flory (1939) particles in cell membranes Finegold and Donnell (1979) and ion implantation in semiconductors. Roman and Majlis (1983)

In its simplest rendition, a RSA sphere packing is produced by randomly, irreversibly, and sequentially placing nonoverlapping objects into a large volume in that at some initial time is empty of spheres. If an attempt to add a sphere at some time results in an overlap with an existing sphere in the packing, the attempt is rejected and further attempts are made until it can be added without overlapping existing spheres. One can stop the addition process at any finite time , obtaining an RSA configuration with a time-dependent packing fraction but this value cannot exceed the maximal saturation packing fraction that occurs in the infinite-time and thermodynamic limits. Even at the saturation state, the spheres are never in contact and hence are not jammed in the sense described in Sec. III.3; moreover, the pair correlation function decays to unity for large super-exponentially fast. Bonnier et al. (1994) The latter property implies that the structure factor of RSA packings in must be an even and analytic function at Torquato et al. (2006) as indicated in Eq. (8). It is notable that saturated RSA packings are not hyperuniform for any finite space dimension. Torquato et al. (2006); Zhang and Torquato (2013)

In the one-dimensional case, also known as the “car-parking” problem, the saturation packing fraction was obtained analytically: .Reńyi (1963) However, for , in the thermodynamic limit has only been estimated via numerical simulations. The most precise numerical study to date Zhang and Torquato (2013) has yielded for and for . Estimates of in higher dimensions are discussed in Sec. IV. RSA packings have also been examined for spheres with a size distribution Adamczyk et al. (1997) and other particle shapes, including squares, Feder (1980); Brosilow et al. (1991); Viot et al. (1992) rectangles, Vigil and Ziff (1989) ellipses, Viot et al. (1992) superdisks, Gromenko and Privman (2009) and polygons, Cieśla and Barbasz (2014); Zhang (2018) in and spheroids, Sherwood (1997) spherocylinders Viot et al. (1992) and cubes Bonnier (2001); Cieśla and Kubala (2018) in .

Figure 3: (Color online) Examples of two nonequilibrium packing models in two dimensions under periodic boundary conditions. Left panel: A configuration of the standard RSA packing at saturation with . Right panel: A configuration of a ghost RSA packing at a packing fraction very near its maximal value of 0.25. Note that the packing is clearly unsaturated (as defined in Sec. II) and there are no contacting particles.

iii.2.2 “Ghost” Random Sequential Addition

There is a generalization of the aforementioned standard RSA process, which can be viewed as a “thinning” process of a Poisson distribution of sphere and is parameterized by the positive constant that lies in the closed interval ; see Ref. Torquato and Stillinger, 2006b. When , one recovers the standard RSA process and when , one obtains the “ghost” random sequential addition (RSA) process that enables one to obtain exactly not only the time-dependent packing fractions but all of the -particle correlation functions for any and . The reader is referred to Ref. Torquato and Stillinger, 2006b for details about the general model.

In the ghost RSA process, one attempts to place spheres into an initially empty region of space randomly, irreversibly and sequentially. However, here one keeps track of any rejected sphere, which is called a “ghost” sphere. No additional spheres can be added whenever they overlap an existing sphere or a ghost sphere. The packing fraction at time for spheres of unit diameter is given by , where is the volume of sphere of radius . Thus, we see that as , , which is appreciably smaller than the RSA saturation packing fraction in low dimensions; see Fig. 3 for 2D examples. Nonetheless, it is notable that the ghost RSA process is the only hard-sphere packing model that is exactly solvable for any dimension and all realizable densities, which has implications for high-dimensional packings, as discussed in Sec. IV.

iii.2.3 Random Close Packing

The “random close packing” (RCP) notion was pioneered by BernalBernal (1960, 1965) to model the structure of liquids, and has been one of the more persistent themes with a venerable history. Scott and Kilgour (1969); Anonymous (1972); Gotoh and Finney (1974); Berryman (1983); Jodrey and Tory (1985); Zinchenko (1994); Jullien et al. (1997); Pouliquen et al. (1997) Two decades ago, the prevailing notion of the RCP state was that it is the maximum density that a large, random collection of congruent (identical) spheres can attain and that this density is a well-defined quantity. This traditional view has been summarized as follows: “Ball bearings and similar objects have been shaken, settled in oil, stuck with paint, kneaded inside rubber balloons–and all with no better result than (a packing fraction of) ”; see Ref. Anonymous, 1972.

Torquato, Truskett and Debenedetti Torquato et al. (2000) have argued that this RCP-state concept is actually ill-defined because “randomness” and “closed-packed” were never defined and, even if they were, are at odds with one another. Using the Lubachevsky-Stillinger (LS) Lubachevsky and Stillinger (1990) molecular-dynamics growth algorithm to generate jammed packings, it was shown Torquato et al. (2000) that fastest particle growth rates generated the most disordered sphere packings with , but that by slowing the growth rates larger packing fractions could be continuously achieved up to the densest value such that the degree of order increased monotonically with . These results demonstrated that once can increase the packing fraction by an arbitrarily small amount at the expense of correspondingly small increases in order and thus the notion of RCP is ill-defined as the highest possible density that a random sphere packing can attain. To remedy these flaws, Torquato, Truskett and Debenedetti Torquato et al. (2000) replaced the notion of “close packing” with “jamming” categories (defined precisely in Sec. III.3) and introduced the notion of an “order metric” to quantify the degree of order (or disorder) of a single packing configuration. This led them to supplant the concept of RCP with the maximally random jammed (MRJ) state, which is defined, roughly speaking, to be that jammed state with a minimal value of an order metric (see Sec. III.3.4 for details). This work pointed the way toward a quantitative means of characterizing all packings, namely, the geometric-structure approach.


Figure 4: Typical protocols, used to generate disordered sphere packings in three dimensions, produce highly crystalline packings in two dimensions. Left panel: A highly ordered collectively jammed configuration (Sec. III.3.1) of 1000 disks with produced using the Lubachevsky-Stillinger (LS) algorithm Lubachevsky and Stillinger (1990) with a fast expansion rate. Donev et al. (2004a) Right panel: A 3D MRJ-like configuration of 500 spheres with produced using the LS algorithm with a fast expansion rate. Torquato et al. (2000)

We note that whereas the LS packing protocol with a fast growth rate typically leads to disordered jammed states in three dimensions, it invariably produces highly crystalline “collectively” jammed packings in two dimensions. Figure 4 vividly illustrates the differences between the textures produced in three and in two dimensions (see Section VII for further remarks).

iii.3 Geometric-Structure Approach to Jammed Packings

A “jammed” packing is one in which each particle is in contact with its nearest neighbors such that mechanical stability of a specific type is conferred to the packing, as detailed below. Two conceptual approaches for their study have emerged. One is the “ensemble” approach, Bernal (1960, 1965); Edwards (1994); Edwards and Grinev (2001); Liu and Nagel (1998); Makse and Kurchan (2002); Silbert et al. (2002); O’Hern et al. (2003); Wyart et al. (2005); Silbert et al. (2005); Gao et al. (2006); Song et al. (2008); Parisi and Zamponi (2010) which for a given packing protocol aims to understand typical configurations and their frequency of occurrence. The other, more recently, is the “geometric-structure” approach, Torquato et al. (2000); Torquato and Stillinger (2001); Kansal et al. (2002a); Torquato et al. (2003); Donev et al. (2004a, b, 2007a); Torquato and Stillinger (2007) which emphasizes quantitative characterization of single-packing configurations without regard to their occurrence frequency in the protocol used to produce them. Here we focus on the latter approach, which enables one to enumerate and classify packings with a diversity of order/disorder and packings fractions, including extremal packings, such as the densest sphere packing and MRJ packings.

iii.3.1 Jamming Categories

Three broad and mathematically precise “jamming” categories of sphere packings can be distinguished depending on the nature of their mechanical stability. Torquato and Stillinger (2001, 2003) In order of increasing stringency (stability), for a finite sphere packing, these are the following:

  1. Local jamming: Each particle in the packing is locally trapped by its neighbors (at least contacting particles, not all in the same hemisphere), i.e., it cannot be translated while fixing the positions of all other particles;

  2. Collective jamming: Any locally jammed configuration is collectively jammed if no subset of particles can be simultaneously displaced so that its members move out of contact with one another and with the remainder set; and

  3. Strict jamming: Any collectively jammed configuration that disallows all uniform volume-nonincreasing strains of the system boundary is strictly jammed.

These hierarchical jamming categories are closely related to the concepts of “rigid” and “stable” packings found in the mathematics literature, Connelly et al. (1998) and imply that there can be no “rattlers” (i.e., movable but caged particles) in the packing. The jamming category of a given packing depends on the boundary conditions employed; see Refs. Torquato and Stillinger, 2001, 2010 for specific examples. Rigorous and efficient linear-programming (LP) algorithms have been devised to assess whether a particular sphere packing is locally, collectively, or strictly jammed. Donev et al. (2004a, c) These jamming categories can now be ascertained in real-system experiments by applying the LP tests to configurational coordinates of a packing determined via a variety of imaging techniques, including tomography Aste et al. (2006), confocal microscopy Brujić et al. (2003) and magnetic resonance imaging.  Man et al. (2005)

iii.3.2 Polytope Picture and Pressure in Jamming Limit

A packing of hard spheres of diameter in a jammed framework in is specified by a -dimensional configurational position vector . Isostatic jammed packings possess the minimal number of contacts for a jamming category and boundary conditions. Donev et al. (2005b) The relative differences between isostatic collective and strict jammed packings diminish as becomes large, and since the number of degrees of freedom is essentially equal to , an isostatic packing has a mean contact number per particle, , equal to in the large- limit. Packings having more and fewer contacts than isostatic ones are hyperstatic and hypostatic, respectively. Sphere packings that are hypostatic cannot be collectively or strictly jammed.Donev et al. (2007a)

Consider decreasing the density slightly in a sphere packing that is at least collectively jammed by reducing the particle diameter by , so that the packing fraction is lowered to , where is called the jamming gap or distance to jamming. There is a sufficiently small that does not destroy the jamming confinement property. For fixed and sufficiently small , it can be shown asymptotically, through first order in ), that the set of displacements accessible to the packing approaches a convex limiting polytope (high-dimensional polyhedron) . Salsburg and Wood (1962); Stillinger and Salsburg (1969) This polytope is determined from the linearized impenetrability equations Donev et al. (2004a, c) and is necessarily bounded for a jammed configuration.

Now consider adding thermal kinetic energy to a nearly jammed sphere packing in the absence of rattlers. While the system will not be globally ergodic over the full system configuration space and thus not in thermodynamic equilibrium, one can still define a macroscopic pressure for the trapped but locally ergodic system by considering time averages as the system executes a tightly confined motion around the particular jammed configuration . The probability distribution of the time-averaged interparticle forces has been rigorously linked to the contact value of the pair correlation function . Donev et al. (2005b) Moreover, the available (free) configuration volume scales with the jamming gap such that the reduced pressure is asymptotically given by the free-volume equation of state: Salsburg and Wood (1962); Stillinger and Salsburg (1969); Donev et al. (2005b)


where is the absolute temperature and is the number density. Relation (16) is remarkable, since it enables one to determine accurately the true jamming density of a given packing, even if the actual jamming point has not quite yet been reached, just by measuring the pressure and extrapolating to . This free-volume form has been used to estimate the equation of state along “metastable” extensions of the hard-sphere fluid up to the infinite-pressure endpoint, assumed to be disordered jammed states (see Fig. 1 and Sec. III.1). Torquato (1995b, 2002) Kamien and Liu Kamien and Liu (2007) assumed the same free-volume form to fit data for the pressure of “metastable” hard-sphere states.

iii.3.3 Hard-Particle Jamming Algorithms

For many years, the Lubachevsky-Stillinger (LS) algorithm Lubachevsky and Stillinger (1990) has served as the premier numerical scheme to generate a wide spectrum of dense jammed hard-sphere packings with variable disorder in both two and three dimensions. Torquato et al. (2000); Truskett et al. (2000); Kansal et al. (2002a) This is an event-driven or collision-driven molecular-dynamics algorithm: an initial configuration of spheres of a given size within a periodic cell are given initial random velocities and the motion of the spheres are followed as they collide elastically and also grow uniformly until the spheres can no longer expand. The final density is generally inversely related to the particle growth rate. This algorithm has been generalized by Donev, Torquato, and Stillinger Donev et al. (2005c, d) to generate jammed packings of smoothly shaped nonspherical particles, including ellipsoids, Donev et al. (2005d, 2007a) superdisks, Jiao et al. (2008), and superballs. Jiao et al. (2009) Event-driven packing protocols with growing particles, however, do not guarantee jamming of the final packing configuration, since jamming is not explicitly incorporated as a termination criterion.

The task of generating dense packings of particles of general shape within an adaptive periodic fundamental cell has been posed by Torquato and Jiao Torquato and Jiao (2009a, b) as an optimization problem called the adaptive-shrinking cell (ASC) scheme. The negative of the packing fraction, (which can be viewed as an “energy”) is to be minimized subject to constraints, and the design variable are the shape and size of the deformable periodic fundamental cell as well as the centroids and orientations of the particles. The so-called Torquato-Jiao (TJ) sphere-packing algorithm Torquato and Jiao (2010a) is a sequential linear-programming (SLP) solution of the ASC optimization problem for the special case of packings of spheres with a size distribution for which linearization of the design variables is exact. The deterministic SLP solution in principle always leads to strictly jammed packings up to a high numerical tolerance with a wider range of densities and degree of disorder than packings produced by the LS algorithm. From an initial configuration, the TJ algorithm leads to a mechanically stable local “energy” minimum (local density maximum), which in principle is the inherent structure associated with the starting initial many-particle configuration; see Fig. 5. A broad range of inherent structures can be obtained across space dimensions, including locally maximally dense and mechanically stable packings, such as MRJ states,Torquato and Jiao (2010a); Atkinson et al. (2013, 2014, 2016a) disordered hyperstatic packings Jiao et al. (2011) and the globally maximally dense inherent structures,Torquato and Jiao (2010a); Hopkins et al. (2011a, 2012a); Marcotte and Torquato (2013); Kallus et al. (2013) with very small computational cost, provided that the system sizes are not too large.

It is notable the TJ algorithm creates disordered jammed sphere packings that are closer to the ideal MRJ state than previous algorithms.Atkinson et al. (2013) It was shown that the rattler concentration of these packings converges towards 1.5% in the infinite-system limit, which is markedly lower than previous estimates for the MRJ state using the LS protocol (with about 3% rattlers).

Figure 5: (Color online). A schematic diagram of inherent structures (local density maxima) for sphere packings of spheres, as taken from Ref. Torquato and Jiao, 2010a. The horizontal axis labeled stands for the entire set of centroid positions, and (“energy”) decreases downward. The jagged curve is the boundary between accessible configurations (yellow) and inaccessible ones (blue). The deepest point of the accessible configurations corresponds to the maximal density packings of hard spheres.

iii.3.4 Order Metrics

The enumeration and classification of both ordered and disordered sphere packings is an outstanding problem. Since the difficulty of the complete enumeration of packing configurations rises exponentially with the number of particles, it is desirable to devise a small set of intensive parameters that can characterize packings well. One obvious property of a sphere packing is the packing fraction . Another important characteristic of a packing is some measure of its “randomness” or degree of disorder. Devising such measures is a highly nontrivial challenge, but even the tentative solutions that have been put forth during the last two decades have been fruitfully applied to characterize not only sphere packings Torquato et al. (2000); Truskett et al. (2000); Kansal et al. (2002a); Torquato and Stillinger (2003) but also simple liquids,Rintoul et al. (1996); Rintoul and Torquato (1996); Truskett et al. (2000); Errington et al. (2003) glasses,Truskett et al. (2000); Zhang et al. (2016) water, Errington and Debenedetti (2001); Errington et al. (2002) disordered ground states, Torquato et al. (2015) random media, DiStasio et al. (2018) and biological systems.Jiao et al. (2014)

It is quite reasonable to consider “entropic” measures to characterize the randomness of packings. However, as pointed out by Kansal et al. Kansal et al. (2002a), a substantial hurdle to overcome in implementing such an order metric is the necessity to generate all possible jammed states or, at least, a representative sample of such states in an unbiased fashion using a “universal” protocol in the large-system limit, each of which is an intractable problem. Even if such a universal protocol could be developed, however, the issue of what weights to assign the resulting configurations remains. Moreover, there are other basic problems with the use of entropic measures as order metrics, as we will discuss in Sec. III.5, including its significance for certain 2D hard-disk packings.

We know that a many-body system of particles is completely characterized statistically by its -body probability density function that is associated with finding the -particle system with configuration at some time . Such complete information is virtually never available for large and, in practice, one must settle for reduced information, such as a scalar order metric . Any order metric conventionally possesses the following three properties: (1) it is a well-defined scalar function of a configuration ; (2) it is subject typically to the normalization ; and, (3) for any two configurations and , implies that configuration is to be considered as more ordered than configuration . The set of order metrics that one selects is unavoidably subjective, given that there appears to be no single universally applicable scalar metric capable of describing order across all lengths scales. However, one can construct order metrics that lead to consistent results, as we will be discussed below.

Many relevant order metrics have been devised. While a comprehensive discussion of this topic is beyond the scope of this article, it is useful to here to note some order metrics that have been identified, including bond-orientational order metrics in two Kansal et al. (2000) and three dimensions, Steinhardt et al. (1983) translational order metrics, Torquato et al. (2000); Truskett et al. (2000); Torquato et al. (2015) and hyperuniformity order metrics. Torquato and Stillinger (2003); Zachary and Torquato (2009) These specific order metrics have both strengths and weaknesses. This raises the question of what are the characteristics of a good order metric? It has been suggested that a good scalar order metric should have the following additional properties: Kansal et al. (2002a); Torquato and Stillinger (2010) (1) sensitivity to any type of ordering without bias toward any reference system; (2) ability to reflect the hierarchy of ordering between prototypical systems given by common physical intuition (e.g., perfect crystals with high symmetry should be highly ordered, followed by quasicrystals, correlated disordered packings without long-range order, and finally spatially uncorrelated or Poisson distributed particles); (3) incorporation of both the variety of local coordination patterns as well as the spatial distribution of such patterns should be included; and (4) the capacity to detect translational and orientational order at any length scale. Moreover, any useful set of order metrics should consistently produce results that are positively correlated with one another. Torquato et al. (2000); Torquato (2002) The development of improved order metrics deserves continued research attention.

Order metrics and maps have been fruitfully extended to characterize the degree of structural order in condensed phases in which the constituent particles (jammed or not) possess both attractive and repulsive interactions. This includes the determination of the order maps of models of simple liquids, glasses and crystals with isotropic interactions, Rintoul et al. (1996); Rintoul and Torquato (1996); Truskett et al. (2000); Errington et al. (2003) models of water, Errington and Debenedetti (2001); Errington et al. (2002) disordered ground states of long-ranged isotropic pair potentials,Torquato et al. (2015) and models of amorphous polymers. Stachurski (2003)

iii.4 Order Maps and Extremal Packings

The geometric-structure classification naturally emphasizes that there is a great diversity in the types of attainable packings with varying magnitudes of overall order, packing fraction , mean contact number per particle, , among many other intensive parameters. Any attainable hard-sphere configuration, jammed or not, can be mapped to a point in this multidimensional parameter space that we call an order map. The use of “order maps” in combination with the mathematically precise “jamming categories” enable one to view and characterize individual packings, including the densest sphere packing (e.g., Kepler’s conjecture in 3D) and MRJ packings as extremal states in the order map for a given jamming category. Indeed, this picture encompasses not only these special jammed states, but an uncountably infinite number of other packings, some of which have only recently been identified as physically significant, e.g., the jamming-threshold states Torquato and Stillinger (2007) (least dense jammed packings) as well as states between these and MRJ.

Figure 6: Schematic order map of sphere packings in in the density-order (-) plane. White and blue regions contain the attainable packings, blue regions represent the jammed subspaces, and dark shaded regions contain no packings. The boundaries of the jammed region delineate extremal structures. The locus of points - correspond to the lowest-density jammed packings. The locus of points - correspond to the densest jammed packings. Points MRJ represent the maximally random jammed states, i.e., the most disordered states subject to the jamming constraint.


Figure 7: Three different extremal strictly jammed packings in identified in Fig. 6, as taken from Ref. Torquato and Stillinger, 2010. Left panel: (A) A tunneled crystal with . Middle panel: MRJ packing with (isostatic). Right panel: (B) FCC packing with .

The simplest order map is the one that maps any hard-sphere configuration to a point in the - plane. This two-parameter description is but a very small subset of the relevant parameters that are necessary to fully characterize a configuration, but it nonetheless enables one to draw important conclusions. Figure 6 shows such an order map that delineates the set of strictly jammed packings Torquato et al. (2000); Kansal et al. (2002a); Donev et al. (2004a); Torquato and Stillinger (2007) from non-jammed packings in three dimensions. The boundaries of the jammed region delineate extremal structures (see Fig. 7). The densest sphere packings, Hales (2005) which lie along the locus - with , are strictly jammed. Torquato and Stillinger (2001); Donev et al. (2004a) Point B corresponds to the face-centered-cubic (fcc) packing, i.e., it is the most ordered and symmetric densest packing, implying that their shear moduli are infinitely large.Torquato et al. (2003) The most disordered subset of the stacking variants of the fcc packing is denoted by point . In two dimensions, the strictly jammed triangular lattice is the unique densest packing Fejes Tóth (1964) and so the line - in collapses to a single point in .

The least dense strictly jammed packings are conjectured to be the “tunneled crystals” in with , corresponding to the locus of points -.Torquato and Stillinger (2007) These infinitely-degenerate sparse structures were analytically determined by appropriate stackings of planar “honeycomb” layers of spheres. These constitute a set of zero measure among the possible packings with the same density and thus are virtually impossible to discover using packing algorithms, illustrating the importance of analyzing individual packings, regardless of their frequency of occurrence in the space of jammed configurations, as advocated by the geometric-structure approach. The tunneled crystals are subsets of the densest packings but are permeated with infinitely long chains of particle vacancies. Torquato and Stillinger (2007) Every sphere in a tunneled crystal contacts 7 immediate neighbors. Interestingly, the tunneled crystals exist at the edge of mechanical stability, since removal of any one sphere from the interior would cause the entire packing to collapse. The tunneled crystals are magnetically frustrated structures and Burnell and SondhiBurnell and Sondhi (2008) showed that their underlying topologies greatly simplify the determination of their antiferromagnetic properties. In , the “reinforced” kagomé packing with exactly four contacts per particle (isostatic value) is evidently the lowest density strictly jammed packing Donev et al. (2004a) with and so the line - in collapses to a single point in .

iii.5 Maximally Random Jammed States

Among all strictly jammed sphere packings in , the one that exhibits maximal disorder (minimizes some given order metric ) is of special interest. This is called the maximally random jammed (MRJ) state Torquato et al. (2000); see Fig. 6. The MRJ state is automatically compromised by passing either to the maximal packing fraction (fcc and its stacking variants) or the minimal possible density for strict jamming (tunneled crystals), thereby causing any reasonable order metric to rise on either side. A variety of sensible order metrics produce a 3D MRJ state with a packing fraction (see Ref. Kansal et al., 2002a), close to the traditionally advocated density of the RCP state, and with an isostatic mean contact number (see Ref. Donev et al., 2005b). This consistency among the different order metrics speaks to their utility, even if a perfect order metric has not yet been identified. However, the density of the MRJ state is not sufficient to fully specify it. It is possible to have a rather ordered strictly jammed packing at this very same density, Kansal et al. (2002a) as indicated in Fig. 6. The MRJ state refers to a single packing that is maximally disordered subject to the strict jamming constraint, regardless of its probability of occurrence in some packing protocol. Thus, the MRJ state is conceptually and quantitatively different from random close packed (RCP) packings Bernal (1960, 1965), which, more recently, have been suggested to be the most probable jammed configurations within an ensemble.O’Hern et al. (2003) The differences between these states are even starker in two dimensions, e.g., MRJ packings of identical circular disks in have been shown to be dramatically different from their RCP counterparts, including their respective densities, average contact numbers, and degree of order, Atkinson et al. (2014) as detailed below.

Figure 8: Pair statistics for packings of spheres of diameter in the immediate neighborhood of the 3D MRJ state with . Top panel: The pair correlation function versus . Donev et al. (2005b) The split second peak, the discontinuity at twice the sphere diameter, and the divergence near contact are clearly visible. Bottom panel: The corresponding structure factor as a function of the dimensionless wavenumber for a million-particle packing.Donev et al. (2005a) The inset shows that tends to zero (within numerical error) linearly in hence an MRJ packing is effectively hyperuniform.

The MRJ state under the strict-jamming constraint is a prototypical glass Torquato and Stillinger (2007) in that it is maximally disordered (according to a variety of order metrics) without any long-range order (Bragg peaks) and perfectly rigid (i.e., the elastic moduli are unbounded. Torquato et al. (2003); Torquato and Stillinger (2010)) Its pair correlation function can be decomposed into a Dirac-delta-function contribution from particle contacts and a continuous-function contribution :


where is given by (11), is a radial Dirac delta function and . The corresponding structure factor in the long-wavenumber limit is


For , possesses the well-known feature of a split second peak,Zallen (1983) with a prominent discontinuity at twice the sphere diameter, as shown in Fig. 8. Interestingly, an integrable power-law divergence ( with ) exists for near contacts. Donev et al. (2005b, a) Moreover, an MRJ packing in has a structure factor that tends to zero linearly in (within numerical error) in the limit and hence is hyperuniform (see Fig. 8), belonging to the same class as Fermi-sphere point processes Torquato et al. (2008) and superfluid helium in its ground state. Feynman and Cohen (1956); Reatto and Chester (1967) This nonanalytic behavior at implies that MRJ packings are characterized by negative “quasi-long-range” (QLR) pair correlations in which the total correlation function decays to zero asymptotically like ; see Refs. Donev et al., 2005a; Skoge et al., 2006; Hopkins et al., 2012b; Atkinson et al., 2016a, b; Torquato, 2018. The QLR hyperuniform behavior distinguishes the MRJ state strongly from that of the equilibrium hard-sphere fluid Haberko et al. (2013), which possesses a structure factor that is analytic at [cf. (8)], and thus its decays to zero exponentially fast for large ; see Ref. Torquato and Stillinger, 2010. Interestingly, the hyperuniformity concept enables one to identify a static nonequilibrium growing length scale in overcompressed (rapidly compressed) hard-sphere systems as a function of up to the jammed glassy state.Hopkins et al. (2012b); Atkinson et al. (2016b) This led to the identification of static nonequilibrium growing length scales in supercooled liquids on approaching the glass transition.Marcotte et al. (2013)

The following is a general conjecture due to Torquato and Stillinger Torquato and Stillinger (2003) concerning the conditions under which certain strictly jammed packings are hyperuniform:

Conjecture: Any strictly jammed saturated infinite packing of identical spheres in is hyperuniform.

To date, there is no rigorously known counterexample to this conjecture. Its justification and rigorously known examples are discussed in Refs. Torquato, 2018 and Atkinson et al., 2016a.

A recent numerical study of necessarily finite disordered packings calls into question the connection between jamming and perfect hyperuniformity.Ikeda and Berthier (2015) It is problematic to try to test this conjecture using current packing protocols to construct supposedly disordered jammed states for a variety of reasons. First, we stress that the conjecture eliminates packings that may have a rigid backbone but possess “rattlers” - – a subtle point that has not been fully appreciated in numerical investigations. Ikeda and Berthier (2015); Wu et al. (2015); Ikeda et al. (2017) Current numerically generated disordered packings that are putatively jammed tend to contain a small concentration of rattlers;Lubachevsky et al. (1991); O’Hern et al. (2003); Donev et al. (2005a); Torquato and Jiao (2010a); Atkinson et al. (2013) because of these movable particles, the entire packing cannot be deemed to be “jammed.” Moreover, it has been recently shown that various standard packing protocols struggle to reliably create packings that are jammed for even modest system sizes. Atkinson et al. (2016a) Although these packings appear to be jammed by conventional tests, rigorous linear-programming jamming tests Donev et al. (2004a, c) reveal that they are not. Recent evidence has emerged that suggest that deviations from hyperuniformity in putative MRJ packings also can in part be explained by a shortcoming of the numerical protocols to generate exactly-jammed configurations as a result of a type of “critical slowing down” Atkinson et al. (2016a) as the packing’s collective rearrangements in configuration space become locally confined by high-dimensional “bottlenecks” through which escape is a rare event. Thus, a critical slowing down implies that it becomes increasingly difficult numerically to drive the value of exactly down to its minimum value of zero if a true jammed critical state could be attained; typically,Donev et al. (2005a) . Moreover, the inevitable presence of even a small fraction of rattlers generated by current packing algorithms destroys perfect hyperuniformity. In summary, the difficulty of ensuring jamming in disordered packings as the particle number becomes sufficiently large to access the small-wavenumber regime in the structure factor as well as the presence of rattlers that degrade hyperuniformity makes it extremely challenging impossible to test the Torquato-Stillinger jamming-hyperuniformity conjecture on disordered jammed packings via current numerical protocols. This raises the possibility that the ideal MRJ packings have no rattlers, and provides a challenge to develop packing algorithms that produce large disordered strictly jammed packings that are rattler-free.

A variety of different attributes of MRJ packings generated via the TJ packing algorithm have been investigated in separate studies. In the first paper of a three-part series, Klatt and Torquato Klatt and Torquato (2014) ascertained Minkowski correlation functions associated with the Voronoi cells of such MRJ packings and found that they exhibited even stronger anti-correlations than those shown in the standard pair-correlation function. Klatt and Torquato (2014) In the second paper of this series,Klatt and Torquato (2016) a variety of different correlation functions that arise in rigorous expressions for the effective physical properties of MRJ sphere packings were computed. In the third paper of this series,Klatt and Torquato (2018) these structural descriptors were used to estimate effective transport and electromagnetic properties of composites composed of MRJ sphere packings dispersed throughout a matrix phase and showed, among other things, that electromagnetic waves in the long-wavelength limit can propagate without dissipation, as generally predicted in Ref. Rechtsman and Torquato, 2008. In a separate study, Ziff and Torquato Ziff and Torquato (2017) determined the site and bond percolation thresholds of TJ generated MRJ packings to be and , respectively.

Not surprisingly, ensemble methods that produce “most probable” configurations typically miss interesting extremal points in the order map, such as the locus of points - and the rest of the jamming-region boundary. However, numerical protocols can be devised to yield unusual extremal jammed states. For example, disordered jammed packings can be created in the entire non-trivial range of packing fraction Torquato et al. (2000); Kansal et al. (2002a); Jiao et al. (2011) Thus, in Fig. 6, the locus of points along the boundary of the jammed set to the left and right of the MRJ state are achievable. The TJ algorithm Torquato and Jiao (2010a) was applied to yield disordered strictly jammed packingsJiao et al. (2011) with as low as , which are overconstrained with , and hence are more ordered than the MRJ state.

It is well known that lack of “frustration” Jullien et al. (1997); Torquato (2002) in 2D analogs of 3D computational and experimental protocols that lead to putative RCP states result in packings of identical disks that are highly crystalline, forming rather large triangular coordination domains (grains). Such a 1000-particle packing with is depicted in the right panel of Fig. 4, and is only collectively jammed at this high density. Because such highly ordered packings are the most probable outcomes for these typical protocols, “entropic measures” of disorder would identify these as the most disordered, a misleading conclusion. Recently, Atkinson et al.Atkinson et al. (2014) applied the TJ algorithm to generate relatively large MRJ disk packings with exactly isostatic () jammed backbones and a packing fraction (including rattlers) of . Uncovering such disordered jammed packings of identical hard disks challenges the traditional notion that the most probable distribution is necessarily correlated with randomness and hence the RCP state. An analytical formula was derived for MRJ packing fractions of more general 2D packings,Tian et al. (2015) yielding the prediction in the monodisperse-disk limit, which is in very good agreement with the aforementioned recent numerical estimate.Atkinson et al. (2014)

iii.6 Packings of Spheres With a Size Distribution

Sphere packings with a size distribution (polydispersity), sometimes called hard-sphere mixtures, exhibit intriguing structural features, some of which are only beginning to be understood. Our knowledge of sphere packings with a size distribution is very limited due in part to the infinite-dimensional parameter space, i.e., all particle size ratios and relative concentrations. It is known, for example, that a relatively small degree of polydispersity can suppress the disorder-order phase transition seen in monodisperse hard-sphere systems (see Fig. 1).Henderson et al. (1996) Size polydispersity constitutes a fundamental feature of the microstructure of a wide class of dispersions of technological importance, including those involved in composite solid propellant combustion, Kerstein (1987) flow in packed beds, Scheidegger (1974) sintering of powders,Rahaman (1995) colloids, Russel et al. (1989) transport and mechanical properties of particulate composite materials. Torquato (2002) Packings of biological cells in tissue are better modeled by assuming that the spheres have a size distribution.Jiao et al. (2014); Chen et al. (2016)

Generally, we allow the spheres to possess a continuous distribution in radius , which is characterized by a probability density function . The average of any function is defined by . The overall packing fraction is then defined as


where is the total number density, is given by (4). Two continuous probability densities that have been widely used are the Schulz Schulz (1939) and log-normal Cramer (1954) distributions. One can obtain corresponding results for spheres with discrete different sizes from Eq. (19) by letting


where and are number density and radius of type- particles, respectively, and is the total number density. Therefore, the overall volume fraction using (19) is given by , where is the packing fraction of the th component.

iii.6.1 Equilibrium and Metastable Behavior

The problem of determining the equilibrium phase behavior of hard sphere mixtures is substantially more challenging and richer than that for identical hard spheres, including the possibilities of metastable or stable fluid-fluid and/or solid-solid phase transitions (apart from stable fluid or crystal phases). While much research remains to be done, considerable progress has been made over the years,Lebowitz (1964); Widom and Rowlinson (1970); Lebowitz and Zomick (1971); Mansoori et al. (1971); Salacuse and Stell (1982); Nixon and Silbert (1984); Rosenfeld (1989); Dickman and Stell (1995); Shew and Yethiraj (1996); Lomba et al. (1996); Santos et al. (1999); Dijkstra et al. (1999); Góźdź (2003); Roth (2010); Dijkstra (2014) which we only briefly touch upon here for both additive and nonadditive cases. In additive hard-sphere mixtures, the distance of closest approach between the centers of any two spheres is the arithmetic mean of their diameters. By contrast, in nonadditive hard-sphere mixtures, the distance of closest approach between any two spheres is generally no longer the arithmetic mean.

Additive mixtures have been studied both theoretically and computationally. Lebowitz exactly obtained the pair correlation functions of such systems with components within the Percus-Yevick approximation. Lebowitz (1964) Accurate approximate equations of state under liquid conditions have been found for both discrete Mansoori et al. (1971); Santos et al. (1999) and continuous Salacuse and Stell (1982); Santos et al. (2017) size distributions with additive hard cores. Fundamental-measure theory can provide useful insights about the phase behavior of hard-sphere mixtures.Rosenfeld (1989); Roth (2010) This theory predicts the existence of stable fluid-fluid coexistence for sufficiently large size ratios in additive binary hard-sphere systems, while numerical simulations Dijkstra et al. (1999) indicate that such phase separation is metastable with respect to fluid-crystal coexistence and also shows stable solid-solid coexistence. Nonetheless, it remains an open theoretical question whether a fluid-fluid demixing transition exists in a binary mixture of additive hard spheres for any size ratio and relative composition. In order to quantify fluid-crystal phase coexistence, one must know the densest crystal structure, which is highly nontrivial, especially for large size ratios, although recent progress has been made; see Sec. III.6.2.

The Widom-Rowlinson model is an extreme case of nonadditivity in which like species do not interact and unlike species interact via a hard-core repulsion.Widom and Rowlinson (1970) As the density is increased, this model exhibits a fluid-fluid demixing transition in low dimensions and possesses a critical point that is in the Ising universality class.Dickman and Stell (1995); Shew and Yethiraj (1996) More general nonadditive hard-sphere mixtures in which all spheres interact have been studied. Mixtures of hard spheres with positive nonadditivity (unlike-particle diameters greater than the arithmetic mean of the corresponding like-particle diameters) can exhibit fluid-fluid demixing transition Lomba et al. (1996); Dijkstra et al. (1999) that belongs to the Ising universality class,Góźdź (2003) while those with negative nonadditivity have tendencies to form alloyed (mixed) fluid phases.Nixon and Silbert (1984) For certain binary mixtures of hard sphere fluids with nonadditive diameters, the pair correlation function has been determined in the Percus-Yevick approximation.Lebowitz and Zomick (1971); Nixon and Silbert (1984)

Fluid phases of hard sphere mixtures, like their monodisperse counterparts, are not hyperuniform. However, multicomponent equilibrium plasmas consisting of nonadditive hard spheres with Coulombic interactions enable one to generate a very wide class of disordered hyperuniform as well as “multihyperuniform” Jiao et al. (2014) systems at positive temperatures. Lomba et al. (2017, 2018)

iii.6.2 Maximally Random Jammed States

The study of dense disordered packings of 3D polydisperse additive spheres has received considerable attention. Visscher and Bolsterli (1972); Clarke and Wiley (1987); Yang et al. (1996); Richard et al. (1998); Santiso and Müller (2002); de Lange Kristiansen et al. (2005) However, these investigations did not consider their mechanical stability via our modern understanding of jamming. Not surprisingly, the determination of the MRJ state for an arbitrary polydisperse sphere packing is a challenging open question, but some progress has been made recently, as described below.

Jammed states of polydisperse spheres, whether disordered or ordered, will generally have higher packing fractions when “alloyed” than their monodisperse counterparts. The first investigation that attempted to generate 3D amorphous jammed sphere packings with a polydispersity in size was carried out by Kansal et al. Kansal et al. (2002b) using the LS packing algorithm. It was applied to show that disordered binary jammed packings with a small-to-large size ratio and relative concentration can be obtained whose packing fractions exceeds 0.64 and indeed can attain for (the smallest value considered). Chaudhuri et al. Chaudhuri et al. (2010) numerically generated amorphous 50-50 binary packings with packing fractions in the range for . It is notable Clusel et al. Clusel et al. (2009) carried out a series of experiments to visualize and characterize 3D random packings of frictionless polydisperse emulsion droplets using confocal microscopy.

Until recently, packing protocols that have attempted to produce disordered binary sphere packings have been limited in producing mechanically stable, isostatic packings across a broad spectrum of packing fractions. Many previous simulation studies of disordered binary sphere packings have produced packings that are not mechanically stable Clarke and Wiley (1987); Sakai et al. (2002); de Lange Kristiansen et al. (2005) and report coordination numbers as opposed to contact numbers. Whereas a coordination number indicates only proximity of two spheres, a contact number reflects mechanical stability and is derived from a jammed network.Donev et al. (2004c); Torquato and Stillinger (2010)

Using an the TJ linear programming packing algorithm,Torquato and Jiao (2010a) Hopkins et al. Hopkins et al. (2013) were recently able to generate high-fidelity strictly jammed, isostatic, disordered binary packings with an anomalously large range of packing fractions, , with the size ratio restriction . These packings are MRJ like due to the nature of the TJ algorithm and the use of RSA initial conditions. Additionally, the packing fractions for certain values of and approach those of the corresponding densest known ordered packings.Hopkins et al. (2011a, 2012a) These findings suggest that these high-density disordered packings should be good glass formers for entropic reasons and hence may be easy to prepare experimentally. The identification and explicit construction of binary packings with such high packing fractions could have important practical implications for the design of improved solid propellants, concrete, and ceramics. In this connection, a recent study of MRJ binary sphere packings under confinement is particularly relevant.Chen and Torquato (2015)

The LS algorithm has been used successfully to generate disordered strictly jammed packings of binary disksDonev et al. (2006) with and . By explicitly constructing an exponential number of jammed packings of binary disks with densities spanning the spectrum from the accepted amorphous glassy state to the phase-separated crystal, it has been argued Donev et al. (2006, 2007b) that there is no “ideal glass transition”. Parisi and Zamponi (2005) The existence of an ideal glass transition remains a hotly debated topic of research.

Simulational Zachary et al. (2011a, b); Berthier et al. (2011) as well as experimental Dreyfus et al. (2015); Ricouvier et al. (2017) studies of disordered jammed spheres with a size distribution reveal that they are effectively hyperuniform. In such cases, it has been demonstrated that the proper means of investigating hyperuniformity is through the spectral density ,Zachary et al. (2011a, b) which is defined by condition (14).

iii.6.3 Densest Packings

The densest packings of spheres with a size distribution is of great interest in crystallography, chemistry and materials science. It is notable that the densest packings of hard-sphere mixtures are intimately related to high-pressure phases of molecular systems, including intermetallic compoundsDemchyna et al. (2006) and solid rare-gas compoundsCazorla et al. (2009) for a range of temperatures.

Except for trivial space-filling structures, there are no provably densest packings when the spheres possess a size distribution, which is a testament to the mathematical challenges that they present. We begin by noting some rigorous bounds on the maximal packing fraction of packings of spheres with different radii in . Specifically, the overall maximal packing fraction of such a general mixture in [where is defined by (19) with (20)] is bounded from above and below in terms of the maximal packing fraction for a monodisperse packing in the infinite-volume limit:Torquato (2002)


The lower bound corresponds to the case when the components completely demix, each at the density . The upper bound corresponds to an ideal sequential packing process for arbitrary in which one takes the limits , .Torquato (2002) Specific nonsequential protocols (algorithmic or otherwise) that can generate structures that approach the upper bound (21) for arbitrary values of are currently unknown and thus the development of such protocols is an open area of research. We see that in the limit , the upper bound approaches unity, corresponding to space-filling polydisperse spheres with an infinitely wide separation in sizes Herrmann et al. (1990) or a continuous size distribution with sizes ranging to the infinitesimally small.Torquato (2002)

Even the characterization of the densest jammed binary sphere packings offers great challenges and is very far from completion. Here we briefly review some work concerning such packings in two and three dimensions. Ideally, it is desired to obtain as a function of and . In practice, we have a sketchy understanding of the surface defined by .

Among the 2D and 3D cases, we know most about the determination of the maximally dense binary packings in . Fejes Tóth Fejes Tóth (1964) reported a number of candidate maximally dense packings for certain values of the radii ratio in the range . Maximally dense binary disk packings have been also studied to determine the stable crystal phase diagram of such alloys.Likos and Henley (1993) The determination of for sufficiently small amounts to finding the optimal arrangement of the small disks within a tricusp: the nonconvex cavity between three close-packed large disks. A particle-growth Monte Carlo algorithm was used to generate the densest arrangements of a fixed number of small identical disks within such a tricusp. Uche et al. (2004) All of these results can be compared to a relatively tight upper bound on given by Florian (1960)


Inequality (22) also applies to general multicomponent packings, where is taken to be the ratio of the smallest disk radius to the largest one.

There is great interest in finding the densest binary sphere packings in in the physical sciences because of their relationship to binary crystal alloys found in molecular systems; see Refs. Hopkins et al., 2012a and Hudson and Harrowell, 2008 as well as references therein for details and some history. Past efforts to identify such optimal packings have employed simple crystallographic techniques (filling the interstices in uniform 3D tilings of space with spheres of different sizes) and algorithmic methods, e.g., Monte Carlo calculations and a genetic algorithm.Kummerfeld et al. (2008); Filion and Dijkstra (2009) However, these methods have achieved only limited success, in part due to the infinite parameter space that is involved. When using traditional algorithms, difficulties result from the enormous number of steps required to escape from local minima in the “energy” (negative of the packing fraction). Hopkins et al. Hopkins et al. (2011a, 2012a) have presented the most comprehensive determination to date of the “phase diagram” for the densest binary sphere packings via the TJ linear-programming algorithm.Torquato and Jiao (2010a) In Ref. Hopkins et al., 2012a, 19 distinct crystal alloys (compositions of large and small spheres spatially mixed within a fundamental cell) were identified, including 8 that were unknown at the time. Using the TJ algorithm, they were always able to obtain either the densest previously known alloy or denser ones. These structures may correspond to currently unidentified stable phases of certain binary atomic and molecular systems, particularly at high temperatures and pressures.Demchyna et al. (2006); Cazorla et al. (2009) Reference Hopkins et al., 2012a provides details about the structural characteristics of these densest-known binary sphere packings.

Iv Packing Spheres in High Dimensions

Sphere packings in four and higher-dimensional Euclidean spaces are of great interest in the physical and mathematical sciences; see Refs. Conway and Sloane, 1995; Frisch and Percus, 1999; Parisi and Slanina, 2000; Elkies, 2000; Parisi and Zamponi, 2006; Cohn, 2002; Cohn and Elkies, 2003; Torquato and Stillinger, 2006b, a; Torquato, 2006; Skoge et al., 2006; Rohrmann and Santos, 2007; Adda-Bedia et al., 2008; Scardicchio et al., 2008; Cohn and Kumar, 2009; van Meel et al., 2009, 2009; Parisi and Zamponi, 2010; Lue et al., 2010; Torquato and Stillinger, 2010; Marcotte and Torquato, 2013; Kallus et al., 2013; Kallus and Torquato, 2014; Andreanov et al., 2016. Physicists have studied high-dimensional packings to gain insight into liquid, crystal and glassy states of matter in lower dimensions.Frisch and Percus (1999); Parisi and Slanina (2000); Parisi and Zamponi (2006); Skoge et al. (2006); Rohrmann and Santos (2007); Adda-Bedia et al. (2008); van Meel et al. (2009, 2009); Parisi and Zamponi (2010) Finding the densest packings in arbitrary dimension in Euclidean and compact spaces is a problem of long-standing interest in discrete geometry.Conway and Sloane (1998); Cohn and Elkies (2003); Viazovska (2017); Cohn et al. (2017) Remarkably, the optimal way of sending digital signals over noisy channels corresponds to the densest sphere packing in a high-dimensional space.Shannon (1948); Conway and Sloane (1998) These “error-correcting” codes underlie a variety of systems in digital communications and storage, including compact disks, cell phones and the Internet.

iv.1 Equilibrium and Metastable Phase Behavior

The properties of equilibrium and metastable states of hard spheres have been studied both theoretically and computationally in spatial dimensions greater than three. Luban and Baram (1982); Frisch and Percus (1999); Skoge et al. (2006); Clisby and McCoy (2006); van Meel et al. (2009, 2009); Lue et al. (2010) This includes the evaluation of various virial coefficients across dimensions, Luban and Baram (1982); Clisby and McCoy (2006) as well as the pressure along the liquid, metastable and crystal branches,Skoge et al. (2006); Rohrmann and Santos (2007); van Meel et al. (2009, 2009); Lue et al. (2010) especially for and 7. Using the LS packing algorithm, Skoge et al. Skoge et al. (2006) numerically estimated the freezing and melting packing fractions for the equilibrium hard-sphere fluid-solid transition, and , respectively, for , and and , respectively, for . These authors showed that nucleation appears to be strongly suppressed with increasing dimension. The same conclusion was subsequently reached in a study by van Meel et al..van Meel et al. (2009) Finken, Schmidt and Löwen Finken et al. (2001) used a variety of approximate theoretical methods to show that equilibrium hard spheres have a first-order freezing transition for dimensions as high as .

Any disordered packing in which the pair correlation function at contact, , is bounded (such as equilibrium hard spheres) induces a power-law decay in the structure factor in the limit for any dimension given by


There is a remarkable duality between the equilibrium hard-hypersphere (hypercube) fluid system in and the continuum percolation model of overlapping hyperspheres (hypercubes) in . In particular, the pair connectedness function of the latter is to a good approximation equal to the negative of the total correlation function of the former evaluated at negative density.Torquato (2012) This mapping becomes exact for and in the large- limit.

iv.2 Nonequilibrium Disordered Packings Via Sequential Addition

In Sec. III.2.2, we noted that the “ghost” RSA packing Torquato and Stillinger (2006b) is a disordered but unsaturated packing construction whose -particle correlation functions are known exactly and rigorously achieves the infinite-time packing fraction for any ; see Fig. 3 for a 2D realization of such a packing.

Saturated RSA sphere packings have been numerically generated and structurally characterized for dimensions up through (Ref. Torquato et al., 2006). A more efficient numerical procedure was devised to produce such packings for dimensions up through (Ref. Zhang and Torquato, 2013). The current best estimates of the maximal saturation packing fraction for , , , and are , , , , and , respectively.Zhang and Torquato (2013) These are lower than the corresponding MRJ packing fractions in those dimensions (see Sec. IV.3). The quantity apparently scales as or possibly for large ; see Refs. Torquato et al., 2006; Zhang and Torquato, 2013. While saturated RSA packings are nearly but not exactly hyperuniform,Torquato et al. (2006) as increases, the degree of hyperuniformity increases and pair correlations markedly diminish,Zhang and Torquato (2013) consistent with the decorrelation principle,Torquato and Stillinger (2006a) which is described in Sec. IV.5.

iv.3 Maximally Random Jammed States

Using the LS algorithm, Skoge et al. Skoge et al. (2006) generated and characterized MRJ packings in four, five and six dimensions. In particular, they estimated the MRJ packing fractions, finding , and for , and , respectively. To a good approximation, the MRJ packing fraction obeys the scaling form , where and , which appears to be consistent with high-dimensional asymptotic limit, albeit with different coefficients. The dominant large- density scaling is supported by theoretical studies.Torquato and Stillinger (2002, 2006a); Jin et al. (2010); Kallus and Torquato (2014) Skoge et al. Skoge et al. (2006) also determined the MRJ pair correlation function and structure factor for these states and found that short-range ordering appreciably decreases with increasing dimension, consistent with the decorrelation principle.Torquato and Stillinger (2006a) This implies that, in the limit , tends to unity for all outside the hard-core, except for a Dirac-delta function at contact due to the jamming constraint.Torquato and Stillinger (2006a) As for (where ), the MRJ packings were found to be isostatic, hyperuniform, and have a power-law divergence in at contact; ; with as tends to . Across dimensions, the cumulative number of neighbors was shown to equal the kissing (contact) number of the conjectured densest packing close to where has its first minimum. Disordered jammed packings were also simulated and studied in dimensions 7-10; see Ref. Charbonneau et al., 2012.

iv.4 Maximally Dense Sphere Packings

The sphere packing problem seeks to answer the following question: Conway and Sloane (1998) Among all packings of congruent spheres in , what is the maximal packing fraction and the corresponding arrangements of the spheres? Until 2017, the optimal solutions were known only for the first three space dimensions. Hales (2005) For and , these are the triangular lattice () with and checkerboard (fcc) lattice () and its stacking variants with , respectively. We now know the optimal solutions in two other space dimensions; namely, the and lattices are the densest packings among all possible packings in and , respectively; see Refs. Viazovska, 2017 and Cohn et al., 2017. For , the densest known packings are (Bravais) lattice packings. Conway and Sloane (1998) The “checkerboard” lattice is believed to be optimal in and . Interestingly, the non-lattice (periodic) packing (with 40 spheres per fundamental cell) is the densest known packing in , which is the lowest dimension in which the best known packing is not a (Bravais) lattice.

Table 1 lists the densest known or optimal sphere packings in for selected . For the first three space dimensions, the optimal sphere-packing solutions (or their “dual” solutions) are directly related to the best known solutions of the number-variance problem Torquato and Stillinger (2003); Zachary and Torquato (2009) as well as of two other well-known problems in discrete geometry: the covering and quantizer problems, Conway and Sloane (1998); Sikirić et al. (2008) but such relationships may or may not exist for , depending on the peculiarities of the dimensions involved. Torquato (2010a)

   Packing  Packing Fraction, Kissing number,
1 2
2 6
3 12
4 24
5 40
6 72
7 126
8 240
9 272
10 372
16 4320
24 196360
Table 1: The densest known or optimal sphere packings in for selected . For each packing, we provide the packing fraction and kissing number . Except for the non-lattice packing in , all of the other densest known packings listed in this table are lattice packings: is the integer lattice, is the triangular lattice, is the checkerboard lattice (a generalization of the fcc lattice), is one the root lattices, and is the laminated lattice. For and , it has recently been proved that and are optimal among all packings, respectively; see Refs. Viazovska, 2017 and Cohn et al., 2017. The reader is referred to Ref. Conway and Sloane, 1998 for additional details.

The TJ linear-programming packing algorithm was adapted by Marcotte and Torquato Marcotte and Torquato (2013) to determine the densest lattice packing (one particle per fundamental cell) in some high dimension. These authors applied it for and showed that it was able to rapidly and reliably discover the densest known lattice packings without a priori knowledge of their existence. It was found to be appreciably faster than previously known algorithms at that time.Gravel et al. (2011); Andreanov and Scardicchio (2012) The TJ algorithm was used to generate an ensemble of isostatic jammed hard-sphere lattices and study the associated pair statistic and force distributions.Kallus et al. (2013) It was shown that this special ensemble of lattice-sphere packings retain many of the crucial structural features of the classical hard-sphere model.

It is noteworthy that for sufficiently large , lattice packings are most likely not the densest (see Fig. 9), but it becomes increasingly difficult to find explicit dense packing constructions as increases. Indeed, the problem of finding the shortest lattice vector in a particular lattice packing grows super-exponentially with and is in the class of NP-hard problems. Ajtai (1998)

Figure 9: (color online) Lattice packings in sufficiently high dimensions are not dense because the “holes” (space exterior to the spheres) eventually dominates the space and hence become unsaturated.Conway and Sloane (1998) For illustration purposes, we consider the hypercubic lattice . Left panel: A fundamental cell of represented in two dimensions. The distance between the point of intersection of the longest diagonal in the hypercube with the hypersphere boundary and the vertex of the cube along this diagonal is given by for a sphere of unit radius. This means that already becomes unsaturated at . Placing an additional sphere in doubles the density of and yields the four-dimensional checkerboard lattice packing , which is believed to be the optimal packing in . Right panel: A schematic “effective” distorted representation of the hypersphere within the hypercubic fundamental cell for large , illustrating that the volume content of the hypersphere relative to the hypercube rapidly diminishes asymptotically. Indeed, the packing fraction of is given by . Indeed, the checkerboard lattice with packing fraction becomes suboptimal in relatively low dimensions because it also becomes dominated by larger and larger holes as increases.

For large , the best that one can do theoretically is to devise upper and lower bounds on . Conway and Sloane (1998) The nonconstructive lower bound of Minkowski Minkowski (1905) established the existence of reasonably dense lattice packings. He found that the maximal packing fraction among all lattice packings for satisfies


where is the Riemann zeta function. Note that for large values of , the asymptotic behavior of the Minkowski lower bound is controlled by . Since 1905, many extensions and generalizations of the lower bound (24) have been derived Davenport and Rogers (1947); Ball (1992); Conway and Sloane (1998); Vance (2009), but none of these investigations have been able to improve upon the dominant exponential term ; they instead only improve on the latter by a factor linear in .

It is trivial to prove that the packing fraction of a saturated packing of congruent spheres in satisfies Torquato and Stillinger (2006b)


for all . This “greedy” bound (25) has the same dominant exponential term as Minkowski’s bound (24).

Nontrivial upper bounds on in for any have been derived.Blichfeldt (1929); Rogers (1958, 1964); Kabatiansky and Levenshtein (1978); Cohn and Elkies (2003) The linear-programming (LP) upper bounds due to Cohn and Elkies Cohn and Elkies (2003) provides the basic framework for proving the best known upper bounds on for dimensions in the range . Recently, these LP bounds have been used to prove that no packings in and have densities that can exceed those of the (Ref. Viazovska, 2017) and (Ref. Cohn et al., 2017) lattices, respectively. Kabatiansky and Levenshtein Kabatiansky and Levenshtein (1978) found the best asymptotic upper bound, which in the limit yields , proving that the maximal packing fraction tends to zero in the limit . This rather counterintuitive high-dimensional property of sphere packings can be understood by recognizing that almost all of the volume of a -dimensional sphere for large is concentrated near the sphere surface.

iv.5 Are Disordered Packings the Densest in High Dimensions?

Since 1905, many extensions and generalizations of Minkowski’s bound have been derived Conway and Sloane (1998), but none of them have improved upon the dominant exponential term . The existence of the unjammed disordered ghost RSA packing Torquato and Stillinger (2006b) (Sec. III.2.2) that rigorously achieves a packing fraction of strongly suggests that Bravais-lattice packings (which are almost surely unsaturated in sufficiently high ) are far from optimal for large .

Torquato and Stillinger Torquato and Stillinger (2006a) employed a plausible conjecture that strongly supports the counterintuitive possibility that the densest sphere packings for sufficiently large may be disordered or at least possess fundamental cells whose size and structural complexity increase with . They did so using the so-called -invariant optimization procedure that maximizes associated with a radial “test” pair correlation function to provide the putative exponential improvement on Minkowski’s 100-year-old bound on . Specifically, a -invariant process Torquato and Stillinger (2002) is one in which the functional form of a “test” pair correlation function remains invariant as density varies, for all , over the range of packing fractions subject to the satisfaction of the nonnegativity of the structure factor and . When there exist sphere packings with a satisfying these conditions in the interval , then one has the lower bound on the maximal packing fraction given by


Torquato and Stillinger Torquato and Stillinger (2006a) conjectured that a test function is a pair-correlation function of a translationally invariant disordered sphere packing in for for sufficiently large if and only if the nonnegativity conditions on and . The decorrelation principle,Torquato and Stillinger (2006a) among other results,Scardicchio et al. (2008); Torquato and Stillinger (2010); Andreanov et al. (2016) provides justification for the conjecture. This principle states that unconstrained correlations in disordered sphere packings vanish asymptotically in high dimensions and that the for any can be inferred entirely (up to small errors) from a knowledge of and . This is vividly exhibited by the exactly solvable ghost RSA packing process Torquato and Stillinger (2006b) as well as by computer simulations of high-dimensional MRJ Skoge et al. (2006) and RSA Torquato (2006) packings. Interestingly, this optimization problem is the dual of the infinite-dimensional linear program (LP) devised by Cohn and Elkies Cohn (2002) to obtain upper bounds on .

Using a particular test pair correlation corresponding to a disordered sphere packing, Torquato and Stillinger Torquato and Stillinger (2006a) found a conjectural lower bound on that is controlled by , thus providing the first putative exponential improvement on Minkowski’s lower bound (24). Scardicchio, Stillinger and Torquato Scardicchio et al. (2008) studied a wider class of test functions (corresponding to disordered packings) that lead to precisely the same putative exponential improvement on Minkowski’s lower bound and therefore the asymptotic form is much more general and robust than previously surmised.

Zachary and Torquato Zachary and Torquato (2011) studied, among other quantities, pair statistics of high-dimensional generalizations of the periodic 2D kagomé and 3D diamond crystal structures. They showed that the decorrelation principle is remarkably already exhibited in these periodic crystals in low dimensions, suggesting that it applies for any sphere packing in high dimensions, whether disordered or not. This conclusion was bolstered in a subsequent study by Andreanov, Scardicchio and TorquatoAndreanov et al. (2016) who showed that strictly jammed lattice sphere packings visibly decorrelate as increases in relatively low dimensions.

iv.6 Remarks on Packing Problems in Non-Euclidean Spaces

Particle packing problems in non-Euclidean (curved) spaces arises in a variety of fields, including physics, Renes (2005); Bowick et al. (2006); Modes and Kamien (2007) biology, Tammes (1930); Goldberg (1967); Prusinkiewicz and Lindenmayer (1990); Torquato et al. (2002); Zandi et al. (2004) communications theory, Conway and Sloane (1998) and geometry. Conway and Sloane (1998); Bowen and Radin (2003); Hardin and Saff (2004); Cohn and Kumar (2007); Cohn et al. (2011) While a comprehensive overview of this topic is beyond the scope of this review, it is useful to highlight here some of the developments for the interested reader. We will limit the discussion to sphere packings on the positively curved unit sphere . The reader is referred to the review by Torquato and Stillinger Torquato and Stillinger (2010) for some discussion of negatively curved hyperbolic space .

Recall that the kissing number is the number of spheres of unit radius that simultaneously touch a unit sphere . Conway and Sloane (1998) The kissing number problem asks for the maximal kissing number in . The determination of the maximal kissing number in spurred a famous debate between Issac Newton and David Gregory in 1694. The former correctly thought the answer was 12, but the latter wrongly believed that it was 13. The maximal kissing number for is only known in dimensions four, Musin (2008) eight,Levenshtein (1979); Odlyzko and Sloane (1979) and twenty four.Levenshtein (1979); Odlyzko and Sloane (1979)

A packing of congruent spherical caps on the unit sphere in yields a spherical code consisting of the centers of the caps. Conway and Sloane (1998) A spherical code is optimal if the minimal distance (smallest angular separation between distinct points in the code) is as large as possible. The reader is referred to the paper by Cohn and Kumar Cohn and Kumar (2007) and references therein for some developments in the mathematics literature. Interestingly, the jamming of spherical codes in a variety of dimensions has been investigated. Cohn et al. (2011) It is noteworthy that some optimal spherical codes Hopkins et al. (2010a) are related to the densest local packing of spheres around a central sphere. Hopkins et al. (2010b, 2011b)

V Packings of Nonspherical Particles

The packing characteristics of equilibrium phases and jammed states of packings of nonspherical particle are considerably richer than their spherical counterparts. Onsager (1949); Hoylman (1970); Frenkel and Maguire (1983); Bezdek and Kuperberg (1991); Allen (1993); Frenkel (1999); Betke and Henk (2000); Williams and Philipse (2003); Donev et al. (2004b, d); Man et al. (2005); Donev et al. (2005c, d); Bettina and Escobedo (2005); Conway and Torquato (2006); Donev et al. (2007a); Bolhuis and Frenkel (1997); Chen (2008); Torquato and Jiao (2009a, b); Haji-Akbari et al. (2009); Jiao et al. (2010); Torquato (2010b); Kallus et al. (2010); Chen et al. (2010); Jaoshvili et al. (2010); Batten et al. (2010); Jiao and Torquato (2011a, b); Agarwal and Escobedo (2011); Gravel et al. (2011); Marechal et al. (2011); Torquato and Jiao (2012); Henzie et al. (2012); Damasceno et al. (2012a, b); Atkinson et al. (2012); Ni et al. (2012); Gabbrielli et al. (2012); Bautista-Carbajal et al. (2013); Gantapara et al. (2013); Gabbrielli et al. (2014); Chen et al. (2014a); Zong (2014); Chen et al. (2014b, a); Cinacchi and Torquato (2015); Dostert et al. (2017); VanderWerf et al. (2018) This is due to the fact that nonsphericity of the particle shape introduces rotational degrees of freedom not present in sphere packings. Our primary interest is in simple 3D convex shapes, such as ellipsoids, superballs, spherocylinders and polyhedra, although we remark on more complex shapes, such as concave particles as well as congruent ring tori, which are multiply connected nonconvex bodies of genus 1. Organizing principles to characterize and classify very dense, possibly jammed packings of nonspherical particles in terms of shape symmetry of the particles Torquato (2009); Torquato and Jiao (2009a, 2012) are discussed in Sec.V.4.5.

v.1 Simple Nonspherical Convex Shapes

v.1.1 Ellipsoids

One simple generalization of the sphere is an ellipsoid, the family of which is a continuous deformation of a sphere. A -dimensional ellipsoid in is a centrally symmetric body occupying the region


where are Cartesian coordinates and are the semi-axes of the ellipsoid. Thus, we see that an ellipsoid is an affine (linear) transformation of the sphere.

v.1.2 Superballs

A -dimensional superball in is a centrally symmetric body occupying the region


where are Cartesian coordinates and is the deformation parameter (not pressure, as denoted in Sec. III.1), which controls the extent to which the particle shape has deformed from that of a -dimensional sphere (). Thus, superballs constitute a large family of both convex () and concave () particles (see Fig. 10). A “superdisk”,which is the designation in the 2D case, possesses square symmetry. As moves away from unity, two families of superdisks with square symmetry can be obtained depending on whether or . When , the superdisk is concave; see Ref. Jiao et al., 2008.

Figure 10: (color online). Superballs with different values of the deformation parameter . We note that and correspond to a 3D cross, regular octahedron, sphere and cube, respectively.

v.1.3 Spherocylinder

A -dimensional spherocylinder in consists of a cylinder of length and radius capped at both ends by hemispheres of radius , and therefore is a centrally-symmetric convex particle. Its volume for is given by


When , a -dimensional spherocylinder reduces to a -dimensional sphere of radius . Figure 11 shows 3D examples.

Figure 11: (Color online) Three-dimensional spherocylinders composed of cylinder with length , caped at both ends with hemispheres with radius . Left panel: A spherocylinder with aspect ratio . Right panel: A spherocylinder with aspect ratio .

v.1.4 Polyhedra

Figure 12: (color online). The five Platonic solids: tetrahedron (P1), icosahedron (P2), dodecahedron (P3), octahedron (P4), and cube (P5).

The Platonic solids (mentioned in Plato’s Timaeus) are convex polyhedra with faces composed of congruent convex regular polygons. There are exactly five such solids: the tetrahedron (P1), icosahedron (P2), dodecahedron (P3), octahedron (P4), and cube (P5) (see Fig. 12) Note that viral capsids often have icosahedral symmetry; see, for example, Ref. Zandi et al., 2004.

Figure 13: (color online) The thirteen Archimedean solids: truncated tetrahedron (A1), truncated icosahedron (A2), snub cube (A3), snub dodecahedron (A4), rhombicosidodecahdron (A5), truncated icosidodecahedron (A6), truncated cuboctahedron (A7), icosidodecahedron (A8), rhombicuboctahedron (A9), truncated dodecahedron (A10), cuboctahedron (A11), truncated cube (A12), and truncated octahedron (A13).

An Archimedean solid is a highly symmetric, semi-regular convex polyhedron composed of two or more types of regular polygons meeting in identical vertices. The thirteen Archimedean solids are depicted in Fig. 13. This typical enumeration does not count the chiral forms (not shown) of the snub cube (A3) and snub dodecahedron (A4). The remaining 11 Archimedean solids are non-chiral (i.e., each solid is superposable on its mirror image) and the only non-centrally symmetric one among these is the truncated tetrahedron.

It is noteworthy that the tetrahedron (P1) and the truncated tetrahedron (A1) are the only Platonic and non-chiral Archimedean solids, respectively, that are not centrally symmetric. We will see that the central symmetry of the majority of the Platonic and Archimedean solids (P2 – P5, A2 – A13) distinguish their dense packing arrangements from those of the non-centrally symmetric ones (P1 and A1) in a fundamental way.

v.2 Equilibrium and Metastable Phase Behavior

Hard nonspherical particles exhibit a richer phase diagram than that of hard spheres because the former can possess different degrees of translational and orientational order; see Refs. Onsager, 1949; Frenkel and Maguire, 1983; Allen, 1993; Frenkel, 1999; Bautista-Carbajal et al., 2013; Cinacchi and Torquato, 2015; Bolhuis and Frenkel, 1997; Batten et al., 2010; Gabbrielli et al., 2012; Bettina and Escobedo, 2005; Agarwal and Escobedo, 2011; Gantapara et al., 2013; Damasceno et al., 2012a, b; Jiao and Torquato, 2011a; Chen et al., 2014a. Nonspherical-particle systems can form isotropic liquids, a variety of liquid-crystal phases, rotator crystals, and solid crystals. Particles in liquid phases have neither translational nor orientational order. Examples of liquid-crystal states include a nematic phase in which the particles are aligned (i.e., with orientational order) while the system lacks any long-range translational order and a smectic phase in which the particles have ordered orientations and possess translational order in one direction. A rotator (or plastic) phase is one in which particles possess translational order but can rotate freely. Solid crystals are characterized by both translational and orientational order.

Ordering transitions in systems of hard nonspherical particles are entropically driven, i.e., the stable phase is determined by a competition between translational and orientational entropies. This principle was first established in the pioneering work of Onsager,Onsager (1949) where it was shown that needle-like shapes exhibit a liquid-nematic phase transition at low densities because in the nematic phase the drop in orientational entropy is offset by the increase in translational entropy, i.e., the available space for any needle increases as the needle tends to align.

The stable phase formed by systems of hard-nonspherical particles is influenced by its symmetry and surface smoothness, which in turn determines the overall system entropy. Here we briefly highlight 3D numerical studies that use Monte Carlo methods and/or free-energy calculations to determine the phase diagrams. Spheroids exhibit not only fluid, solid crystal, and nematic phases for some aspect ratios, but rotator phases for nearly spherical particle shapes at intermediate densities.Frenkel et al. (1984); Bautista-Carbajal et al. (2013) Hard “lenses” (common volume to two overlapping identical spheres) have a qualitatively similar phase diagram to hard oblate spheroids, but differences between them are more pronounced in the high-density crystal phase up to the densest-known packings. Cinacchi and Torquato (2015) Spherocylinders exhibit five different possible phases, depending on the packing fraction and aspect ratio: isotropic fluid, smectic, nematic, rotator, and solid crystal phases.Bolhuis and Frenkel (1997) Apart from fluid and crystal phases, superballs can form rotator phases at intermediate densities.Batten et al. (2010); Gabbrielli et al. (2012) Systems of tetragonal parallelepipeds can exhibit liquid-crystalline and cubatic phases, depending on the aspect ratio and packing fraction.Bettina and Escobedo (2005) Various convex space-filling polyhedra (including truncated octahedron, rhombic dodecahedron, two types of prisms and cube) possess unusual liquid-crystalline and rotator-crystalline phases at intermediate packing fractions.Agarwal and Escobedo (2011) Truncated cubes exhibit a rich diversity in crystal structures that depend sensitively on the amount of truncation.Gantapara et al. (2013) A study of a wide class of polyhedra revealed that entropy maximization favors mutual alignment of particles along their facets. Damasceno et al. (2012b) By analytically constructing the densest known packings of congruent Archimedean truncated tetrahedra (which nearly fill all of space), the melting properties of such systems were examined by decompressing this densest packing and equilibrating them.Jiao and Torquato (2011a) A study of the entire phase diagram of truncated tetrahedra showed that the system undergoes two first-order phase transitions as the density increases: first a liquid-solid transition and then a solid-solid transition.Chen et al. (2014a)

v.3 Maximally Random Jammed States

The fact that M&M candies (spheroidal particles with aspect ratio ) were shown experimentally to randomly pack more densely than spheres () motivated the development of a modified LS algorithm Donev et al. (2005c, d) to obtain frictionless MRJ-like packings with even higher densities for other aspect ratios. This included nearly-spherical ellipsoids with ,Donev et al. (2004b, 2007a) i.e., packing fractions approaching those of the densest 3D sphere packings. Note that these other densest MRJ packings are realizable experimentally. Man et al. (2005)

Figure 14: Packing fraction versus aspect ratio for MRJ packings of 10,000 ellipsoids as obtained in Ref. Donev et al., 2007a. The semi-axes here are . The inset shows the mean contact number as a function of . Neither the spheroid (oblate or prolate) or general ellipsoids cases attain their isostatic values of or , respectively.

Figure 14 shows separate plots of and mean contact number as a function of as predicted by the more refined simulations of Donev et al.. Donev et al. (2007a) Each plot shows a cusp (i.e., non-differentiable) minimum at the sphere point, and versus aspect ratio possesses a density maximum. The existence of a cusp at the sphere point runs counter to the conventional expectation that for “generic” (disordered) jammed frictionless particles, the total number of (independent) constraints equals the total number of degrees of freedom , which has been referred to as the isostatic conjecture. Alexander (1998) This conjecture implies a mean contact number , where for disks, for ellipses, for spheres, for spheroids, and for general ellipsoids. Since increases discontinuously with the introduction of rotational degrees of freedom as one makes the particles nonspherical, the isostatic conjecture predicts that should have a jump increase at aspect ratio to a value of for a general ellipsoid. Such a discontinuity was not originally observed by Donev et al. Donev et al. (2004b), rather, they found that jammed ellipsoid packings are hypostatic (or sub-isostatic), , near the sphere point, and only become nearly isostatic for large aspect ratios. In fact, the isostatic conjecture is only rigorously true for amorphous sphere packings after removal of rattlers; generic nonspherical-particle packings should generally be hypostatic. Roux (2000); Donev et al. (2007a) It has been rigorously shown that packings of nonspherical particles are generally not jammed to first order in the “jamming gap” [see Sec. III.3.2], but are jammed to second order in due to curvature deviations from the sphere.Donev et al. (2007a) In striking contrast with MRJ sphere packings, the rattler concentrations of the MRJ ellipsoid packings appear practically to vanish outside of some small neighborhood of the sphere point.Donev et al. (2007a) The reader is referred to recent work on hypostatic jammed 2D packings of noncircular particles.VanderWerf et al. (2018)

Jiao, Stillinger and Torquato Jiao et al. (2010) have computed the packing fraction of MRJ packings of binary superdisks in and identical superballs in . They found that increases dramatically and nonanalytically as one moves away from the circular-disk or sphere point (). Moreover, these disordered packings were shown to be hypostatic. Hence, the local particle arrangements are necessarily nontrivially correlated to achieve strict jamming and hence “nongeneric,” the degree of which was quantified. MRJ packings of binary superdisks and of ellipses are effectively hyperuniform.Zachary et al. (2011c) Note that the geometric-structure approach was used to derive a highly accurate formula for the packing fraction of MRJ binary packings of convex superdisksTian et al. (2015) that is valid for almost all size ratios, relative concentrations and deformation parameter . For the special limit of monodisperse circular disks, this formula predicts , which is in very good agreement with the recently numerically discovered MRJ isostatic stateAtkinson et al. (2014) with .

3D MRJ packings of the four non-tiling Platonic solids (tetrahedra, octahedra, dodecahedra, and icosahedra) were generated using the ASC optimization scheme. Jiao and Torquato (2011b) The MRJ packing fractions for tetrahedra, octahedra, dodecahedra, and icosahedra are , , , and , respectively. It was shown that as the number of facets of the particles increases, the translational order in the packings increases while the orientational order decreases. Moreover, such MRJ packings were found to be hyperuniform with a total correlation function that decays to zero asymptotically with same power-law as MRJ spheres, i.e, like . These results suggest that hyperuniform quasi-long-range correlations are a universal feature of MRJ packings of frictionless particles of general shape. However, unlike MRJ packings of ellipsoids, superballs, and superellipsoids, which are hypostatic, MRJ packings of the non-tiling Platonic solids are isostatic.Jiao and Torquato (2011b) In addition, 3D MRJ packings of truncated tetrahedra with an average packing fraction of 0.770 were also generated. Chen et al. (2014a)

v.4 Maximally Dense Packings

We focus here mainly on exact constructions of the densest known packings of nonspherical particles in which the lattice vectors as well as particle positions and orientations are expressible analytically. Where appropriate, we also cite some work concerning packings obtained via computer simulations and laboratory experiments.

Rigorous upper bounds on the maximal packing fraction of packings of nonspherical particles of general shape in can be used to assess their packing efficiency. However, it has been highly challenging to formulate upper bounds for non-tiling particle packings that are nontrivially less than unity. It has recently been shown that of a packing of congruent nonspherical particles of volume in is bounded from above according to


where is the volume of the largest sphere that can be inscribed in the nonspherical particle and is the maximal packing fraction of a packing of -dimensional of identical spheres Torquato and Jiao (2009a, b) (e.g., for and for ). The upper bound (30) will be relatively tight provided that the asphericity (equal to the ratio of the circumradius to the inradius) of the particle is not large. Since bound (30) cannot generally be sharp (i.e., exact) for a non-tiling particle, any packing whose density is close to the upper bound (30) is nearly optimal, if not optimal. Interestingly, a majority of the centrally symmetric Platonic and Archimedean solids have relatively small asphericities, explaining the corresponding small differences between and , the packing fraction of the densest lattice packing. Minkowski (1905); Betke and Henk (2000); Torquato and Jiao (2009a, b) As discussed below, these densest lattice packings are conjectured to be the densest among all packings.Torquato and Jiao (2009a, b) Upper bound (30) will also be relatively tight for superballs (superdisks) for deformation parameters near the sphere (circle) point (). Dostert et al. Dostert et al. (2017) recently obtained upper bounds on of translative packings of superballs and 3D convex bodies with tetrahedral symmetry.

v.4.1 Ellipsoids

Figure 15: The packing fraction of the “laminated” non-Bravais lattice packing of ellipsoids (with a two-particle basis) as a function of the aspect ratio . The point corresponding to the face-centered cubic lattice sphere packing is shown, along with the two sharp maxima in the packing fraction for prolate ellipsoids with and oblate ellipsoids with , as shown in the insets. For both and